Lagrange Multipliers with Maple

Copyright © 2001, 2002 by James F. Hurley, University of Connecticut Department of Mathematics, Unit 3009, Storrs CT 06269-3009. All rights reserved.

The interactive version of this worksheet is in MSB 203 in the Math Lab folder inside the Server1 volume.

The method of Lagrange multipliers applies to constrained optimization problems, in which the object is find the maximum and minimum of a function
*f *
:
** **
®
** **
on a level surface
*g*
(
*x*
,
*y*
,
*z*
) =
*k*
in
. Maple's built-in routine for solving systems of equations is often helpful for such problems, because Lagrange's method involves solving a system of equations

(*)
** Ñ**
*f*
(
*x*
,
*y*
,
*z*
) =
**Ñ**
*g*
(
*x*
,
*y*
,
*z*
),

*g*
(
*x*
,
*y*
,
*z*
) =
*k*

in the variables
*x*
,
*y*
,
*z*
, and
. Since both gradient vectors have three coordinates, this is a system of four equations in four unknowns. Hand solution of such a system usually requires some insight in eliminating variables, and the complexity of the associated algebra can be daunting. Maple's
solve
command can't solve
*every*
such system, but it can handle virtually all exercises in texts, whose authors designed them to be solvable by hand.

**Example 1**
. As of June, 2001, FedEx Corporation's size limitations on rectangular boxes are maximum dimension (called the
*length*
of the box) at most 119 inches, and girth (length plus the perimeter of the ends) at most 165 inches. What is the box of maximum volume that FedEx will accept for shipment?

Solution.
Let
*z*
be the maximum dimension of the box in inches, and let
*x*
and
*y*
denote the other two dimensions in inches. Then the problem is to

maximize
*V*
=
*f*
(
*x*
,
*y*
,
*z*
) =
*xyz*

subject to

0
* ² z*
² 119, 0 ²
*x*
²
*z*
, 0 ²
*y*
²
*z*
, and 2
*x*
+ 2
*y*
+
*z*
² 165.

To get the most out of the constrained variables
*x*
,
*y*
and
*z*
, replace the last inequality by the corresponding
*equation*
. That defines the level surface
*g*
(
*x*
,
*y*
,
*z*
)
= 165 for the function
*g*
, where

.

This problem is feasible for hand solution, but Maple provides a handy way to check the solution
*x*
=
*y*
= 27.5,
*z*
= 55 from hand implementation of the Lagrange-multiplier method. The routine below uses Maple's partial differentiation commands
diff(f(x, y), x)
,
diff(f(x, y), y)
and
diff(f(x, y), z)
to
produce the system of equations (*) above. (This capability is valuable in case of a discrepancy between a hand-computed answer and an answer book: you can check all the partial differentiation in the hand solution.) Note that since this routine doesn't put the partial derivatives together into the gradient vector, invoking the
linalg
package is not necessary. Execute the code to see Maple's answers.

`> `
**f := (x, y, z) -> x*y*z;
g := (x, y, z) -> 2*x + 2*y + z;
f_x := (x, y, z) -> diff(f(x, y, z), x):
f_y := (x, y, z) -> diff(f(x, y, z), y):
f_z := (x, y, z) -> diff(f(x, y, z), z):
evalf(f_x(x, y, z));
evalf(f_y(x, y, z));
evalf(f_z(x, y, z));
g_x := (x, y, z) -> diff(g(x, y, z), x):
g_y := (x, y, z) -> diff(g(x, y, z), y):
g_z := (x, y, z) -> diff(g(x, y, z), z):
evalf(g_x(x, y, z));
evalf(g_y(x, y, z));
evalf(g_z(x, y, z));
**

The next routine asks Maple to solve the system of equations (*) for the
*f*
and
*g*
in this problem.

`> `
**solve( { evalf(f_x(x, y, z))= lambda*evalf(g_x(x, y, z)),
evalf(f_y(x, y, z))= lambda*evalf(g_y(x, y, z)),
evalf(f_z(x, y, z))= lambda*evalf(g_z(x, y, z)),
g(x, y, z) = 165});**

Notice that Maple finds all solutions to the system, even those for which a variable is 0, which of course are not relevant in designing an actual box. Only the last solution has all positive dimensions, and it does agree with the solution obtained by hand calculation. The following example checks the hand solution to an example in the forthcoming second edition of
*Multivariable Calculus*
by the author of this worksheet.

**Example 2**
. Find the maximum and minimum of
*f*
(
*x*
,
*y*
) =
subject to the constraint
.

**Solution. **
First, find the extreme values on the bounding circle
of the feasible region, which is the circular disk of radius 4 centered at the origin. To finish the problem, investigate the behavior of
*f*
inside the feasible region. Proceeding as before, re-code the formulas for
*f*
and
*g*
,
and then execute the edited routine, which of course involves just the two variables
*x*
and
*y*
.

`> `
**f := (x, y) -> 2*x - y ;
g := (x, y) -> x^2 + y^2 - 16;
f_x := (x, y, z) -> diff(f(x, y), x):
f_y := (x, y, z) -> diff(f(x, y), y):
evalf(f_x(x, y));
evalf(f_y(x, y));
g_x := (x, y) -> diff(g(x, y), x):
g_y := (x, y) -> diff(g(x, y), y):
evalf(g_x(x, y));
evalf(g_y(x, y));
**

