%% Problem Set G, Problem 10, Part (b)
% This is a problem illustrating the Lagrange multiplier method. We are trying to maximize and minimize a function f(x,y) on the ellipse g(x,y)=x^2+3y^2=1.
%% (i)
% First, we reduce to a one-variable problem and use fminbnd.

syms x y t real
f = exp(x) - sin(y);
g = x^2 + 3*y^2;
onevarfunct = subs(f, [x, y], [cos(t), sin(t)/sqrt(3)])

%
% Let's plot this one-variable function to see where the maximum and minimum seem to be located.

fplot(onevarfunct, [0, 2*pi])

%
% It seems from this that the maximum occurs near t=6 and the minimum occurs near t=2.3. Indeed, we see:

F = matlabFunction(onevarfunct);
xmin = fminbnd(F, 0, 2*pi)
F(xmin)
xmax = fminbnd(matlabFunction(-onevarfunct), 0, 2*pi)
F(xmax)

%
%% (ii)
% Now we use the Lagrange multiplier method. We will plot the constraint curve in a heavy line style and the Lagrange multiplier equation with a dashed curve. After eliminating \lambda, we want to set to zero the expression

lagr = diff(f,x)*diff(g,y) - diff(f,y)*diff(g,x)
fimplicit(g - 1, [-1, 1, -1, 1], 'LineWidth', 2), hold on
fimplicit(lagr, [-1, 1, -1, 1], '--'), hold off, axis equal

%
% We can see two intersection points, one near (-.6, .4) and one near (.95, -.1). Indeed, using the script newton2d from Problem 7.5, we can get fairly accurate values for x and y, and compute the values of the objective function f at the two points.

pt1 = newton2d(g-1, lagr, x, y, -.6, .4); pt1 
minval = double(subs(f, [x,y], pt1))
pt2 = newton2d(g-1, lagr, x, y, .95, -.1); pt2 
maxval = double(subs(f, [x,y], pt2))

%
% Clearly pt2 gives the maximum value, around 2.7793, and pt1 gives the minimum value, around 0.0957. These match what we found by the other method. Finally, we show that the maximum and minimum occur where the level curves of f are tangent to the ellipse. The dark dots indicate the maximum and minimum points, with the maximum in the lower right and the minimum in the upper left.

list = [minval, maxval, -1:0.5:4];
fimplicit(g - 1, [-1, 1, -1, 1], 'LineWidth', 2), hold on
plot(pt1(1), pt1(2), '.k', 'MarkerSize', 20)
plot(pt2(1), pt2(2), '.k', 'MarkerSize', 20)
fcontour(f, [-1, 1, -1, 1], 'LevelList', list), hold off, axis equal