%% Linear Programming
% MATLAB is ideally suited to handle so-called linear programming
% problems.  These are problems in which you have a quantity,
% depending linearly on several variables, that you want to maximize
% or minimize subject to several constraints that are expressed as
% linear inequalities in the same variables.  If the number of
% variables and the number of constraints are small, then there are
% numerous mathematical techniques for solving a linear programming
% problem -- indeed, these techniques are often taught in high-school
% or university-level courses in finite mathematics.  But sometimes
% these numbers are high, or, even if they are low, the constants in the linear
% inequalities or the object expression for the quantity to be
% optimized may be numerically complicated -- in which case a
% software package like MATLAB is required to effect a solution.  We
% shall illustrate the method of linear programming by means of a
% simple example, giving a combined graphical/numerical solution,
% and then solve both a slightly as well as a substantially more
% complicated problem.
%
% Suppose that a farmer has $75$ acres on which to plant two crops: wheat
% and barley.  To produce these crops, it costs the farmer (for seed,
% fertilizer, etc.) \$120 per acre for the wheat and \$210 per
% acre for the barley.  The farmer has \$15,000 available for
% expenses.  But after the harvest, the farmer must store the crops
% while awaiting favorable market conditions.  The farmer has storage
% space for $4000$ bushels.  Each acre yields an average of $110$
% bushels of wheat or $30$ bushels of barley.  If the net profit per
% bushel of wheat (after all expenses have been subtracted) is
% \$1.30 and for barley is \$2.00, how should the farmer plant the
% 75 acres to maximize profit?
%
% We begin by formulating the problem mathematically.  First, we
% express the objective, that is the profit, and the constraints
% algebraically, then we graph them, and lastly we arrive at the
% solution by graphical inspection and a minor arithmetic
% calculation.
%
% Let $x$ denote the number of acres allotted to wheat and $y$ the
% number of acres allotted to barley.  Then the expression to be
% maximized, that is the profit, is clearly
%
% $$P = (110)(1.30)x + (30)(2.00)y = 143x + 60y.$$
%
% There are three constraint inequalities, specified by the limits on
% expenses, storage, and acreage.  They are, respectively,
%
% $$120x + 210y \leq 15000,$$
% $$110x + 30y \leq 4000,$$
% $$x + y \leq 75.$$
%
% Strictly speaking there are two more constraint inequalities forced
% by the fact that the farmer cannot plant a negative number of acres,
% namely
%
% $$x \geq 0, \quad y \geq 0.$$
%
% Next we graph the regions specified by the constraints.  The last
% two say that we need to consider only the first quadrant in the
% $x$-$y$ plane.  Here's a graph delineating the triangular region in
% the first quadrant determined by the first inequality:

X = 0:125;
Y1 = (15000 - 120.*X)./210; 
area(X, Y1)  

%%
% Now let's put in the other two constraint inequalities:

Y2 = max((4000 - 110.*X)./30, 0); 
Y3 = max(75 - X, 0); 
Ytop = min([Y1; Y2; Y3]); 
area(X, Ytop)  

%%
% It's a little hard to see the polygonal boundary of the region
% clearly.  Let's home in a bit:

area(X, Ytop); axis([0 40 40 75])  

%%
% Now let's superimpose on top of this picture a contour plot of the
% objective function $P$:

hold on
[U V] = meshgrid(0:40, 40:75);
contour(U, V, 143.*U + 60.*V); hold off  

%%
% It seems apparent that the maximum value of $P$ will occur on the
% level curve (that is, level line) that passes through the vertex of
% the polygon that lies near $(22,53)$.  In fact, we can compute:

[x, y] = solve('x + y = 75', '110*x + 30*y = 4000')  

%%

double([x, y])  

%%
% The acreage that results in the maximum profit is $21.875$ for wheat
% and $53.125$ for barley.  In that case, the profit is:

P = 143*x + 60*y  

%%

format bank; double(P)  

%%
% that is, \$6,315.63. 

%%
% This problem illustrates and is governed by the ``Fundamental Theorem
% of Linear Programming,'' which is stated here for two variables: a linear
% expression $ax + by$, defined over a closed bounded convex set $S$
% whose sides are line segments, takes on its maximum
% and minimum values at vertices of $S$.  If $S$ is
% unbounded, there might but need not be an optimum value, but, if there
% is, it occurs at a vertex.  (A convex set is one for which any line
% segment joining two points of the set lies entirely inside the set.)
%
% In fact, the Optimization Toolbox has a built-in function, *|linprog|*, 
% which implements the solution of a linear programming problem.  You can 
% learn about it from the online help. We will use *|linprog|* on the 
% above problem.  After that, we will use it to solve two more complicated 
% problems involving more variables and constraints.  Here is the 
% beginning of its help text:

