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.