Computer assignment in MATLAB
Critical Points of Functions of Two and Three Variables
Assignment Participants
Andreas Karlsson
Daniel Lehtinen
Contents
Problem 1 .............................................................................. 1
Solution procedure .......................................................................
MATLAB calculations .................................................................. 2
MATLAB code ......................................................................
Calculation results ............................................................. 3
Plots ........................................................................................ 4
Conclusion ............................................................................... 5
Problem 2 .............................................................................. 6
Solution procedure .......................................................................
MATLAB calculations .....................................................................
MATLAB code ......................................................................
Calculation results ............................................................. 7
Conclusion ..................................................................................
Problem 3 .............................................................................. 8
Solution procedure .......................................................................
MATLAB calculations .................................................................. 9
MATLAB code .....................................................................
Calculation results ..............................................................
Plots ..................................................................................... 10
Conclusion ............................................................................. 11
Problem 4 ............................................................................ 12
Solution procedure .......................................................................
MATLAB calculations .....................................................................
MATLAB code ......................................................................
Calculation results ........................................................... 13
Conclusion ..................................................................................
Problem 5 ............................................................................ 14
Solution procedure .......................................................................
MATLAB calculations .....................................................................
MATLAB code ......................................................................
Calculation results ........................................................... 15
Plots ...................................................................................... 16
Conclusion ............................................................................. 17
Problem 6 ............................................................................ 18
Solution procedure .......................................................................
MATLAB calculations .....................................................................
MATLAB code ......................................................................
Calculation results ........................................................... 19
Plots ...................................................................................... 20
Conclusion ..................................................................................
Problem 1
Find and classify all critical points of the function
g x , y=
x
4
2x
3
y−6x
2
y
2
y
4
x
4
y
4
1
,
which you studied in the previous lesson, and compare your results with the contour and
gradient plot of g that you have already made.
Solution procedure
To find critical points of g, we must set the partial derivatives equal to 0 and solve for x and y.
With the coordinates in hand, our next step is to calculate the second derivatives. These
second derivatives will make the Hessian matrix. Finally, we calculate the Hessian determinant
and this value will classify the critical points for us.
MATLAB calculations
MATLAB code
syms x y z
f=(x^4 + 2*(x^3)*y - 6*(x^2)*(y^2) + y^4)/(x^4 + y^4 +1)
%figure(1);
%ezsurf(f, [-3, 3, -3, 3])
fx=diff(f,x);
sfx = simplify(fx);
fy=diff(f,y);
sfy = simplify(fy);
[xcr,ycr]=solve(fx,fy);
[xcr,ycr] %punkterna i vilka derivatorna är noll
fxx=diff(fx,x);
sfxx=simplify(fxx);
fxy=diff(fx,y);
sfxy=simplify(fxy);
fyy=diff(fy,y);
sfyy=simplify(fyy);
%hessdetf=fxx*fyy-fxy^2 %hessiandeterminanten
hessdetf_simplified = simplify(hessdetf)
hessdetf_simplified_inserted = subs(hessdetf_simplified, [x,y], [0,0])
gradf = jacobian(f, [x, y])
hessmatf = jacobian(gradf, [x, y])
hessmatf_inserted = subs(hessmatf, [x,y], [0,0])
hessdetf1 = simplify(det(hessmatf))
%simplify(hessdetf - hessdetf1)
[xcr(1), ycr(1), subs(hessdetf, [x,y], [xcr(1), ycr(1)]), subs(fxx, [x,y],
[xcr(1), ycr(1)])]
%syms h k
%q = fxx*h^2 + 2*fxy*h*k + fyy*k^2
%q1 = simplify(q)
%q2 = subs(q, [x,y], [0,0])
%subs(fx, [x,y], [0,0])
%subs(fy, [x,y], [0,0])
[xx, yy] = meshgrid(-3:.1:3,-3:.1:3);
ffun = inline(vectorize(f));
fxfun = inline(vectorize(gradf(1)));
fyfun = inline(vectorize(gradf(2)));
%figure(2);
%contour(xx, yy, ffun(xx,yy), 30)
%hold on
%[xx, yy] = meshgrid(-3:.25:3,-3:.25:3);
%quiver(xx, yy, fxfun(xx,yy), fyfun(xx,yy), 0.6)
%axis equal tight, hold off
Calculation results
Partial derivatives of g
Simplified using MATLABs built-in function simplify().
f '
x
=−2x
−2x
2
x
5
y−3xy
5
−3xy−6y
2
x
4
6y
6
6y
2
x
4
y
4
1
2
f '
y
=2
x
7
−3x
3
y
4
x
3
−6x
6
y6x
2
y
5
−6x
2
y2y
3
x
4
y
4
1
2
Which resulted in one point (0, 0).
Second derivatives of g
Due to the length of the derivatives, we did not reproduce them here.
f '
xx
=...
f '
xy
= f '
yx
=...
f '
yy
=...
Hessian matrix (the square matrix of second partial derivatives)
H x , y =
f '
xx
f '
xy
f '
yx
f '
yy
H 0,0=
0 0
0 0
det H 0,0 =0
Plots
Figure 1. Contour plot
Figure 2. Surf plot
Conclusion
det(H) > 0 and any second derivate is positive => local minimum.
det(H) > 0 and any second derivate is negative => local maximum.
det(H) < 0 => saddlepoint.
Since the Hessian determinant equals 0, the second derivatives test to evaluate if the critical
points is either max, min or saddle points fails. The figure in the plot looks like a pair of
crossed saddle points with some kind of plane in the intersection.
Problem 2
Find and classify the critical points of the function
h x , y , z =3x
2
4y
2
z
2
−9xyz.
You may have trouble substituting the symbolic values of the coordinates of the critical points
into the Hessian matrix; if so, convert them to numerical values, using double.
Solution procedure
We calculate the first partial derivatives and set them to 0 to receive the critical points.
We evaluate the Hessian matrix at the critical points and compute the eigenvalues of the
matrix. The values of the eigenvalues will determine the type of critical points.
MATLAB calculations
MATLAB code
syms x y z
f=(3*x^2 + 4*y^2 + z^2 - 9*x*y*z)
gradf2 = jacobian(f,[x,y,z])
hessmat = jacobian(gradf2, [x,y,z])
[xcr2, ycr2, zcr2] = solve(gradf2(1),gradf2(2),gradf2(3));
[xcr2, ycr2, zcr2]
H1=subs(hessmat,[x,y,z],[xcr2(1), ycr2(1), zcr2(1)])
H2=subs(hessmat,[x,y,z],[xcr2(2), ycr2(2), zcr2(2)])
H3=subs(hessmat,[x,y,z],[xcr2(3), ycr2(3), zcr2(3)])
H4=subs(hessmat,[x,y,z],[xcr2(4), ycr2(4), zcr2(4)])
H5=subs(hessmat,[x,y,z],[xcr2(5), ycr2(5), zcr2(5)])
double(eig(H1))
double(eig(H2))
double(eig(H3))
double(eig(H4))
double(eig(H5))
Calculation results
Partial derivatives
f '
x
=6x−9yz
f '
y
=8y−9xz
f '
z
=2z−9xy
Hessian matrix
H x , y , z=
6 −9z −9y
−9z 8 −9x
−9y −9x 2
Cordinates for critical points
0 0 0
4/9 2/9⋅
3 4/9⋅
3
−4 /9 2/9⋅
3 −4/9⋅
3
−4 /9 −2 /9⋅
3 4/9⋅
3
4/9 −2 /9⋅
3 −4/9⋅
3
Eigenvalues
H1 = (6, 8, 2)
H2 = (14.0565, -4.3445, 6.2880)
H3 = (14.0565, -4.3445, 6.2880)
H4 = (14.0565, -4.3445, 6.2880)
H5 = (14.0565, -4.3445, 6.2880)
Conclusion
The significance of the eigenvalues of the Hessian matrix is that if all of them are positive at a
critical point, the function has a local minimum there; if all are negative, the function has a
local maximum; if they have mixed signs, the function has a saddle point; and if at least one
of them is 0, the critical point is degenerate.
According to this, our function results in a single local minimum and four saddle points.
Problem 3
Find and classify all critical points of the function
x
4
y
4
xy −x
2
− y
2
.
Relate your results to a simultaneous contour and gradient plot of the function.
Solution procedure
First we plot the f'x=0 and f'y=0 to find the points where f'x=f'y=0. We estimated the
coordinates for when the curves intersected in the plot. We sent those estimated values into
the newton2d function to get the exact values.
Then we calculate the second derivatives and use them to calculate the Hessian determinant,
in which we use our coordinates.
MATLAB calculations
MATLAB code
syms x y
k=(x^4 + y^4 + x*y -x^2 -y^2)
kx=diff(k,x)
ky=diff(k,y)
[x1,y1]=meshgrid(-8:.1:8,-8:.1:8);
kxfun=inline(vectorize(kx));
kyfun=inline(vectorize(ky));
figure(1);
contour(x1,y1,kxfun(x1,y1),[0,0],'r'), hold on
contour(x1,y1,kyfun(x1,y1),[0,0],'b'), hold off
[xsol,ysol]=newton2d(kx, ky, -0.8, 0.8);
[xsol,ysol]
[xsol,ysol]=newton2d(kx, ky, 0.8, -0.8);
[xsol,ysol]
[xsol,ysol]=newton2d(kx, ky, 0, 0);
[xsol,ysol]
[xsol,ysol]=newton2d(kx, ky, -0.5, -0.5);
[xsol,ysol]
[xsol,ysol]=newton2d(kx, ky, 0.5, 0.5);
[xsol,ysol]
hessdetk = simplify(diff(kx, x)*diff(ky, y) - diff(kx, y)^2)
double(subs(hessdetk, [x,y], [-3^(1/2)/2, 3^(1/2)/2]))
double(subs(hessdetk, [x,y], [3^(1/2)/2, -3^(1/2)/2]))
double(subs(hessdetk, [x,y], [0, 0]))
double(subs(hessdetk, [x,y], [-0.5, -0.5]))
double(subs(hessdetk, [x,y], [0.5, 0.5]))
double(subs(diff(kx, x), [x,y], [-3^(1/2)/2, 3^(1/2)/2]))
double(subs(diff(kx, x), [x,y], [3^(1/2)/2, -3^(1/2)/2]))
double(subs(diff(kx, x), [x,y], [0,0]))
figure(2);
ezcontour(k, [0.75, 0, 0.75, 0])
figure(3);
ezsurf(k, [-1.5, 1.5, -1.5, 1.5])
Calculation results
Critical points (coordinates where kx = ky = 0)
(-0.8, 0.8)
−
3 /2 ,
3/2
(0.8, -0.8)
3/2 ,−
3/2
(0,0) (0,0)
(-0.5, -0.5) (-0.5, -0.5)
(0.5, 0.5) (0.5, 0.5)
Left values are the values we estimated in the plot, and the right are the exact values according to newton2d function.
Hessian determinant
144x
2
y
2
−24x
2
−24y
2
3
Hessian determinants calculated for critical points
−
3 /2 ,
3/2
48
3/2 ,−
3/2
48
(0,0) 3
(-0.5, -0.5) 0
(0.5, 0.5) 0
We have three local min- or maxpoints and two other points. To find out whether it is min or max, we have to look at
the second derivatives.
Classified critical points
We put the first three points into any second derivate, to find out if they are minimum or
maximum points.
We find the first two to be positive indicating that they're local minimum points. The third one
has a negative secondary derivative, which means that it's a local maximum.
Plots
Figure 3. ezsurf
Figure 4. contour (kx=0: red ky=0: blue)
Figure 5. Ezcontour of one of the two undefined points (the other one is mirrored).
Conclusion
To define the two last points, we look at the contour plot above. We can see that they're saddle
points.
−
3 /2 ,
3/2
min
3/2 ,−
3/2
min
(0,0) max
(-0.5, -0.5) saddle
(0.5, 0.5) saddle
Problem 4
Find and classify all critical points of the function
x
4
y
4
z
4
xyz.
MATLAB will report many critical points, but only a few of them are real.
Solution procedure
To find critical points of the function, we must set the partial derivatives equal to 0 and solve
for x, y and z. Then we calculate the Hessian matrix for the second derivatives.
We discarded the non-real critical points.
MATLAB calculations
MATLAB code
syms x y z
f=(x^4 + y^4 + z^4 + x*y*z)
gradf = jacobian(f,[x,y,z])
hessmat = jacobian(gradf, [x,y,z])
[xcr, ycr, zcr] = solve(gradf(1),gradf(2),gradf(3));
[xcr, ycr, zcr]
H1=subs(hessmat,[x,y,z],[xcr(1), ycr(1), zcr(1)])
H2=subs(hessmat,[x,y,z],[xcr(2), ycr(2), zcr(2)])
H3=subs(hessmat,[x,y,z],[xcr(3), ycr(3), zcr(3)])
H6=subs(hessmat,[x,y,z],[xcr(6), ycr(6), zcr(6)])
H7=subs(hessmat,[x,y,z],[xcr(7), ycr(7), zcr(7)])
eigH1 = double(eig(H1))
eigH2 = double(eig(H2))
eigH3 = double(eig(H3))
eigH6 = double(eig(H6))
eigH7 = double(eig(H7))
Calculation results
Partial derivatives
f '
x
=4x
3
yz
f '
y
=4y
3
xz
f '
z
=4z
3
xy
Hessian matrix
H x , y , z=
12x
2
z y
z 12y
2
x
y x 12z
2
Critical points
x y z
0 0 0
-¼ ¼ ¼
-¼ -¼ -¼
¼ -¼ ¼
¼ ¼ -¼
Eigenvalues
0 0 0
0.250 1 1
0.250 1 1
0.250 1 1
0.250 1 1
Conclusion
Our eigenvalues determine the critical points as being minimum except the first one (0, 0, 0)
which is degenerate.
Problem 5
Find and classify all critical points of the function
h x , y =y
2
⋅exp x
2
−x −3y.
You will need the graphical/numerical method to find the critical points.
Solution procedure
We plot f'x = 0 and f'y = 0 and estimated the critical points where f'x = f'y = 0. We calculate
the Hessian matrix, and puts in the values returned from the solve function. Then we calculate
the eigenvalues to see which type of point it is.
MATLAB calculations
MATLAB code
syms x y
h = y^2*exp(x^2) - x - 3*y
hx = diff(h, x)
hy = diff(h, y)
[x1,y1]=meshgrid(0:.05:5,-.5:.05:1.5);
hxfun=inline(vectorize(hx));
hyfun=inline(vectorize(hy));
figure(1);
contour(x1,y1,hxfun(x1,y1),[0,0],'r'), hold on
contour(x1,y1,hyfun(x1,y1),[0,0],'b'), hold off
gradh = jacobian(h,[x,y])
hessmat = jacobian(gradh, [x,y])
[xcr, ycr] = solve(gradh(1),gradh(2));
[xcr, ycr]
H1=subs(hessmat,[x,y],[xcr(1), ycr(1)])
eigH1 = double(eig(H1))
H2=subs(hessmat,[x,y],[xcr(2), ycr(2)])
eigH2 = double(eig(H2))
figure(2);
ezsurf(h, [0.23, 0.26, 1.4, 1.45])
figure(3);
ezsurf(h, [1.33, 1.35, 0.24, 0.26])
Calculation results
Partial derivatives
h '
x
=2y
2
x⋅exp x
2
−1
h '
y
=2y⋅exp x
2
−3
Estimated coordinates for critical points
x y
0.25 1.42
1.34 0.25
Eigenvalues
H1 = (5.3429, 1.4987)
H2 = (16.8761, -1.3836)
Where the first one is a minimum and the other a saddle point.
Plots
Figure 6. Contour (f'x = 0: red; f'y = 0: blue)
Figure 7. Minimum point
Figure 8. Saddle point.
Conclusion
The plots verifies that there is a minimum and a saddle point.
Problem 6
Find and classify all critical points of the function
h x , y = y
2
⋅exp x− x−3y.
Use both the analytical and the graphical/numerical methods to find the critical points, and
compare the results.
Solution procedure
Analytical
To find critical points of h, we must set the partial derivatives equal to 0 and solve for x and y.
We calculate the second derivatives and put them into the Hessian matrix. Then we calculate
the eigenvalues which will classify our critical points.
Graphical/numerical
With a contour plot of the derivatives, we are able to narrow down the input range close to the
cutting points. Then we calculate the exact values with the newton2d function.
MATLAB calculations
MATLAB code
syms x y
h = y^2*exp(x) - x - 3*y
hx = diff(h, x)
hy = diff(h, y)
[x1,y1]=meshgrid(-5:.01:10,-5:.01:5);
hxfun=inline(vectorize(hx));
hyfun=inline(vectorize(hy));
figure(1);
contour(x1,y1,hxfun(x1,y1),[0,0],'r'), hold on
contour(x1,y1,hyfun(x1,y1),[0,0],'b'), hold off
[xsol,ysol]=newton2d(hx, hy, 0.81, 0.67);
[xsol,ysol]
gradh = jacobian(h,[x,y])
hessmat = jacobian(gradh, [x,y])
[xcr, ycr] = solve(gradh(1),gradh(2));
[xcr, ycr]
H1=subs(hessmat,[x,y],[xcr(1), ycr(1)])
eigH1 = double(eig(H1))
figure(2);
ezsurf(h, [log(9/4)-0.1, log(9/4)+0.1, 2/3-0.1, 2/3+0.1])
Calculation results
Critical point
(log(9/4), 2/3)
Hessian matrix
H x , y=
y
2
⋅e
x
2y⋅e
x
2y⋅e
x
2⋅e
x
H log9/4 ,2/3=
1 3
3 9/2
Eigenvalues
6.2231
-0.7231
Plots
Figure 9. Saddle point.
Conclusion
The graphical (newton2d) and the analytical (solve) methods returns the same results. Mixed signs in the
eigenvalues results in a saddle point at the critical point, just as the plot show us.
