工程优化设计方法研究报告
一、三种约束方法概述
1.1 Zoutendijk’s Method概述
1.1.1算法思想
Zoutendijk’s Method称为最速可行方向法。可行方向法是求解如下不等式约束问题：
（1-1）
的一种直接算法，这种算法的思路是：从可行域内的一个可行点出发，选择一个合适的搜索方向
和步长因子
，使产生的下一个迭代点
（1-2）
既不超出可行域，又使目标函数的值下降得尽可能多。也就是说，使新的迭代点同时满足
和 （1-3）
可以断定，从可行域内的任意初始点出发，只要始终沿着下降可行方向进行考虑约束限制的一维搜索和迭代运算，就能保证迭代点逐步逼近约束最优化问题的最优点。
在一个点的所有下降可行方向中，使目标函数取得最大下降量的方向称为最佳下降可行方向。首先找到最佳下降可行方向，然后在这个方向上进行一维搜索，通过迭代，就能找出最优点。
1.1.2理论基础
函数在一点的梯度方向是函数在该点函数值上升得最快的方向，与梯度值成锐角的方向是函数值上升的方向，与梯度成钝角的方向是函数值下降的方向。因此，函数在点的下降方向就是满足以下关系的方向S：
（1-4） 可行方向是那些指向可行域内的方向。当点
位于可行域内时，从该点出发的任意方向
S上都必然存在满足：
（1-5） 的可行点，因此所有方向都是可行方向。
当点
位于某一起作用约束的边界
，为了满足上述条件，必须使方向S
和对应的约束函数在该点的梯度相交成钝角，即使
（1-6）
因为约束条件定义为
的形式，所以约束函数的梯度
是指向可行
域之外的。并且在该梯度方向上约束函数的值增大得最快。故与该梯度成钝角的方向必定是指向可行域内的可行方向。
最佳下降可行方向可以在满足可行条件的前提下，通过极小化目标函数的方向导数得到。由此构成如下寻求最佳下降可行方向S的最优化问题：
（1-7）
上式中，
由于点
上的梯度
。
和
都是已知的常数向量，所以上式中
的线性函数。因此这个问题变成了线性规划
。
的目标函数和约束函数都是变量
问题。利用求解线性规划问题的单纯形算法可以方便得求出最佳下降可行方向
1.1.3算法实现
1：选取初始数据，选取初始可行点x0，允许误差ε>0，令k=0；
2：确定起作用约束，确定点xk外起作用约束的指标集，I(xk)=|i|gi(xk)≥0; i=1,2⋯m
3：求解辅助线性规划问题，求解线性规划问题：
minz
s.t ∇f(xk)Td−z≤0
的最优解（zk,dk） T
∇g(x)d+z≥0,i∈I()ikxk
{−1≤dj≤1,j=1,2⋯n
4：构造可行下降方向，若|zk|≤ε，停止迭代，输出xk；否则，得到可行下降方向dk，进入下一步；
5：进行一维搜索，求解min
0≤λ≤λmax
f(xk+λdk)，得到最优解dk。其中λmax按照下式计
算：
λmax=sup|λ|gi(xk+λdk)≥0,i=1,2⋯m;
令xk+1=xk+λdk，k=k+1，返回第二步。
注意：计算实践和理论分析表明，该算法可能失效或者出现锯齿现象，使算法可能收敛很慢甚至不收敛到最优点或K—T点。
1.1.4适用范围
这种方法主要对线性不等式约束较为有效，对非线性不等式约束优化问题效果不佳。
1.2Rosen’s Gradient Projection Method方法概述
1.2.1算法思想
Rosen’s Gradient Projection Method是利用梯度的投影技巧求约束非线性规划问题最优解的一种方法。该方法可以处理一般的非线性规划问题，但在实际上它最适宜于处理约束条件为线性函数的非线性规划问题。在一般工程中，功能函数呈现线性状态，采用Rosen’s Gradient Projection Method方法能够有效的处理非线性功能函数，而且计算效率和精度高，通用性强。
Rosen’s Gradient Projection Method的基本思想为：当迭代点在可行域内部时，取该点处的负梯度方向为可行下降方向；当迭代点在可行域边界上时，取该点处负梯度方向在可行域边界上的投影产生一个可行下降方向。
1.2.2理论基础
梯度投影算法并不通过进行线性规划求解，而是通过将目标函数的负梯度投影到当前的有效约束上。问题描述为：
满足： （1-8）
设有效约束为
，那么有效约束的梯度即为：
（1-9）
定义一个n x p的矩阵 （1-10）
那么寻找可行方向S的问题可以被转化为以下问题：
找到一个S，使得最小。并满足以下条件：
（1-11）
为了解决上述的等式约束问题，我们构造Lagrangian函数如下：
（1-12） 通过上述转换，将优化条件变为下式：
（1-13）
求解S得到：
（1 （1-14）
如果S被标准化，
不为0，由此：
（1-15） 所以：
（1-16） 其中P就是投影矩阵。忽略常数束向量的正交平面的相交面上。如果则
，可以看出，投影矩阵将向量
投影到约
是独立的，那么矩阵N的列向量线性无关，
是非奇异矩阵并可以被转置。向量S可以被标准化如下：
（1-17） 设
是第i次迭代的初始点。(
严格满足)，我们可以找到
：
（1-18） 是：
在矩阵P下的投影。如果那么从出发，沿着方向可以找到点
（1-19） 如果
那么：
（1-20） 其中
（1-21） 上式表明，目标函数的负梯度可以通过在
处的有效约束的线性组合给出。如果所有的
都是非负的，那么Kuhn-Tucker条件满足，迭代可以被终止。 如果一些
是负的，且
则一些约束梯度
方向与
在
处的角度
为负
是钝角，这表明这些约束虽然是有效约束但是并不应该包括在内。在使用时，当
时，并不需要忽略所有的有效约束，新的矩阵N为：
（1-22）
由此计算投影矩阵P：
（1-23） 搜索方向为：
（1-24）
1.2.3算法实现
1：选取初始数据，选取初始可行点x0，允许误差允许误差ε>0，令k=0； Ab
2：分解出起作用的约束，在xk处分解A=[1]，b=[1]，使得A1xk=b1，A2xk=
A2b2b2
A
3：确定投影矩阵，令M=[1]，若M为空，则令P=In，否则令P=In−
E
MT(MMT)−1M
4：构造可行下降方向，令dk=−P∇f(xk)，若‖dk‖≤ε，则进行第5步，否则，进行第6步
5：检查起作用的约束数，判定迭代点xk是否在可行域内部，若M为空则停止迭代，u
输出xk，否则令W=(MMT)−1M，∇f(xk)=[]，其中u和v分别对应于A1和E，如果u≥
v0，则停止迭代，输出xk；如果，u
̂=b2−A2xk，d̂=A2dk，求解min6：进行一维搜索，计算b
0≤λ≤λmax
f(xk+λdk)
，得到
最优解λk，其中
̂ib
̂i
令xk+1
̂≥0时∞，当d
=xk+λdk，k=k+1，转到第2步。
1.2.4 适用范围
这种方法主要对线性不等式约束较为有效，对非线性不等式约束优化问题效果不佳。
1.3Generalized Reduced Gradient Method方法概述
1.3.1算法思想
1963年,Wolfe将线性规划的单纯形法推广到具有非线性目标函数的问题，提出了产生可行下降方向的另一类方法，称为即约梯度法。随后推广到求解带非线性等式约束的情形，即现在的Generalized Reduced Gradient Method，该方法是目前求解约束非线性最优化问题的最有效的方法之一。
Generalized Reduced Gradient Method的基本思想：利用约束条件将所考察问题中某些变量用其它一组独立变量来表示，如对线性约束，将非基变量当作独立变量，把基变量用非基变量表示，消去目标函数中的基变量，那么目标函数在低维（n－m维）空间中的负梯度方向便是原问题的一个下降方向，再考虑到可行性，就可确定一个下降可行方向。
1.3.2理论基础
广义简约梯度法（GRG）是简约梯度法的延伸，要求
的最小值，约束为：
（1-25）
通过在不等式约束上添加一个非负的松弛变量，将约束转化为：
（1-26）
由此可以将问题写做有n+m个变量
的形式：
（1-27）
GRG方法的算法思想是对等式约束中的变量进行消元。从理论上来说，一个方程可以解出一个变量，可以将n+m个变量分为两个系列：
（1-28）
其中Y是设计变量，X是状态变量。设计变量是完全独立的，状态变量是非独立的，利用状态变量去满足约束
对于第一个目标函数和约束函数而言：
其中
（1-29）
即那个初始X满足约束，对于向量X的任何改变，必须满足，这样，X+dX处
约束才是满足的。所以
（1-30） 进而有：
（1-31）
上式被称为广义简约梯度，从几何学上说，简约梯度可以描述为将原始的n维梯度到n-m维可行域上的一个投影。
要使
（1-32） 必须有
（1-33） 所以
（1-34）
1.3.3算法实现
1：选取初始数据，选取初始可行点x0，允许误差允许误差ε>0，令k=0；
2：进行矩阵分解，从xk中选取m个最大分量，其下标集记作Jk，令B是由下标属于Jk的A构成的m阶方阵，余下的列构成矩阵N。
∇xf(xk)(k)
3：求即约梯度，计算∇f(xk)=[B]，求r{xN}=∇xNf(xk)−
∇xNf(xk)(B−1N)T∇xBf(xk)
4：构造可行下降方向 dNj
(k)
kk
−xNkr{x}，当r{x}>0时jjNNj={ (k)(k)
−rj{xN}，当rj{xN}≤0时
(k)()
()
()
dB=−B−1NdN
5：检查终止条件，若‖dk‖≤ε，则停止计算，输出xk；否则，转到第6步；
6：进行一维搜索，求解min
0≤λ≤λmax
(k)
f(xk+λdk)得到最优解λk，其中λmax按照下式计算： xjkdj
()
λmax=
min{ {
|dk
∞，当dk≥0时
()
1.3.4适用范围
这种方法对具有非线性不等式约束和等式约束效果都很好。 （1）、适合具有非线性约束的问题； （2）、不太适合具有非光滑约束的问题（需要线性逼近）； （3）、对于中小型，大型约束优化问题均适用，且计算稳定性好。
二、实例分析
选取第一种和第三种方法作为题示示例的解答
2.1模型的规范化形式比较
一般约束优化问题模型为：
（2-1）
Zoutendijk’s Method是将优化问题转化为了标准线性规划的问题，规范化模型中不含等式约束；Rosen’s Gradient Projection Method最终得到的是一个投影矩阵，是利用投影矩阵去找最优解；Generalized Reduced Gradient Method是将不等式约束问题转化为等式约束问题求解，规范化模型中只有等式约束。
2.２ Zoutendijk’s Method与Generalized Reduced Gradient Method计算比较
2.2.1Zoutendijk’s Method
见第二次作业。
2.2.2Generalized Reduced Gradient Method
初始点的选取：
function [x0,fx0,gx]=Random_Initial_Point(f,g_cons,xl,xu) %输入 f—目标函数 % g_cons—约束条件 % xl—最小极限 % xu—最大极限
%输出 x0—-随即采样生成初始点 % fx0—-初始点目标函数值 % gx—-初始点约束条件值 N=length(xl); M=size(g_cons); M=length(M(:,1)); gx=ones(M,1); while max(gx)>=0 dir0=rand(N,1); x0=xl+dir0.*xu; gx=feval(g_cons,x0); end
fx0=feval(f,x0); end
目标函数如下：
function f=func(x)
%设置目标函数
f=x(1)+x(2)+x(3);
% f=x(1)^2+x(1)*x(2)+2*x(2)^2-6*x(1)-2*x(2)-12*x(3); End
约束函数如下：
function [c ceq]=fun_cons(x) %·非线性／线性约束条件
c=[0.0025*(x(4)+x(6))-1 0.0025*(-x(4)+x(5)+x(7))-1 0.01*(-x(5)+x(8))-1
100*x(1)-x(1)*x(6)+833.33252*x(4)-83333.333 x(2)*x(4)-x(2)*x(7)-1250*x(4)+1250*x(5) x(3)*x(5)-x(3)*x(8)-2500*x(5)+1250000 100-x(1) x(1)-10000 1000-x(2) x(2)-10000 1000-x(3) x(3)-10000 10-x(4) x(4)-1000 10-x(5) x(5)-1000 10-x(6) x(6)-1000 10-x(7) x(7)-1000 10-x(8)
x(8)-1000];
ceq=[];
end
一阶偏导计算如下：
function [g]=fun_diff(fun,x,dh)
%计算一阶偏导
h2=2*dh;
M=length(fun);
N=length(x);
I=eye(N);
for n=1:N
g(:,n)=(feval(fun,x+I(n,:)*dh)-feval(fun,x-I(n,:)*dh))/h2; end
end
利用黄金分割法计算步长
function [X,Step,Error]=GoldenSectionMethod(f,a,b,minError,maxStep)
x1=a+0.382*(b-a);
x2=a+0.618*(b-a);
y1=feval(f,x1);
y2=feval(f,x2);
Error=b-a;
Step=0;
while (Error>minError)&&(maxStep>Step)
Step=Step+1;
if y1
b=x2;
x2=x1;
y2=y1;
x1=a+0.382*(b-a);
y1=feval(f,x1);
else
a=x1;
x1=x2;
y1=y2;
x2=a+0.618*(b-a);
y2=feval(f,x2);
end
Error=b-a;
end
X=(x1+x2)/2;
end
function lamda=fmin(func,fun_cons,x0,d,lamdmin,lamdmax)
%
lamda=GoldenSectionMethod(@(lamda)func(x0+lamda*d),lamdmin,lamdmax,1.0e-3,20);
lamda=fminbnd(@(lamda)func(x0+lamda*d),lamdmin,lamdmax);%求解初始步长
x=x0+lamda*d;
g_cons1=fun_cons(x0);
g_cons2=fun_cons(x);
g=max(g_cons2);% 判断初始步长是否满足约束条件
max_x=find(g_cons2==g);
while (g>1.0e-1)%对不满足约束条件的步长进行折算
lamda=(0-g_cons1(max_x))*lamda/(g_cons2(max_x)-
g_cons1(max_x));
x=x0+lamda*d;
g_cons2=fun_cons(x);
g=max(g_cons2);
max_x=find(g_cons2==g);
end
end
Generalized Reduced Gradient Method的算法文件：
function gradient()
clear
clc
syms d;
delta=0.001;
x0=init();xm=x0;
[yf,zf,C,D,gr]=gx(x0);normest(gr)
[d0,xx,S,T]=initd(C,D,gr,x0);
xkt=[];fxt=[];
for t=0:1:100
fmin=fff(x0); %目标函数值
xkt=[xkt x0']; %将搜索点输入至xkt
fxt=[fxt fmin]; %将目标函数值输入至fxt
fprintf('搜索点X： ');vpa(x0',6)
fprintf('搜索方向： ');[S;T]
fprintf('搜索点序号： ');t
fprintf('min f(x)： ');vpa(fmin,6)
if normest(gr)
fprintf('找到最优解');
break;
end %若满足终止条件，结束循环
dm=step(d0,xx);
fprintf('搜索步长λ： ');vpa(dm,6)
x=subs(xx,d,dm);
x0=x;
[yf,zf,C,D,gr]=gx(x0);
[d0,xx,S,T]=initd(C,D,gr,x0);
end
%绘图表示搜索过程
figure(1)
plot(xkt);
xlabel('搜索点序号');
ylabel('搜索点');
figure(2)
plot(fxt);
xlabel('搜索点序号');
ylabel('目标函数');
运行结果如下：
初始计算X0结果：
X0=
1.0e+03 *
2.[**************]
7.[**************]
7.[**************]
0.[**************]
0.[**************]
0.[**************]
0.[**************]
0.[**************]
最终计算结果方向：
-1.[**************]
0.[**************]
-1.[**************]
-1.[**************]
-1.[**************]
-1.[**************]
1.[**************]
-1.[**************]
最终计算结果步长：
1.[**************]e-08
最终计算F结果：
F=
1.[**************]e+04
最终计算X0结果：
X0=
1.0e+03 *
0.[**************]
6.[**************]
6.[**************]
0.[**************]
0.[**************]
0.[**************]
0.[**************]
0.[**************]
最终计算结果约束集：
5
在当前所选的初始点（随机生成的初始点），不能取得一个较好的结果，取迭代200次后
的值作为最后的结果。
X*=[100,1276.3,7948.4,115.5,225.4,261.5,238.8,316.4]' ，f*=9324.64 。
2.3结果比较与分析
Zoutendijk’s Method Generalized Reduced Gradient Method
Zoutendijk’s Method Generalized Reduced Gradient Method
通过上图比较可知，Zoutendijk’s Method计算过程中，目标函数值呈阶梯形下降，通过增加循环上限次数，最终结果将越来越逼近最优解；而Generalized Reduced Gradient Method计算过程中，目标函数值在经过几次迭代之后基本就保持不变了。
总而言之，Zoutendijk’s Method可能出现锯齿现象，当某一个迭代点在约束边界上时，如果可行方向选取不当，会使迭代点跑出边界，甚至使得算法不收敛。Rosen’s
Gradient Projection Method求解带线性约束的非线性规划问题十分有效。从一个基本可行解开始，由约束条件确定出凸约束集边界上梯度的投影，以便求出下次的搜索方向和步长。每次搜索后，都要进行检验，直到满足精度要求为止。Generalized Reduced Gradient Method收敛快，精度高，被认为求解带线性约束的非线性规划问题时有显著效果。
三、致谢
通过工程优化设计这门课，自己学习了关于最优化设计的一些算法与基础知识，也让
我了解到自己在设计上的不足，同时黄老师的授课与严谨的教学，也让我懂得了认真严谨
是作为一个科研工作者最基本的素质和态度，所以在此对黄老师表示真挚的感谢。
我也会在今后的学习中，脚踏实地，严谨治学，积极向老师和同学们请教与讨论，争取早日取得优秀的成果。
工程优化设计方法研究报告
一、三种约束方法概述
1.1 Zoutendijk’s Method概述
1.1.1算法思想
Zoutendijk’s Method称为最速可行方向法。可行方向法是求解如下不等式约束问题：
（1-1）
的一种直接算法，这种算法的思路是：从可行域内的一个可行点出发，选择一个合适的搜索方向
和步长因子
，使产生的下一个迭代点
（1-2）
既不超出可行域，又使目标函数的值下降得尽可能多。也就是说，使新的迭代点同时满足
和 （1-3）
可以断定，从可行域内的任意初始点出发，只要始终沿着下降可行方向进行考虑约束限制的一维搜索和迭代运算，就能保证迭代点逐步逼近约束最优化问题的最优点。
在一个点的所有下降可行方向中，使目标函数取得最大下降量的方向称为最佳下降可行方向。首先找到最佳下降可行方向，然后在这个方向上进行一维搜索，通过迭代，就能找出最优点。
1.1.2理论基础
函数在一点的梯度方向是函数在该点函数值上升得最快的方向，与梯度值成锐角的方向是函数值上升的方向，与梯度成钝角的方向是函数值下降的方向。因此，函数在点的下降方向就是满足以下关系的方向S：
（1-4） 可行方向是那些指向可行域内的方向。当点
位于可行域内时，从该点出发的任意方向
S上都必然存在满足：
（1-5） 的可行点，因此所有方向都是可行方向。
当点
位于某一起作用约束的边界
，为了满足上述条件，必须使方向S
和对应的约束函数在该点的梯度相交成钝角，即使
（1-6）
因为约束条件定义为
的形式，所以约束函数的梯度
是指向可行
域之外的。并且在该梯度方向上约束函数的值增大得最快。故与该梯度成钝角的方向必定是指向可行域内的可行方向。
最佳下降可行方向可以在满足可行条件的前提下，通过极小化目标函数的方向导数得到。由此构成如下寻求最佳下降可行方向S的最优化问题：
（1-7）
上式中，
由于点
上的梯度
。
和
都是已知的常数向量，所以上式中
的线性函数。因此这个问题变成了线性规划
。
的目标函数和约束函数都是变量
问题。利用求解线性规划问题的单纯形算法可以方便得求出最佳下降可行方向
1.1.3算法实现
1：选取初始数据，选取初始可行点x0，允许误差ε>0，令k=0；
2：确定起作用约束，确定点xk外起作用约束的指标集，I(xk)=|i|gi(xk)≥0; i=1,2⋯m
3：求解辅助线性规划问题，求解线性规划问题：
minz
s.t ∇f(xk)Td−z≤0
的最优解（zk,dk） T
∇g(x)d+z≥0,i∈I()ikxk
{−1≤dj≤1,j=1,2⋯n
4：构造可行下降方向，若|zk|≤ε，停止迭代，输出xk；否则，得到可行下降方向dk，进入下一步；
5：进行一维搜索，求解min
0≤λ≤λmax
f(xk+λdk)，得到最优解dk。其中λmax按照下式计
算：
λmax=sup|λ|gi(xk+λdk)≥0,i=1,2⋯m;
令xk+1=xk+λdk，k=k+1，返回第二步。
注意：计算实践和理论分析表明，该算法可能失效或者出现锯齿现象，使算法可能收敛很慢甚至不收敛到最优点或K—T点。
1.1.4适用范围
这种方法主要对线性不等式约束较为有效，对非线性不等式约束优化问题效果不佳。
1.2Rosen’s Gradient Projection Method方法概述
1.2.1算法思想
Rosen’s Gradient Projection Method是利用梯度的投影技巧求约束非线性规划问题最优解的一种方法。该方法可以处理一般的非线性规划问题，但在实际上它最适宜于处理约束条件为线性函数的非线性规划问题。在一般工程中，功能函数呈现线性状态，采用Rosen’s Gradient Projection Method方法能够有效的处理非线性功能函数，而且计算效率和精度高，通用性强。
Rosen’s Gradient Projection Method的基本思想为：当迭代点在可行域内部时，取该点处的负梯度方向为可行下降方向；当迭代点在可行域边界上时，取该点处负梯度方向在可行域边界上的投影产生一个可行下降方向。
1.2.2理论基础
梯度投影算法并不通过进行线性规划求解，而是通过将目标函数的负梯度投影到当前的有效约束上。问题描述为：
满足： （1-8）
设有效约束为
，那么有效约束的梯度即为：
（1-9）
定义一个n x p的矩阵 （1-10）
那么寻找可行方向S的问题可以被转化为以下问题：
找到一个S，使得最小。并满足以下条件：
（1-11）
为了解决上述的等式约束问题，我们构造Lagrangian函数如下：
（1-12） 通过上述转换，将优化条件变为下式：
（1-13）
求解S得到：
（1 （1-14）
如果S被标准化，
不为0，由此：
（1-15） 所以：
（1-16） 其中P就是投影矩阵。忽略常数束向量的正交平面的相交面上。如果则
，可以看出，投影矩阵将向量
投影到约
是独立的，那么矩阵N的列向量线性无关，
是非奇异矩阵并可以被转置。向量S可以被标准化如下：
（1-17） 设
是第i次迭代的初始点。(
严格满足)，我们可以找到
：
（1-18） 是：
在矩阵P下的投影。如果那么从出发，沿着方向可以找到点
（1-19） 如果
那么：
（1-20） 其中
（1-21） 上式表明，目标函数的负梯度可以通过在
处的有效约束的线性组合给出。如果所有的
都是非负的，那么Kuhn-Tucker条件满足，迭代可以被终止。 如果一些
是负的，且
则一些约束梯度
方向与
在
处的角度
为负
是钝角，这表明这些约束虽然是有效约束但是并不应该包括在内。在使用时，当
时，并不需要忽略所有的有效约束，新的矩阵N为：
（1-22）
由此计算投影矩阵P：
（1-23） 搜索方向为：
（1-24）
1.2.3算法实现
1：选取初始数据，选取初始可行点x0，允许误差允许误差ε>0，令k=0； Ab
2：分解出起作用的约束，在xk处分解A=[1]，b=[1]，使得A1xk=b1，A2xk=
A2b2b2
A
3：确定投影矩阵，令M=[1]，若M为空，则令P=In，否则令P=In−
E
MT(MMT)−1M
4：构造可行下降方向，令dk=−P∇f(xk)，若‖dk‖≤ε，则进行第5步，否则，进行第6步
5：检查起作用的约束数，判定迭代点xk是否在可行域内部，若M为空则停止迭代，u
输出xk，否则令W=(MMT)−1M，∇f(xk)=[]，其中u和v分别对应于A1和E，如果u≥
v0，则停止迭代，输出xk；如果，u
̂=b2−A2xk，d̂=A2dk，求解min6：进行一维搜索，计算b
0≤λ≤λmax
f(xk+λdk)
，得到
最优解λk，其中
̂ib
̂i
令xk+1
̂≥0时∞，当d
=xk+λdk，k=k+1，转到第2步。
1.2.4 适用范围
这种方法主要对线性不等式约束较为有效，对非线性不等式约束优化问题效果不佳。
1.3Generalized Reduced Gradient Method方法概述
1.3.1算法思想
1963年,Wolfe将线性规划的单纯形法推广到具有非线性目标函数的问题，提出了产生可行下降方向的另一类方法，称为即约梯度法。随后推广到求解带非线性等式约束的情形，即现在的Generalized Reduced Gradient Method，该方法是目前求解约束非线性最优化问题的最有效的方法之一。
Generalized Reduced Gradient Method的基本思想：利用约束条件将所考察问题中某些变量用其它一组独立变量来表示，如对线性约束，将非基变量当作独立变量，把基变量用非基变量表示，消去目标函数中的基变量，那么目标函数在低维（n－m维）空间中的负梯度方向便是原问题的一个下降方向，再考虑到可行性，就可确定一个下降可行方向。
1.3.2理论基础
广义简约梯度法（GRG）是简约梯度法的延伸，要求
的最小值，约束为：
（1-25）
通过在不等式约束上添加一个非负的松弛变量，将约束转化为：
（1-26）
由此可以将问题写做有n+m个变量
的形式：
（1-27）
GRG方法的算法思想是对等式约束中的变量进行消元。从理论上来说，一个方程可以解出一个变量，可以将n+m个变量分为两个系列：
（1-28）
其中Y是设计变量，X是状态变量。设计变量是完全独立的，状态变量是非独立的，利用状态变量去满足约束
对于第一个目标函数和约束函数而言：
其中
（1-29）
即那个初始X满足约束，对于向量X的任何改变，必须满足，这样，X+dX处
约束才是满足的。所以
（1-30） 进而有：
（1-31）
上式被称为广义简约梯度，从几何学上说，简约梯度可以描述为将原始的n维梯度到n-m维可行域上的一个投影。
要使
（1-32） 必须有
（1-33） 所以
（1-34）
1.3.3算法实现
1：选取初始数据，选取初始可行点x0，允许误差允许误差ε>0，令k=0；
2：进行矩阵分解，从xk中选取m个最大分量，其下标集记作Jk，令B是由下标属于Jk的A构成的m阶方阵，余下的列构成矩阵N。
∇xf(xk)(k)
3：求即约梯度，计算∇f(xk)=[B]，求r{xN}=∇xNf(xk)−
∇xNf(xk)(B−1N)T∇xBf(xk)
4：构造可行下降方向 dNj
(k)
kk
−xNkr{x}，当r{x}>0时jjNNj={ (k)(k)
−rj{xN}，当rj{xN}≤0时
(k)()
()
()
dB=−B−1NdN
5：检查终止条件，若‖dk‖≤ε，则停止计算，输出xk；否则，转到第6步；
6：进行一维搜索，求解min
0≤λ≤λmax
(k)
f(xk+λdk)得到最优解λk，其中λmax按照下式计算： xjkdj
()
λmax=
min{ {
|dk
∞，当dk≥0时
()
1.3.4适用范围
这种方法对具有非线性不等式约束和等式约束效果都很好。 （1）、适合具有非线性约束的问题； （2）、不太适合具有非光滑约束的问题（需要线性逼近）； （3）、对于中小型，大型约束优化问题均适用，且计算稳定性好。
二、实例分析
选取第一种和第三种方法作为题示示例的解答
2.1模型的规范化形式比较
一般约束优化问题模型为：
（2-1）
Zoutendijk’s Method是将优化问题转化为了标准线性规划的问题，规范化模型中不含等式约束；Rosen’s Gradient Projection Method最终得到的是一个投影矩阵，是利用投影矩阵去找最优解；Generalized Reduced Gradient Method是将不等式约束问题转化为等式约束问题求解，规范化模型中只有等式约束。
2.２ Zoutendijk’s Method与Generalized Reduced Gradient Method计算比较
2.2.1Zoutendijk’s Method
见第二次作业。
2.2.2Generalized Reduced Gradient Method
初始点的选取：
function [x0,fx0,gx]=Random_Initial_Point(f,g_cons,xl,xu) %输入 f—目标函数 % g_cons—约束条件 % xl—最小极限 % xu—最大极限
%输出 x0—-随即采样生成初始点 % fx0—-初始点目标函数值 % gx—-初始点约束条件值 N=length(xl); M=size(g_cons); M=length(M(:,1)); gx=ones(M,1); while max(gx)>=0 dir0=rand(N,1); x0=xl+dir0.*xu; gx=feval(g_cons,x0); end
fx0=feval(f,x0); end
目标函数如下：
function f=func(x)
%设置目标函数
f=x(1)+x(2)+x(3);
% f=x(1)^2+x(1)*x(2)+2*x(2)^2-6*x(1)-2*x(2)-12*x(3); End
约束函数如下：
function [c ceq]=fun_cons(x) %·非线性／线性约束条件
c=[0.0025*(x(4)+x(6))-1 0.0025*(-x(4)+x(5)+x(7))-1 0.01*(-x(5)+x(8))-1
100*x(1)-x(1)*x(6)+833.33252*x(4)-83333.333 x(2)*x(4)-x(2)*x(7)-1250*x(4)+1250*x(5) x(3)*x(5)-x(3)*x(8)-2500*x(5)+1250000 100-x(1) x(1)-10000 1000-x(2) x(2)-10000 1000-x(3) x(3)-10000 10-x(4) x(4)-1000 10-x(5) x(5)-1000 10-x(6) x(6)-1000 10-x(7) x(7)-1000 10-x(8)
x(8)-1000];
ceq=[];
end
一阶偏导计算如下：
function [g]=fun_diff(fun,x,dh)
%计算一阶偏导
h2=2*dh;
M=length(fun);
N=length(x);
I=eye(N);
for n=1:N
g(:,n)=(feval(fun,x+I(n,:)*dh)-feval(fun,x-I(n,:)*dh))/h2; end
end
利用黄金分割法计算步长
function [X,Step,Error]=GoldenSectionMethod(f,a,b,minError,maxStep)
x1=a+0.382*(b-a);
x2=a+0.618*(b-a);
y1=feval(f,x1);
y2=feval(f,x2);
Error=b-a;
Step=0;
while (Error>minError)&&(maxStep>Step)
Step=Step+1;
if y1
b=x2;
x2=x1;
y2=y1;
x1=a+0.382*(b-a);
y1=feval(f,x1);
else
a=x1;
x1=x2;
y1=y2;
x2=a+0.618*(b-a);
y2=feval(f,x2);
end
Error=b-a;
end
X=(x1+x2)/2;
end
function lamda=fmin(func,fun_cons,x0,d,lamdmin,lamdmax)
%
lamda=GoldenSectionMethod(@(lamda)func(x0+lamda*d),lamdmin,lamdmax,1.0e-3,20);
lamda=fminbnd(@(lamda)func(x0+lamda*d),lamdmin,lamdmax);%求解初始步长
x=x0+lamda*d;
g_cons1=fun_cons(x0);
g_cons2=fun_cons(x);
g=max(g_cons2);% 判断初始步长是否满足约束条件
max_x=find(g_cons2==g);
while (g>1.0e-1)%对不满足约束条件的步长进行折算
lamda=(0-g_cons1(max_x))*lamda/(g_cons2(max_x)-
g_cons1(max_x));
x=x0+lamda*d;
g_cons2=fun_cons(x);
g=max(g_cons2);
max_x=find(g_cons2==g);
end
end
Generalized Reduced Gradient Method的算法文件：
function gradient()
clear
clc
syms d;
delta=0.001;
x0=init();xm=x0;
[yf,zf,C,D,gr]=gx(x0);normest(gr)
[d0,xx,S,T]=initd(C,D,gr,x0);
xkt=[];fxt=[];
for t=0:1:100
fmin=fff(x0); %目标函数值
xkt=[xkt x0']; %将搜索点输入至xkt
fxt=[fxt fmin]; %将目标函数值输入至fxt
fprintf('搜索点X： ');vpa(x0',6)
fprintf('搜索方向： ');[S;T]
fprintf('搜索点序号： ');t
fprintf('min f(x)： ');vpa(fmin,6)
if normest(gr)
fprintf('找到最优解');
break;
end %若满足终止条件，结束循环
dm=step(d0,xx);
fprintf('搜索步长λ： ');vpa(dm,6)
x=subs(xx,d,dm);
x0=x;
[yf,zf,C,D,gr]=gx(x0);
[d0,xx,S,T]=initd(C,D,gr,x0);
end
%绘图表示搜索过程
figure(1)
plot(xkt);
xlabel('搜索点序号');
ylabel('搜索点');
figure(2)
plot(fxt);
xlabel('搜索点序号');
ylabel('目标函数');
运行结果如下：
初始计算X0结果：
X0=
1.0e+03 *
2.[**************]
7.[**************]
7.[**************]
0.[**************]
0.[**************]
0.[**************]
0.[**************]
0.[**************]
最终计算结果方向：
-1.[**************]
0.[**************]
-1.[**************]
-1.[**************]
-1.[**************]
-1.[**************]
1.[**************]
-1.[**************]
最终计算结果步长：
1.[**************]e-08
最终计算F结果：
F=
1.[**************]e+04
最终计算X0结果：
X0=
1.0e+03 *
0.[**************]
6.[**************]
6.[**************]
0.[**************]
0.[**************]
0.[**************]
0.[**************]
0.[**************]
最终计算结果约束集：
5
在当前所选的初始点（随机生成的初始点），不能取得一个较好的结果，取迭代200次后
的值作为最后的结果。
X*=[100,1276.3,7948.4,115.5,225.4,261.5,238.8,316.4]' ，f*=9324.64 。
2.3结果比较与分析
Zoutendijk’s Method Generalized Reduced Gradient Method
Zoutendijk’s Method Generalized Reduced Gradient Method
通过上图比较可知，Zoutendijk’s Method计算过程中，目标函数值呈阶梯形下降，通过增加循环上限次数，最终结果将越来越逼近最优解；而Generalized Reduced Gradient Method计算过程中，目标函数值在经过几次迭代之后基本就保持不变了。
总而言之，Zoutendijk’s Method可能出现锯齿现象，当某一个迭代点在约束边界上时，如果可行方向选取不当，会使迭代点跑出边界，甚至使得算法不收敛。Rosen’s
Gradient Projection Method求解带线性约束的非线性规划问题十分有效。从一个基本可行解开始，由约束条件确定出凸约束集边界上梯度的投影，以便求出下次的搜索方向和步长。每次搜索后，都要进行检验，直到满足精度要求为止。Generalized Reduced Gradient Method收敛快，精度高，被认为求解带线性约束的非线性规划问题时有显著效果。
三、致谢
通过工程优化设计这门课，自己学习了关于最优化设计的一些算法与基础知识，也让
我了解到自己在设计上的不足，同时黄老师的授课与严谨的教学，也让我懂得了认真严谨
是作为一个科研工作者最基本的素质和态度，所以在此对黄老师表示真挚的感谢。
我也会在今后的学习中，脚踏实地，严谨治学，积极向老师和同学们请教与讨论，争取早日取得优秀的成果。