helptext = help('linprog'); helptext(1:180)

%%
% So:

f = [-143 -60];  
A = [120 210; 110 30; 1 1; -1 0; 0 -1]; 
b = [15000; 4000; 75; 0; 0];  

format short; linprog(f, A, b)

%%
% This is the same answer that we obtained before.  Note that we entered
% the negative of the coefficients for the objective function $P$
% into the vector *|f|* because *|linprog|* searches for a minimum rather 
% than a maximum.  Note also that the non-negativity constraints are
% accounted for in the last two rows of *|A|* and *|b|*.

%%
% Well, we could have done this problem by hand.  But suppose that the
% farmer is dealing with a third crop, say corn, and that the
% corresponding data are:
%
% * Cost per acre: \$150.75,
% * Yield per acre: $125$ bushels,
% * Profit per bushel: \$1.56.
%
% If we denote the number of acres allotted to corn by $z$, then the
% objective function becomes
%
% $$P = (110)(1.30)x + (30)(2.00)y  + (125)(1.56)= 143x + 60y + 195z,$$
%
% and the constraint inequalities are
%
% $$120x + 210y + 150.75z \leq 15000,$$
% $$110x + 30y + 125z \leq 4000,$$
% $$x + y + z \leq 75,$$
% $$x \geq 0, \quad y \geq 0, \quad z \geq 0.$$
%
% The problem is solved with *|linprog|* as follows:

f = [-143 -60 -195];
A = [120 210 150.75; 110 30 125; 1 1 1; -1 0 0; ...
    0 -1 0; 0 0 -1]; 
b = [15000; 4000; 75; 0; 0; 0];  
linprog(f, A, b) 

%%
% So the farmer should ditch the wheat and plant $56.5789$ acres of
% barley and $18.4211$ acres of corn.

%%
% There is no practical limit on the number of variables and
% constraints that MATLAB can handle -- certainly none that the
% relatively unsophisticated user will encounter.  Indeed, in many
% real applications of the technique of linear programming, you may need
% to deal with many variables and constraints.  The solution of such a
% problem by hand is not feasible, and software like MATLAB is crucial
% to success.  For example, in the farming problem with which we have
% been working, you could have more crops than two or three -- think
% agribusiness instead of family farmer.  Also, you could have
% constraints that arise from other things beside expenses, storage,
% and acreage limitations -- for example the following.
%
% * Availability of seed.  This might lead to constraint inequalities
% such as $x_j \leq k$. 
% * Personal preferences.  The farmer's spouse might have a
% preference for one variety over another and insist on a
% corresponding planting, or something similar with a collection of
% crops; thus, constraint inequalities such as $x_i \leq x_j$ or $x_1 +
% x_2 \geq x_3$.
% * Government subsidies.  It may take a moment's reflection on the
% reader's part, but this could lead to inequalities such as $x_j \geq
% k$.
%
% Below is a sequence of commands that solves exactly such a problem.
% You should be able to recognize the objective expression and the
% constraints from the data that are entered.  But as an aid, you might
% answer the following questions.
%
% * How many crops are under consideration?
% * What are the corresponding expenses? How much is available for expenses?
% * What are the yields in each case? What is the storage capacity?
% * How many acres are available?
% * What crops are constrained by seed limitations? To what extent?
% * What about preferences?
% * What are the minimum acreages for each crop?

f = [-110*1.3 -30*2.0 -125*1.56 -75*1.8 -95*.95 ...
    -100*2.25 -50*1.35];  
A = [120 210 150.75 115 186 140 85; ...
    110 30 125 75 95 100 50; 1 1 1 1 1 1 1; 1 0 0 0 0 0 0; ...
    0 0 1 0 0 0 0; 0 0 0 0 0 1 0; 1 -1 0 0 0 0 0; ...
    0 0 1 0 -2 0 0; 0 0 0 -1 0 -1 1; -1 0 0 0 0 0 0; ...
    0 -1 0 0 0 0 0; 0 0 -1 0 0 0 0; 0 0 0 -1 0 0 0; ...
    0 0 0 0 -1 0 0; 0 0 0 0 0 -1 0; 0 0 0 0 0 0 -1];
b = [55000; 40000; 400; 100; 50; 250; 0; 0; 0; -10; -10; ...
    -10; -10; -20; -20; -20]; 
linprog(f, A, b)  

%%
% Note that, despite the complexity of the problem, MATLAB solves it
% almost instantaneously.  It seems the farmer should bet the farm on
% crop number 6.  We suggest that you alter the expense and/or
% the storage limit in the problem and see the effect on the
% answer.