We are fortunate to have the optimization toolbox with MATLAB 5.3 in the computer labs.
In the last lecture, we demonstrated the simplex algorithm on the following simple problem
max 3*x1 + 5*x2 min 4*y1 + 12*y2 + 18*y3
subject to subject to
x1 <= 4 y1 + 3*y3 >= 3
2*x2 <= 12 2*y2 + 2*y3 >= 5
3*x1 + 2*x2 <= 18 y1 >= 0, y2 >= 0
x1 >= 0, x2 >= 0
At the right we see the dual version of this problem.
We can interpret the original problem on the left as
the determination of the optimal mix for a production
of two goods, x1 (sells at $3/unit) and x2 (sells at $5/unit),
subject to production capacity constraints.
The duality theory is needed to solve this problem in MATLAB, as the linprog command only treats the minimization.
% L-16 MCS 494 Wed 2 Oct 2002 Using the command linprog in MATLAB
diary h:\lec16.txt % start logging to file
help linprog % read documentation about linprog
f = [4 12 18]' % coefficients of objective
A = [-1 0 -3 % coefficients of constraints
0 -2 -2
-1 0 0 % -y(1) <= 0
0 -1 0
0 0 -1]
b = [-3 -5 0 0 0]' % right-hand side vector
[y,fval] = linprog(f,A,b) % solves the dual problem
% now we wish to see the solution to the primal problem
B = A(1:2,2:3) % basis in y is y(2) and y(3)
-B'\f(2:3) % solve for x
The last two steps may require some explanation...
From the solution to the dual problem, we see that the
nonzero variables are y(2) and y(3). With the basis we
select those linear equations that determine the vertex
at the optimum. So we have the linear system
B*[y(2) y(3)]' = [b(1) b(2)]'To determine the vertex in the original x coordinates, we transpose the matrix B, flip the signs and solve
-B'*[x(1) x(2)]' = [f(2) f(3)]'knowing that the right-hand vector of the dual problem contains the coefficients of the objective function of the primal problem, and vice versa. (By the way, the backslash is the shortcut for solving a linear system in MATLAB.)
We continue our session, introducing Chapter 6 on regression, with the command polyfit of MATLAB: The experiment we conduct is an illustration on how polynomial fitting can be used to remove noise from data.
% L-16 MCS 494 Wed 2 Oct 2002 Fitting with Polynomials in MATLAB help polyval % see how MATLAB treats polynomials c = [3 -2.23 -5.1 9.8]; % 3*x^3 - 2.23*x^2 - 5.1*x + 9.8 x = -1:0.1:1; % generate range to plot the polynomial y = polyval(c,x); % evaluate the polynomial plot(x,y) % we will generate some noise on the vector y noise = randn(1,size(y,2)); % random vector of same size as y noise = noise/norm(noise); % make vector of length one ey = y + noise; % perturbed data % now we have a data set (x,ey) and we wish to find a polynomial figure(1) % make first figure current plot(x,ey,'r+'); hold on plot(x,y) % we wish to reconstruct the polynomial ec = polyfit(x,ey,3); % cubic fit yy = polyval(ec,x); % sample fitted polynomial plot(x,yy,'r--'); diary off % stop logging to file