Then, have Maple solve the system of equations

** Ñ**
*f*
(
*x*
,
*y*
) =
**Ñ**
*g*
(
*x*
,
*y*
),

.

Execute the code below to see Maple's solutions to this system.

`> `
**(ans1, ans2) := solve( {
evalf(f_x(x, y))= lambda*evalf(g_x(x, y)),
evalf(f_y(x, y))= lambda*evalf(g_y(x, y)),
g(x, y) = 0
} );
**

`> `
**subs(ans1, f(x, y));
**

`> `
**subs(ans2, f(x, y));**

From this, Maple's answers for extreme-value points and the corresponding maximum and minimum of
*f*
subject to the two constraints are clear. How do they compare to the values found by hand:

¥ minimum
at
*x*
=
,
*y*
=
,

¥ maximum
at the point with coordinates that are the negatives of those values?

Maple converts those numbers to floating-point numbers that agree with the output from its
solve
command except in the last decimal place. (That discrepancy is the result of floating-point roundoff error.)

`> `
**evalf([-8/sqrt(5), 4/sqrt(5), -4*sqrt(5)]);
**

To identify the maximum and minimum of
*f*
on the entire closed disk , note that the gradient is never 0, so there can be no local maximum or minimum inside the circle. Thus, the values found above are the actual maximum and minimum of
*f*
on the entire feasible region.

The concluding example below (Example 10.7, Section 4.10 of the current edition of
*Multivariable Calculus*
) illustrates the Lagrange-multiplier method for optimizing a function subject to more than one constraint. .

**Example 3**
. Find the maximum and minimum of
on the curve of intersection of the cylinders
and
.

**Solution**
.
Proceeding as before, re-code the formulas in Exampl 1 for
and
, add the second constraint function
*h*
(
*x*
,
*y*
,
*z*
) =
to the first routine, and then execute the edited routine:

`> `
**f:= (x, y, z) -> x*y + x*z;
g := (x, y, z) -> x^2 + y^2 - 1;
h := (x, y, z) -> x*z - 1;
f_x := (x, y, z) -> diff(f(x, y, z), x):
f_y := (x, y, z) -> diff(f(x, y, z), y):
f_z := (x, y, z) -> diff(f(x, y, z), z):
evalf(f_x(x, y, z));
evalf(f_y(x, y, z));
evalf(f_z(x, y, z));
g_x := (x, y, z) -> diff(g(x, y, z), x):
g_y := (x, y, z) -> diff(g(x, y, z), y):
g_z := (x, y, z) -> diff(g(x, y, z), z):
evalf(g_x(x, y, z));
evalf(g_y(x, y, z));
evalf(g_z(x, y, z));
h_x := (x, y, z) -> diff(h(x, y, z), x):
h_y := (x, y, z) -> diff(h(x, y, z), y):
h_z := (x, y, z) -> diff(h(x, y, z), z):
evalf(h_x(x, y, z));
evalf(h_y(x, y, z));
evalf(h_z(x, y, z));**

Then, use the following routine to have Maple solve the system of four equations in the four variables
*x*
,
*y*
,
*z*
and
.
(If Maple gives an error message about multiple assignments, try adding or deleting answers: here we initially guess that there may be as many as four, which turns out to be accurate.)

** Ñ**
*f*
(
*x*
,
*y*
,
*z*
) =
**Ñ**
*g*
(
*x*
,
*y*
,
*z*
) +
**Ñ**
*g*
(
*x*
,
*y*
,
*z*
)

,
= 1.

`> `
**(ans1, ans2, ans3, ans4) := solve( {
evalf(f_x(x, y, z))= lambda*evalf(g_x(x, y, z)) + mu*evalf(h_x(x, y, z)),
evalf(f_y(x, y, z))= lambda*evalf(g_y(x, y, z)) + mu*evalf(h_y(x, y, z)),
evalf(f_z(x, y, z))= lambda*evalf(g_z(x, y, z)) + mu*evalf(h_z(x, y, z)), g(x, y, z) = 0, h(x, y, z) = 0
} );**

`> `
**subs(ans1, f(x, y, z));**

`> `
**subs(ans2, f(x, y, z));
**

`> `
**subs(ans3, f(x, y, z));
**

`> `
**subs(ans4, f(x, y, z));**

From this, Maple's answers for extreme-value points and the corresponding maximum
*f*
subject to the two constraints are clear. The minimum value appears so close to 0.5 that the discrepancy appears due to floating-point roundoff error.

How do these answers compare to the values found by hand? As you can check (or read in the text), hand computation gives minimum
at
*x*
=
,
*y*
=
,
*z*
=
and also at
; the same hand work yields maximum value
at the points
and
. As the following shows, Maple converts those numbers to floating-point numbers that substantially agree with the output from the
solve
command.

`> `
**evalf([-1/sqrt(2), 1/sqrt(2), -sqrt(2), 1/2]);**