Python实现拉格朗日求解最优问题【椭圆内接长方体的最大体积】
拉格朗日乘数法的基本思想
作为一种优化算法，拉格朗日乘子法主要用于解决约束优化问题，它的基本思想就是通过引入拉格朗日乘子来将含有n个变量和k个约束条件的约束优化问题转化为含有（n+k）个变量的无约束优化问题。拉格朗日乘子背后的数学意义是其为约束方程梯度线性组合中每个向量的系数。
如何将一个含有n个变量和k个约束条件的约束优化问题转化为含有（n+k）个变量的无约束优化问题？拉格朗日乘数法从数学意义入手，通过引入拉格朗日乘子建立极值条件，对n个变量分别求偏导对应了n个方程，然后加上k个约束条件（对应k个拉格朗日乘子）一起构成包含了（n+k）变量的（n+k）个方程的方程组问题，这样就能根据求方程组的方法对其进行求解。
更多详细说明请参考：Math & Algorithm] 拉格朗日乘数法
一、问题描述
给定椭圆 x 2 a 2 + y 2 b 2 + z 2 c 2 = 1 \frac{x^2}{a^2} + \frac{y^2}{b^2} + \frac{z^2}{c^2} = 1 a2x2+b2y2+c2z2=1 ，求这个椭球的内接长方体的最大体积。
这个问题实际上就是条件极值问题，即在条件 x 2 a 2 + y 2 b 2 + z 2 c 2 = 1 \frac{x^2}{a^2} + \frac{y^2}{b^2} + \frac{z^2}{c^2} = 1 a2x2+b2y2+c2z2=1下，求 f ( x , y , z ) = 8 x y z f(x,y,z)=8xyz f(x,y,z)=8xyz 的最大值。
二、拉格朗日手工求解问题
① 带入拉格朗日乘子法中
F
(
x
,
y
,
z
,
λ
)
=
f
(
x
,
y
,
z
)
+
λ
ψ
(
x
,
y
,
z
)
=
8
x
y
z
+
λ
(
x
2
a
2
+
y
2
b
2
+
z
2
c
2
−
1
)
F(x,y,z,\lambda) = f(x,y,z) + \lambda\psi(x,y,z) \\ =8xyz + \lambda(\frac{x^2}{a^2} + \frac{y^2}{b^2} + \frac{z^2}{c^2} -1)
F(x,y,z,λ)=f(x,y,z)+λψ(x,y,z)=8xyz+λ(a2x2+b2y2+c2z2−1)
② 对 $F(x,y,z,\lambda) $求偏导
∂
F
(
x
,
y
,
z
,
λ
)
∂
x
=
8
y
z
+
2
λ
x
a
2
=
0
∂
F
(
x
,
y
,
z
,
λ
)
∂
y
=
8
x
z
+
2
λ
y
b
2
=
0
∂
F
(
x
,
y
,
z
,
λ
)
∂
z
=
8
x
y
+
2
λ
z
c
2
=
0
∂
F
(
x
,
y
,
z
,
λ
)
∂
λ
=
x
2
a
2
+
y
2
b
2
+
z
2
c
2
−
1
=
0
\frac{\partial F(x,y,z,\lambda)}{\partial x} = 8yz + \frac{2\lambda x}{a^2} = 0 \\ \frac{\partial F(x,y,z,\lambda)}{\partial y} = 8xz + \frac{2\lambda y}{b^2} = 0 \\ \frac{\partial F(x,y,z,\lambda)}{\partial z} = 8xy + \frac{2\lambda z}{c^2} = 0 \\ \frac{\partial F(x,y,z,\lambda)}{\partial \lambda} = \frac{x^2}{a^2} + \frac{y^2}{b^2} + \frac{z^2}{c^2} -1 = 0
∂x∂F(x,y,z,λ)=8yz+a22λx=0∂y∂F(x,y,z,λ)=8xz+b22λy=0∂z∂F(x,y,z,λ)=8xy+c22λz=0∂λ∂F(x,y,z,λ)=a2x2+b2y2+c2z2−1=0
解方程可得：
x
=
3
3
a
y
=
3
3
b
z
=
3
3
c
x=\frac{\sqrt{3}}{3}a \qquad y=\frac{\sqrt{3}}{3}b \qquad z=\frac{\sqrt{3}}{3}c
x=33
ay=33
bz=33
c
③ 带入原始目标函数，求解最大体积
V
m
a
x
=
f
(
3
3
a
,
3
3
b
,
3
3
c
)
=
8
3
9
a
b
c
V_{max} = f(\frac{\sqrt{3}}{3}a, \frac{\sqrt{3}}{3}b, \frac{\sqrt{3}}{3}c)=\frac{8\sqrt{3}}{9}abc
Vmax=f(33
a,33
b,33
c)=983
abc
三、Python 编程实现
实验环境
Anaconda + python3.6 + jupyter
① 导入 sympy 库
from sympy import *
② 设置变量
x,y,z,k = symbols('x,y,z,k')
a,b,c=symbols('a,b,c')
f = 8*x*y*z
g = x**2/a**2+y**2/b**2+z**2/c**2-1
③ 构造拉格朗日等式
L=f+k*g
④ 求偏导
dx = diff(L, x) # 对x求偏导
print("dx=",dx)
dy = diff(L,y) #对y求偏导
print("dy=",dy)
dz = diff(L,z) #对z求偏导
print("dz=",dz)
dk = diff(L,k) #对k求偏导
print("dk=",dk)
④ 求解
m= solve([dx,dy,dz,dk],[x,y,z,k])
print(m)
根据题目要求，最终求得的 x , y , z x,y,z x,y,z 应该为正数，所以标注出来的应为对应变量的解
⑤ 求最大体积
#变量赋值
x=sqrt(3)*a/3
y=sqrt(3)*b/3
z=sqrt(3)*c/3
#计算方程的值
f = 8*x*y*z
print("内接长方体的最大体积：",f)