The patrol problem is given by\begin{align*} u_{i} & \geq 0\\ 2u_{1}+2u_{2} & \leq 0\\ 2u_{1}+2u_{2} & \geq 0\\ u_{2} & \geq 1.5u_{1} \end{align*}
And the objective function which we want to minimize is \[ J\left ( u\right ) =\frac{1}{30}u_{1}+\frac{1}{15}u_{2}\] The above is the raw LP. We convert it to standard LP by introducing slack and surplus variable and rename the variables to \(x_{i}\) from \(u_{i}\). Therefore the first table of the initial phase is
| \(x_{1}\) | \(x_{2}\) | \(x_{3}\) | \(x_{4}\) | \(x_{5}\) | \(b\) | |
| row 1 | \(2\)
| \(2\)
| \(1\)
| \(0\)
| \(0\)
| \(10\) |
| row 2 | \(2\) | \(2\) | \(0\) | \(-1\) | \(0\) | \(4\)
|
| row 3 | \(-1.5\)
| \(1\)
| \(0\)
| \(0\)
| \(-1\)
| \(0\)
|
| \(J\left ( x\right ) \) | \(\frac{1}{30}\) | \(\frac{1}{15}\) | \(0\) | \(0\) | \(0\) | \(0\) |
Since there are two surplus variables (these are the ones associated with \(-1\) entries), we have to start with phase one LP. If there were no surplus variables (i.e. only slack variables), then we go directly to phase two. Phase one is only needed when there are surplus variables. We introduce two new artificial variables \(y_{1},y_{2}\) and an artificial objective \(J\left ( y\right ) \) function which we want to minimize over \(y\) to zero, \[ \min _{y\geq 0,x\geq 0}J\left ( y\right ) =y_{1}+y_{2}\] The first table of first phase is (where we now use the artificial objective function \(J\left ( y\right ) \) in place of \(J\left ( x\right ) \).)
| \(x_1\) | \(x_2\) | \(x_3\) | \(x_4\) | \(x_5\) | \(y_1\) | \(y_2\) | \(b\) | |
| row 1 | \(2\)
| \(2\)
| \(1\)
| \(0\)
| \(0\)
| \(0\)
| \(0\)
| \(10\) |
| row 2 | \(2\) | \(2\) | \(0\) | \(-1\) | \(0\) | \(1\) | \(0\) | \(4\)
|
| row 3 | \(-1.5\)
| \(1\)
| \(0\)
| \(0\)
| \(-1\)
| \(0\)
| \(1\)
| \(0\)
|
| row 4, \(J\left ( y\right ) \) | \(0\) | \(0\) | \(0\) | \(0\) | \(0\) | \(1\) | \(1\) | \(0\) |
We start by making last row canonical (this means we need to zero out last row entries under \(y_{1},y_{2}\) columns). Doing row(4) = row(4)-row(2) gives
| \(x_{1}\) | \(x_{2}\) | \(x_3\) | \(x_{4}\) | \(x_{5}\) | \(y_{1}\) | \(y_{2}\) | \(b\) | |
| row 1 | \(2\)
| \(2\)
| \(1\)
| \(0\)
| \(0\)
| \(0\)
| \(0\)
| \(10\) |
| row 2 | \(2\) | \(2\) | \(0\) | \(-1\) | \(0\) | \(1\) | \(0\) | \(4\)
|
| row 3 | \(-1.5\)
| \(1\)
| \(0\)
| \(0\)
| \(-1\)
| \(0\)
| \(1\)
| \(0\)
|
| row 4, \(J\left ( y\right ) \) | \(-2\) | \(-2\) | \(0\) | \(1\) | \(0\) | \(0\) | \(1\) | \(-4\) |
To zero out last row under \(y_{2}\), row(4)=row(4)-row(3) applied to the above gives
| \(x_{1}\) | \(x_{2}\) | \(x_{3}\) | \(x_{4}\) | \(x_{5}\) | \(y_{1}\) | \(y_{2}\) | \(b\) | |
| row 1 | \(2\)
| \(2\)
| \(1\)
| \(0\)
| \(0\)
| \(0\)
| \(0\)
| \(10\) |
| row 2 | \(2\) | \(2\) | \(0\) | \(-1\) | \(0\) | \(1\) | \(0\) | \(4\)
|
| row 3 | \(-1.5\)
| \(1\)
| \(0\)
| \(0\)
| \(-1\)
| \(0\)
| \(1\)
| \(0\)
|
| row 4, \(J\left ( y\right ) \) | \(-0.5\) | \(-3\) | \(0\) | \(1\) | \(1\) | \(0\) | \(0\) | \(-4\) |
We see that the artificial basic feasible solution \(J\left ( y\right ) \) is not optimal, since there are negative values on the last row. second column has the largest negative value in the last row, at \(-3\). So we need to move this column in the basis vectors. To decide on the pivot row we look at ratio of \(\frac{b}{\text{second\ column}}\) which gives \(\min \left \{ \frac{10}{2},\frac{4}{2},\frac{0}{1}\right \} =0\) which is associated with the third row. So the third row is the pivot row. Now we need to zero out all other entries in the second column. To make (1,2) zero: row 1 = row(1) - \(\left ( 2\times \text{row(3)}\right ) \) gives
| \(x_{1}\) | \(x_{2}\) | \(x_{3}\) | \(x_{4}\) | \(x_{5}\) | \(y_{1}\) | \(y_{2}\) | \(b\) | |
| row 1 | \(5\)
| \(0\)
| \(1\)
| \(0\)
| \(2\)
| \(0\)
| \(-2\)
| \(10\) |
| row 2 | \(2\) | \(2\) | \(0\) | \(-1\) | \(0\) | \(1\) | \(0\) | \(4\)
|
| row 3 (pivot) | \(-1.5\)
| \(1\)
| \(0\)
| \(0\)
| \(-1\)
| \(0\)
| \(1\)
| \(0\)
|
| row 4, \(J\left ( y\right ) \) | \(-0.5\) | \(-3\) | \(0\) | \(1\) | \(1\) | \(0\) | \(0\) | \(-4\) |
To make (2,2) zero, row (2) = row(2) - 2 row(3) gives
| \(x_{1}\) | \(x_{2}\) | \(x_{3}\) | \(x_{4}\) | \(x_{5}\) | \(y_{1}\) | \(y_{2}\) | \(b\) | |
| row 1 | \(5\)
| \(0\)
| \(1\)
| \(0\)
| \(2\)
| \(0\)
| \(-2\)
| \(10\) |
| row 2 | \(5\) | \(0\) | \(0\) | \(-1\) | \(2\) | \(1\) | \(-2\) | \(4\)
|
| row 3 (pivot) | \(-1.5\)
| \(1\)
| \(0\)
| \(0\)
| \(-1\)
| \(0\)
| \(1\)
| \(0\)
|
| row 4, \(J\left ( y\right ) \) | \(-0.5\) | \(-3\) | \(0\) | \(1\) | \(1\) | \(0\) | \(0\) | \(-4\) |
Finally to make (4,2) entry zero, row(4)=row(4)+3 row(3) gives
| \(x_{1}\) | \(x_{2}\) | \(x_{3}\) | \(x_{4}\) | \(x_{5}\) | \(y_{1}\) | \(y_{2}\) | \(b\) | |
| row 1 | \(5\)
| \(0\)
| \(1\)
| \(0\)
| \(2\)
| \(0\)
| \(-2\)
| \(10\) |
| row 2 | \(5\) | \(0\) | \(0\) | \(-1\) | \(2\) | \(1\) | \(-2\) | \(4\)
|
| row 3 (pivot) | \(-1.5\)
| \(1\)
| \(0\)
| \(0\)
| \(-1\)
| \(0\)
| \(1\)
| \(0\)
|
| row 4, \(J\left ( y\right ) \) | \(-5\) | \(0\) | \(0\) | \(1\) | \(-2\) | \(0\) | \(3\) | \(-4\) |
We see that the artificial basic feasible solution is still not optimal since there is negative values on last row. The most negative is in first column. This is the column to move in. Taking the ratio of \(\frac{b}{\text{first\ column}}\) gives \(\min \left \{ \frac{10}{5},\frac{4}{5}\right \} =\frac{4}{5}\) (we do not divide by negative entries). This minimum is associated with second row. So the second row is the pivot row. We start by normalizing the second row (the pivot row) so that entry \(\left ( 2,1\right ) \) is one (it is 5 now). Row(2) = row(2)/5
| \(x_{1}\) | \(x_{2}\) | \(x_{3}\) | \(x_{4}\) | \(x_{5}\) | \(y_{1}\) | \(y_{2}\) | \(b\) | |
| row 1 | \(5\)
| \(0\)
| \(1\)
| \(0\)
| \(2\)
| \(0\)
| \(-2\)
| \(10\) |
| row 2 (pivot) | \(1\) | \(0\) | \(0\) | \(-\frac{1}{5}\) | \(\frac{2}{5}\) | \(\frac{1}{5}\) | \(\frac{-2}{5}\) | \(\frac{4}{5}\)
|
| row 3 | \(-1.5\)
| \(1\)
| \(0\)
| \(0\)
| \(-1\)
| \(0\)
| \(1\)
| \(0\)
|
| row 4, \(J\left ( y\right ) \) | \(-5\) | \(0\) | \(0\) | \(1\) | \(-2\) | \(0\) | \(3\) | \(-4\) |
To make \((1,1)\) entry zero, then row(1)=row(1)-5 row(2)
| \(x_{1}\) | \(x_{2}\) | \(x_{3}\) | \(x_{4}\) | \(x_{5}\) | \(y_{1}\) | \(y_{2}\) | \(b\) | |
| row 1 | \(0\)
| \(0\)
| \(1\)
| \(1\)
| \(0\)
| \(-1\)
| \(0\)
| \(6\) |
| row 2 (pivot) | \(1\) | \(0\) | \(0\) | \(-\frac{1}{5}\) | \(\frac{2}{5}\) | \(\frac{1}{5}\) | \(\frac{-2}{5}\) | \(\frac{4}{5}\)
|
| row 3 | \(-1.5\)
| \(1\)
| \(0\)
| \(0\)
| \(-1\)
| \(0\)
| \(1\)
| \(0\)
|
| row 4, \(J\left ( y\right ) \) | \(-5\) | \(0\) | \(0\) | \(1\) | \(-2\) | \(0\) | \(3\) | \(-4\) |
To make \((1,3)\) zero, then row(3)=row(3)+1.5 row(2)
| \(x_{1}\) | \(x_{2}\) | \(x_{3}\) | \(x_{4}\) | \(x_{5}\) | \(y_{1}\) | \(y_{2}\) | \(b\) | |
| row 1 | \(0\)
| \(0\)
| \(1\)
| \(1\)
| \(0\)
| \(-1\)
| \(0\)
| \(6\) |
| row 2 (pivot) | \(1\) | \(0\) | \(0\) | \(-\frac{1}{5}\) | \(\frac{2}{5}\) | \(\frac{1}{5}\) | \(\frac{-2}{5}\) | \(\frac{4}{5}\)
|
| row 3 | \(0\)
| \(1\)
| \(0\)
| \(-\frac{3}{10}\)
| \(-\frac{2}{5}\)
| \(\frac{3}{10}\)
| \(\frac{2}{5}\)
| \(\frac{6}{5}\)
|
| row 4, \(J\left ( y\right ) \) | \(-5\) | \(0\) | \(0\) | \(1\) | \(-2\) | \(0\) | \(3\) | \(-4\) |
To make \((4,1)\) zero, row(4)=row(4)+5 row(2)
| \(x_{1}\) | \(x_{2}\) | \(x_{3}\) | \(x_{4}\) | \(x_{5}\) | \(y_{1}\) | \(y_{2}\) | \(b\) | |
| row 1 | \(0\)
| \(0\)
| \(1\)
| \(1\)
| \(0\)
| \(-1\)
| \(0\)
| \(6\) |
| row 2 (pivot) | \(1\) | \(0\) | \(0\) | \(-\frac{1}{5}\) | \(\frac{2}{5}\) | \(\frac{1}{5}\) | \(\frac{-2}{5}\) | \(\frac{4}{5}\)
|
| row 3 | \(0\)
| \(1\)
| \(0\)
| \(\frac{-3}{10}\)
| \(\frac{-2}{5}\)
| \(\frac{3}{10}\)
| \(\frac{2}{5}\)
| \(\frac{6}{5}\)
|
| row 4, \(J\left ( y\right ) \) | \(0\) | \(0\) | \(0\) | \(0\) | \(0\) | \(1\) | \(1\) | \(0\) |
We have driven \(J\left ( y\right ) \) to zero with no positive entries in last row. This completes phase one. Now we remove the last row and also remove the \(y_{1},y_{2}\) columns from the above table, and put back the \(J\left ( x\right ) \) in its place in last row. This next tableau is the starting of phase 2.
| \(x_{1}\) | \(x_{2}\) | \(x_{3}\) | \(x_{4}\) | \(x_{5}\) | \(b\) | |
| row 1 | \(0\)
| \(0\)
| \(1\)
| \(1\)
| \(0\)
| \(6\) |
| row 2 | \(1\) | \(0\) | \(0\) | \(-\frac{1}{5}\) | \(\frac{2}{5}\) | \(\frac{4}{5}\)
|
| row 3 | \(0\)
| \(1\)
| \(0\)
| \(\frac{-3}{10}\)
| \(\frac{-2}{5}\)
| \(\frac{6}{5}\)
|
| \(J\left ( x\right ) \) | \(\frac{1}{30}\) | \(\frac{1}{15}\) | \(0\) | \(0\) | \(0\) | \(0\) |
We see that all entries in the last row positive, therefore phase two is now complete. There is nothing to do in phase two. All the hard work was done in phase one. This is special case and we were lucky. Note: The text book says that we should now zero out the last row so that zeros appear under the basis columns. But I find this not needed, since it does not change the optimal \(x^{\ast }\). We can always calculate \(J(x)\) once we know \(x^{\ast }\), and there is no need to have \(J(x)\) show up in the bottom right corner of the tableau really. Therefore, I did not do this extra step as not needed.
Now we read out the solution from the above tableau \[ \boldsymbol{x}=\begin{bmatrix} \frac{4}{5}\\ \frac{6}{5}\\ 6\\ 0\\ 0 \end{bmatrix} =\begin{bmatrix} 0.8\\ 1.2\\ 6\\ 0\\ 0 \end{bmatrix} \] To find the corresponding \(J^{\ast }(u)\), since \(x_1=u_1\) and \(x_2=u_2\) and \(J(u)=\frac{1}{30} u_1+\frac{1}{15}u_2\), then \begin{align*} J^{\ast } & = \frac{1}{30}(0.8) +\frac{1}{15}(1.2) \\ & = 0.10667 \end{align*}
Using Matlab linprog
The output from the above is
Using my own nma_simple.m which prints all the intermediate tableau and solutions \(x\) during the search
The output from the above is
Which gives same answer.
This section shows different views of the problem. In the next section, the solution itself is given. Under each plot, the small code used to generate the plot is shown. First, the feasibility region given by the constraints is plotted.
This plot shows each constraint, superimposed on top of the feasibility region.
The following is 3D plot, showing the four objective functions
The following is same 3D plot, but now with the constraints added
The solution found is now added the above plot., which is \(x_{1}=0,x_{2}=4\) with the min max value of \(J(u)=4\) marked with small red point below.
The first step is to convert this multi-objective minimax problem with constraints, to a pure linear programming problem so that Matlab linprog can be used. Introducing extra variable \(z\) the problem can be written as\[\begin{array} [c]{ccc}\underset{z}{\min } & \tilde{J}\left ( z\right ) =z & \\ s.t. & z\geq J_{i}\left ( x\right ) & i=1\cdots 4\\ & 2x_{1}+2x_{2}\leq 10 & \\ & -2x_{1}-2x_{2}\leq 0 & \\ & x_{2}\geq 1.5x_{1}+4 & \end{array} \] Using the \(J_{i}\left ( x\right ) \) given, the above becomes\[\begin{array} [c]{ccc}\underset{z}{\min } & \tilde{J}\left ( z\right ) =z & \\ s.t. & \frac{1}{30}x_{1}+\frac{1}{15}x_{2}-z\leq 0 & \\ & \frac{3}{30}x_{1}+\frac{1}{5}x_{2}-z\leq 0 & \\ & 2x_{1}+\frac{1}{2}x_{2}-z\leq 0 & \\ & x_{1}+x_{2}-z\leq 0 & \\ & 2x_{1}+2x_{2}\leq 10 & \\ & -2x_{1}-2x_{2}\leq 0 & \\ & x_{2}\geq 1.5x_{1}+4 & \end{array} \] To use Matlab linprog, we first convert the above to the standard form. We do this by introducing slack and surplus variables. The above becomes (added tilde on top of the \(x\) to make it clear which one variables are the original raw LP variables, and which is ones are the new variables).\[\begin{array} [c]{ccc}\underset{z}{\min } & \tilde{J}\left ( z\right ) =z & \\ s.t. & \frac{1}{30}x_{1}+\frac{1}{15}x_{2}-z+\tilde{x}_{3}=0 & \\ & \frac{3}{30}x_{1}+\frac{1}{5}x_{2}-z+\tilde{x}_{4}=0 & \\ & 2x_{1}+\frac{1}{2}x_{2}-z+\tilde{x}_{5}=0 & \\ & x_{1}+x_{2}-z+\tilde{x}_{6}=0 & \\ & 2x_{1}+2x_{2}+\tilde{x}_{7}=10 & \\ & -2x_{1}-2x_{2}+\tilde{x}_{8}=0 & \\ & -1.5x_{1}+x_{2}-\tilde{x}_{9}=4 & \end{array} \] For \(x,\tilde{x}\geq 0\). Hence \(Ax=b\) becomes\[\begin{pmatrix} \frac{1}{30} & \frac{1}{15} & 1 & 0 & 0 & 0 & 0 & 0 & 0 & -1\\ \frac{3}{30} & \frac{1}{5} & 0 & 1 & 0 & 0 & 0 & 0 & 0 & -1\\ 2 & \frac{1}{2} & 0 & 0 & 1 & 0 & 0 & 0 & 0 & -1\\ 1 & 1 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & -1\\ 2 & 2 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0\\ -2 & -2 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0\\ -1.5 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & -1 & 0 \end{pmatrix}\begin{pmatrix} x_{1}\\ x_{2}\\ \tilde{x}_{3}\\ \tilde{x}_{4}\\ \tilde{x}_{5}\\ \tilde{x}_{6}\\ \tilde{x}_{7}\\ \tilde{x}_{8}\\ \tilde{x}_{9}\\ z \end{pmatrix} =\begin{pmatrix} 0\\ 0\\ 0\\ 0\\ 10\\ 0\\ 4 \end{pmatrix} \] And \(\min c^{T}x\) becomes\[ \min \ \begin{pmatrix} 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 \end{pmatrix}\begin{pmatrix} x_{1}\\ x_{2}\\ \tilde{x}_{3}\\ \tilde{x}_{4}\\ \tilde{x}_{5}\\ \tilde{x}_{6}\\ \tilde{x}_{7}\\ \tilde{x}_{8}\\ \tilde{x}_{9}\\ z \end{pmatrix} \] Therefore, using the above, we can write the first tableau table. Then use Matlab to solve it.
\[\begin{tabular} [c]{|c||cccccccccc||c||}\hline & $x_{1}$ & $x_{2}$ & $x_{3}$ & $x_{4}$ & $x_{5}$ & $x_{6}$ & $x_{7}$ & $x_{8}$ & $x_{9}$ & $z$ & $b$\\\hline \hline row 1 & $\frac{1}{30}$ & $\frac{1}{15}$ & $1$ & $0$ & $0$ & $0$ & $0$ & $0$ & $0$ & $-1$ & $0$\\ row 2 & $\frac{3}{30}$ & $\frac{1}{5}$ & $0$ & $1$ & $0$ & $0$ & $0$ & $0$ & $0$ & $-1$ & $0$\\ row 3 & $2$ & $\frac{1}{2}$ & $0$ & $0$ & $1$ & $0$ & $0$ & $0$ & $0$ & $-1$ & $0$\\ row 4 & $1$ & $1$ & $0$ & $0$ & $0$ & $1$ & $0$ & $0$ & $0$ & $-1$ & $0$\\ row 5 & $2$ & $2$ & $0$ & $0$ & $0$ & $0$ & $1$ & $0$ & $0$ & $0$ & $10$\\ row 6 & $-2$ & $-2$ & $0$ & $0$ & $0$ & $0$ & $0$ & $1$ & $0$ & $0$ & $0$\\ row 7 & $-1.5$ & $1$ & $0$ & $0$ & $0$ & $0$ & $0$ & $0$ & $-1$ & $0$ & $4$\\\hline \hline $\tilde{J}\left ( x,z\right ) $ & $0$ & $0$ & $0$ & $0$ & $0$ & $0$ & $0$ & $0$ & $0$ & $1$ & $0$\\\hline \end{tabular} \ \]
Below is the Matlab code used to solve the above. The solution is\[ \fbox{$\begin{pmatrix} x_{1}\\ x_{2}\\ \tilde{x}_{3}\\ \tilde{x}_{4}\\ \tilde{x}_{5}\\ \tilde{x}_{6}\\ \tilde{x}_{7}\\ \tilde{x}_{8}\\ \tilde{x}_{9}\\ z \end{pmatrix} =\begin{pmatrix} 0\\ 4\\ 3.7333\\ 3.2\\ 2\\ 0\\ 2\\ 8\\ 0\\ 4 \end{pmatrix} $}\]
And minimum of the maximum of \(J_{i}\left ( x\right ) \) is\[ \fbox{$\min \max J_i\left ( x\right ) =4$}\]
The output from the above is
This was also solved using Mathematica. Here is the result, which confirms the above as well. This command is from Mathematica help, and used it to apply to this problem:
And the output is
Which is the same.
The variables are given in this table
| variable | description |
cost per unit | calories | potassium (mg) | sodium (mg) | calcium (mg) | phosphorus (mg) |
| \(u_{1}\) | whole milk |
\(0.4\) | \(660\) | \(210\) | \(75\) | \(1140\) | \(930\) |
| \(u_{2}\) | ice cream |
\(0.35\) | \(300\) | \(170\) | \(140\) | \(175\) | \(150\) |
| \(u_{3}\) | Eggs |
\(0.15\) | \(220\) | \(140\) | \(338\) | \(60\) | \(222\) |
| \(u_{4}\) | Cheese |
\(0.05\) | \(70\) | \(30\) | \(180\) | \(133\) | \(128\) |
| \(u_{5}\) | Beef | \(0.25\) | \(185\) | \(340\) | \(110\) | \(10\) | \(158\) |
| \(u_{6}\) | Broiled chicken | \(0.12\) | \(185\) | \(350\) | \(50\) | \(10\) | \(250\) |
| \(u_{7}\) | Baked vegetable | \(0.25\) | \(200\) | \(585\) | \(235\) | \(22\) | \(344\) |
| \(u_{8}\) | French fried |
\(0.07\) | \(155\) | \(510\) | \(6\) | \(9\) | \(6\) |
| \(u_{9}\) | Frozen Grain |
\(0.30\) | \(330\) | \(1315\) | \(4\) | \(69\) | \(115\) |
| \(u_{10}\) | Rice |
\(0.25\) | \(677\) | \(300\) | \(6\) | \(53\) | \(244\) |
Since the goal is to minimize the total cost of daily feeding, then
\[ J\left ( u\right ) =0.4u_{1}+0.35u_{2}+0.15u_{3}+0.05u_{4}+0.25u_{5}+0.12u_{6}+0.25u_{7}+0.07u_{8}+0.30u_{9}+0.25u_{10}\]
To be able to express the constraints, we need to first convert the minerals to same units used in \(J\left ( u\right ) \). For example, one \(u_{1}\) is one qt, which is \(976\) grams. The mineral calcium is \(1.14\) gram, therefore in units of \(u_{1}\) it becomes \(\frac{1.114}{976}u_{1}\), the same for all other minerals. Hence the above table becomes (where now each mineral is fraction of unit)
| variable | description | cost per unit | calories | potassium) | sodium (mg) | calcium (mg) | phosphorus (mg) |
| \(u_{1}\) | whole milk | \(0.4\) | \(660\) | \(\frac{0.210}{976}\) | \(\frac{0.075}{976}\) | \(\frac{1.140}{976}\) | \(\frac{0.930}{976}\) |
| \(u_{2}\) | ice cream | \(0.35\) | \(300\) | \(\frac{0.17}{188}\) | \(\frac{0.14}{188}\) | \(\frac{0.175}{188}\) | \(\frac{0.15}{188}\) |
| \(u_{3}\) | Eggs | \(0.15\) | \(220\) | \(\frac{0.14}{128}\) | \(\frac{0.338}{128}\) | \(\frac{0.06}{128}\) | \(\frac{0.222}{128}\) |
| \(u_{4}\) | Cheese | \(0.05\) | \(70\) | \(\frac{0.03}{17}\) | \(\frac{0.18}{17}\) | \(\frac{0.133}{17}\) | \(\frac{0.128}{17}\) |
| \(u_{5}\) | Beef | \(0.25\) | \(185\) | \(\frac{0.34}{85}\) | \(\frac{0.11}{85}\) | \(\frac{0.01}{85}\) | \(\frac{0.158}{85}\) |
| \(u_{6}\) | Broiled chicken | \(0.12\) | \(185\) | \(\frac{0.35}{85}\) | \(\frac{0.05}{85}\) | \(\frac{0.01}{85}\) | \(\frac{0.25}{85}\) |
| \(u_{7}\) | Baked vegetable | \(0.25\) | \(200\) | \(\frac{0.585}{100}\) | \(\frac{0.235}{100}\) | \(\frac{0.022}{100}\) | \(\frac{0.344}{100}\) |
| \(u_{8}\) | French fried | \(0.07\) | \(155\) | \(\frac{0.51}{60}\) | \(\frac{0.006}{60}\) | \(\frac{0.009}{60}\) | \(\frac{0.006}{60}\) |
| \(u_{9}\) | Frozen Grain | \(0.30\) | \(330\) | \(\frac{1.315}{210}\) | \(\frac{0.004}{210}\) | \(\frac{0.069}{210}\) | \(\frac{0.115}{210}\) |
| \(u_{10}\) | Rice | \(0.25\) | \(677\) | \(\frac{0.3}{187}\) | \(\frac{0.006}{187}\) | \(\frac{0.053}{187}\) | \(\frac{0.244}{187}\) |
Now we formulate each constraint. Let \(A\) be the total daily calories given by
\[ A=660u_{1}+300u_{2}+220u_{3}+70u_{4}+185u_{5}+185u_{6}+200u_{7}+155u_{8}+330u_{9}+677u_{10}\]
Hence constraint (i) in the problem becomes
\begin{align*} A & \geq 2400\\ A & \leq 2800 \end{align*}
For formulating constraint (ii), let \(B\) be the daily potassium and \(C\) be the daily sodium
\begin{align*} B & =\frac{0.210}{976}u_{1}+\frac{0.17}{188}u_{2}+\frac{0.14}{128}u_{3}+\frac{0.03}{17}u_{4}+\frac{0.34}{85}u_{5}+\frac{0.35}{85}u_{6}+\frac{0.585}{100}u_{7}+\frac{0.51}{60}u_{8}+\frac{1.315}{210}u_{9}+\frac{0.3}{187}u_{10}\\ C & =\frac{0.075}{976}u_{1}+\frac{0.14}{188}u_{2}+\frac{0.338}{128}u_{3}+\frac{0.18}{17}u_{4}+\frac{0.11}{85}u_{5}+\frac{0.05}{85}u_{6}+\frac{0.235}{100}u_{7}+\frac{0.006}{60}u_{8}+\frac{0.004}{210}u_{9}+\frac{0.006}{187}u_{10} \end{align*}
Hence constraint (ii) in the problem becomes
\begin{align*} B & \leq 1.15C\\ B & \geq 0.85C \end{align*}
For formulating constraint (iii), let \(D\) be the daily calcium and \(E\) be the daily phosphorus
\begin{align*} D & =\frac{1.140}{976}u_{1}+\frac{0.175}{188}u_{2}+\frac{0.06}{128}u_{3}+\frac{0.133}{17}u_{4}+\frac{0.01}{85}u_{5}+\frac{0.01}{85}u_{6}+\frac{0.022}{100}u_{7}+\frac{0.009}{60}u_{8}+\frac{0.069}{210}u_{9}+\frac{0.053}{187}u_{10}\\ E & =\frac{0.930}{976}u_{1}+\frac{0.15}{188}u_{2}+\frac{0.222}{128}u_{3}+\frac{0.128}{17}u_{4}+\frac{0.158}{85}u_{5}+\frac{0.25}{85}u_{6}+\frac{0.344}{100}u_{7}+\frac{0.006}{60}u_{8}+\frac{0.115}{210}u_{9}+\frac{0.244}{187}u_{10} \end{align*}
Hence constraint (iii) in the problem becomes
\[ D\geq 0.75E \]
Therefore the LP problem is
\(\begin{array} [c]{ccc}\min & & J\left ( u\right ) \\ s.t. & & A\geq 2400\\ & & A\leq 2800\\ & & B\leq 1.15C\\ & & B\geq 0.85C\\ & & D\geq 0.75E \end{array} \)
Converting to standard form and using \(x_{i}\) instead of \(u_{i}\) the above becomes
\(\begin{array} [c]{ccc}\min & & J\left ( x\right ) =0.4x_{1}+0.35x_{2}+0.15x_{3}+0.05x_{4}+0.25x_{5}+0.12x_{6}\\ & & +0.25x_{7}+0.07x_{8}+0.30x_{9}+0.25x_{10}\\ \text{subject to} & & A-\tilde{x}_{11}=2400\\ & & A+\tilde{x}_{12}=2800\\ & & B+\tilde{x}_{13}-1.15C=0\\ & & B-\tilde{x}_{14}-0.85C=0\\ & & D-\tilde{x}_{15}-0.75E=0 \end{array} \)
In full form, the above is
\( \begin{array} [c]{ccc}\min & & J\left ( x\right ) =0.4x_{1}+0.35x_{2}+0.15x_{3}+0.05x_{4}+0.25x_{5}+0.12x_{6}+0.25x_{7}+0.07x_{8}+0.30x_{9}+0.25x_{10}\\ s.t. & & 660x_{1}+300x_{2}+220x_{3}+70x_{4}+185x_{5}+185x_{6}+200x_{7}+155x_{8}+330x_{9}+677x_{10}-\tilde{x}_{11}=2400\\ & & 660x_{1}+300x_{2}+220x_{3}+70x_{4}+185x_{5}+185x_{6}+200x_{7}+155x_{8}+330x_{9}+677x_{10}+\tilde{x}_{12}=2800\\ & & \frac{0.210}{976}x_{1}+\frac{0.17}{188}x_{2}+\frac{0.14}{128}x_{3}+\frac{0.03}{17}x_{4}+\frac{0.34}{85}x_{5}+\frac{0.35}{85}x_{6}+\frac{0.585}{100}x_{7}+\frac{0.51}{60}x_{8}+\frac{1.315}{210}x_{9}+\frac{0.3}{187}x_{10}+\tilde{x}_{13}-1.15\left ( \frac{0.075}{976}x_{1}+\frac{0.14}{188}x_{2}+\frac{0.338}{128}x_{3}+\frac{0.18}{17}x_{4}+\frac{0.11}{85}x_{5}+\frac{0.05}{85}x_{6}+\frac{0.235}{100}x_{7}+\frac{0.006}{60}x_{8}+\frac{0.004}{210}x_{9}+\frac{0.006}{187}x_{10}\right ) =0\\ & & \frac{0.210}{976}x_{1}+\frac{0.17}{188}x_{2}+\frac{0.14}{128}x_{3}+\frac{0.03}{17}x_{4}+\frac{0.34}{85}x_{5}+\frac{0.35}{85}x_{6}+\frac{0.585}{100}x_{7}+\frac{0.51}{60}x_{8}+\frac{1.315}{210}x_{9}+\frac{0.3}{187}x_{10}-\tilde{x}_{14}-0.85\left ( \frac{0.075}{976}x_{1}+\frac{0.14}{188}x_{2}+\frac{0.338}{128}x_{3}+\frac{0.18}{17}x_{4}+\frac{0.11}{85}x_{5}+\frac{0.05}{85}x_{6}+\frac{0.235}{100}x_{7}+\frac{0.006}{60}x_{8}+\frac{0.004}{210}x_{9}+\frac{0.006}{187}x_{10}\right ) =0\\ & & \frac{1.140}{976}x_{1}+\frac{0.175}{188}x_{2}+\frac{0.06}{128}x_{3}+\frac{0.133}{17}x_{4}+\frac{0.01}{85}x_{5}+\frac{0.01}{85}x_{6}+\frac{0.022}{100}x_{7}+\frac{0.009}{60}x_{8}+\frac{0.069}{210}x_{9}+\frac{0.053}{187}x_{10}-\tilde{x}_{15}-0.75\left ( \frac{0.930}{976}x_{1}+\frac{0.15}{188}x_{2}+\frac{0.222}{128}x_{3}+\frac{0.128}{17}x_{4}+\frac{0.158}{85}x_{5}+\frac{0.25}{85}x_{6}+\frac{0.344}{100}x_{7}+\frac{0.006}{60}x_{8}+\frac{0.115}{210}x_{9}+\frac{0.244}{187}x_{10}\right ) =0 \end{array} \)
Simplifying gives
\( \begin{array} [c]{ccc}\min & & J\left ( x\right ) =0.4x_{1}+0.35x_{2}+0.15x_{3}+0.05x_{4}+0.25x_{5}+0.12x_{6}+0.25x_{7}+0.07x_{8}+0.3x_{9}+0.25x_{10}\\ s.t. & & 660x_{1}+300x_{2}+220x_{3}+70x_{4}+185x_{5}+185x_{6}+200x_{7}+155x_{8}+330x_{9}+677x_{10}-\tilde{x}_{11}=2400\\ & & 660x_{1}+300x_{2}+220x_{3}+70x_{4}+185x_{5}+185x_{6}+200x_{7}+155x_{8}+330x_{9}+677x_{10}+\tilde{x}_{12}=2800\\ & & 1.268\times 10^{-4}x_{1}+4.787\times 10^{-5}x_{2}-1.943\times 10^{-3}x_{3}-1.041\times 10^{-2}x_{4}+2.512\times 10^{-3}x_{5}+3.441\times 10^{-3}x_{6}+3.148\times 10^{-3}x_{7}+8.385\times 10^{-3}x_{8}+0.006x_{9}+1.567\times 10^{-3}x_{10}+\tilde{x}_{13}=0\\ & & 1.499\times 10^{-4}x_{1}+2.713\times 10^{-4}x_{2}-1.151\times 10^{-3}x_{3}-7.235\times 10^{-3}x_{4}+0.003x_{5}+3.618\times 10^{-3}x_{6}+3.853\times 10^{-3}x_{7}+8.415\times 10^{-3}x_{8}+6.246\times 10^{-3}x_{9}+1.577\times 10^{-3}x_{10}-\tilde{x}_{14}=0\\ & & 4.534\times 10^{-4}x_{1}+3.325\times 10^{-4}x_{2}-8.32\times 10^{-4}x_{3}+2.176\times 10^{-3}x_{4}-1.277\times 10^{-3}x_{5}-2.088\times 10^{-3}x_{6}-0.002x_{7}+7.5\times 10^{-5}x_{8}-8.214\times 10^{-5}x_{9}-6.952\times 10^{-4}x_{10}-\tilde{x}_{15}=0 \end{array} \)
Therefore \(c^{T}x\) becomes\[ \min c^{T}x=\begin{pmatrix} 0.4 & 0.35 & 0.15 & 0.05 & 0.25 & 0.12 & 0.25 & 0.07 & 0.3 & 0.25 & 0 & 0 & 0 & 0 & 0 \end{pmatrix}\begin{pmatrix} x_{1}\\ x_{2}\\ x_{3}\\ x_{4}\\ x_{5}\\ x_{6}\\ x_{7}\\ x_{8}\\ x_{9}\\ x_{10}\\ \tilde{x}_{11}\\ \tilde{x}_{12}\\ \tilde{x}_{13}\\ \tilde{x}_{14}\\ \tilde{x}_{15}\end{pmatrix} \] subject to \(Ax=b\)
The above is now solved using Matlab linprog. The following is the result, followed by the Matlab code used.\[ x_{optimal}=\begin{pmatrix} x_{1}\\ x_{2}\\ x_{3}\\ x_{4}\\ x_{5}\\ x_{6}\\ x_{7}\\ x_{8}\\ x_{9}\\ x_{10}\\ \tilde{x}_{11}\\ \tilde{x}_{12}\\ \tilde{x}_{13}\\ \tilde{x}_{14}\\ \tilde{x}_{15}\end{pmatrix} =\begin{pmatrix} 0\\ 0\\ 0\\ 1.0792\\ 0\\ 0\\ 0\\ 0.2889\\ 0\\ 3.41\\ 0\\ 400\\ 0.0035\\ 0\\ 0 \end{pmatrix} \] In terms of raw \(u_{i}\) variables, the above says to buy 2
Verification: Calories from the above is \[ \left ( 1.0792\right ) \left ( 70\right ) +\left ( 0.2889\right ) \left ( 155\right ) +\left ( 3.41\right ) \left ( 677\right ) =2428.9 \] Which is within daily allowance OK. Potassium from the above is \[ \left ( 1.0792\right ) \left ( \frac{0.03}{17}\right ) +\left ( 0.2889\right ) \left ( \frac{0.51}{60}\right ) +\left ( 3.41\right ) \left ( \frac{0.3}{187}\right ) =9.831\text{ mg}\] While sodium from the above is \[ \left ( 1.0792\right ) \left ( \frac{0.18}{17}\right ) +\left ( 0.2889\right ) \left ( \frac{0.006}{60}\right ) +\left ( 3.41\right ) \left ( \frac{0.006}{187}\right ) =10.1565\text{ mg}\] Hence Potassium is within \(15\%\) of sodium OK. Daily calcium from above is\[ \left ( 1.0792\right ) \left ( \frac{0.133}{17}\right ) +\left ( 0.2889\right ) \left ( \frac{0.009}{60}\right ) +\left ( 3.41\right ) \left ( \frac{0.053}{187}\right ) =9.453\text{ mg}\] While daily phosphorus is\[ \left ( 1.0792\right ) \left ( \frac{0.128}{17}\right ) +\left ( 0.2889\right ) \left ( \frac{0.006}{60}\right ) +\left ( 3.41\right ) \left ( \frac{0.244}{187}\right ) =12.6\text{ mg}\] Hence Daily calcium is at least \(75\%\) of daily phosphorus OK. All verified OK. The corresponding optimal daily cost is\[ \fbox{$J^\ast =0.9267$ dollars}\]
The output from the above is
Now we add a new constraint, which is \(B_{1}+B_{2}=\gamma (fat)\). In terms of \(u_{i}\), this becomes\[ B_{1}=\frac{0.00032}{976}u_{1}+\frac{0.0006}{210}u_{9}+\frac{0.0003}{187}u_{10}\] And\[ B_{2}=\frac{0.0017}{976}u_{1}+\frac{0.0003}{188}u_{2}+\frac{0.0004}{128}u_{3}+\frac{0.0001}{17}u_{4}+\frac{0.0001}{85}u_{6}+\frac{0.0001}{210}u_{9}\] And\[ fat=\frac{40}{976}u_{1}+\frac{18}{188}u_{2}+\frac{16}{128}u_{3}+\frac{6}{17}u_{4}+\frac{10}{85}u_{5}+\frac{9}{85}u_{6}+\frac{8}{100}u_{7}+\frac{7}{60}u_{8}\] Hence the new constraint is\[ \left ( B_{1}+B_{2}\right ) -\gamma \left ( fat\right ) =0 \] Converting to standard LP form, the above becomes
Now the new \(Ax=b\) becomes (\(c^{T}x\) do not change from part(a) and remain the same).
The above was solve in Matlab for different values of \(\gamma .\)For \(\gamma =0.00001\), the optimal \(x\) was\[ x_{optimal}=\begin{pmatrix} x_{1}\\ x_{2}\\ x_{3}\\ x_{4}\\ x_{5}\\ x_{6}\\ x_{7}\\ x_{8}\\ x_{9}\\ x_{10}\\ \tilde{x}_{11}\\ \tilde{x}_{12}\\ \tilde{x}_{13}\\ \tilde{x}_{14}\\ \tilde{x}_{15}\end{pmatrix} =\begin{pmatrix} 0\\ 0\\ 0\\ 4.2150\\ 6.8174\\ 0\\ 0\\ 3.0043\\ 0\\ 1.0022\\ 0\\ 400\\ 0\\ 0.0161\\ 0 \end{pmatrix} \] With optimal \(J^{\ast }=2.376\) dollars. For \(\gamma =0.00002\), the optimal \(x\) was\[ x_{optimal}=\begin{pmatrix} x_{1}\\ x_{2}\\ x_{3}\\ x_{4}\\ x_{5}\\ x_{6}\\ x_{7}\\ x_{8}\\ x_{9}\\ x_{10}\\ \tilde{x}_{11}\\ \tilde{x}_{12}\\ \tilde{x}_{13}\\ \tilde{x}_{14}\\ \tilde{x}_{15}\end{pmatrix} =\begin{pmatrix} 0\\ 0\\ 0\\ 2.0750\\ 0\\ 0\\ 0\\ 1.1779\\ 0\\ 3.2348\\ 0\\ 400\\ 0.0067\\ 0\\ 0.0024 \end{pmatrix} \] With optimal \(J^{\ast }=0.9949\) dollars. For \(\gamma =0.00003\), the optimal \(x\) was\[ x_{optimal}=\begin{pmatrix} x_{1}\\ x_{2}\\ x_{3}\\ x_{4}\\ x_{5}\\ x_{6}\\ x_{7}\\ x_{8}\\ x_{9}\\ x_{10}\\ \tilde{x}_{11}\\ \tilde{x}_{12}\\ \tilde{x}_{13}\\ \tilde{x}_{14}\\ \tilde{x}_{15}\end{pmatrix} =\begin{pmatrix} 0\\ 0\\ 0\\ 1.0712\\ 0\\ 0\\ 0\\ 0.2078\\ 0.1119\\ 3.3629\\ 0\\ 400\\ 0.0034\\ 0\\ 0 \end{pmatrix} \] With optimal \(J^{\ast }=0.9424\) dollars. For \(\gamma =0.00004\), the optimal \(x\) was\[ x_{optimal}=\begin{pmatrix} x_{1}\\ x_{2}\\ x_{3}\\ x_{4}\\ x_{5}\\ x_{6}\\ x_{7}\\ x_{8}\\ x_{9}\\ x_{10}\\ \tilde{x}_{11}\\ \tilde{x}_{12}\\ \tilde{x}_{13}\\ \tilde{x}_{14}\\ \tilde{x}_{15}\end{pmatrix} =\begin{pmatrix} 0.3431\\ 0\\ 0\\ 0.8519\\ 0\\ 0\\ 0\\ 0.7094\\ 2.8071\\ 0\\ 0\\ 400\\ 0\\ 0.0027\\ 0 \end{pmatrix} \] With optimal \(J^{\ast }=1.0944\) dollars\(.\)When going more than \(\gamma =0.00004\), say \(\gamma =0.00005\), Matlab was not able to obtain solution using simplex method. The message was ”Exiting: The constraints are overly stringent; no feasible starting point found”
When going smaller than \(\gamma =0.00001\) same problem was seen. Out of these \(\gamma \) values, \(\gamma =0.00003\) gave the least cost of \(0.9424\) dollars.
