最优化上机实习报告

1、 用黄金分割法编程求解最优化问题：

minf(x) 2x2 x 1s..t 1 x 1

精度 0.16

①算法（黄金分割法）： 给定a,b(a b)及 0。

Step1 令x2 a 0.618(b a),f2 f(x2)转step2. Step 2 令x1 a 0.382(b a),f1 f(x1)转step3. Step3 若b a ，则x*

a b

，停。否则转step4. 2

Step4 若f1 f2，则b x2,x2 x1,f2 f1，转step2； 若f1 f2，则a x1,b x2转step1； 若f1 f2，则a x1,x1 x2,f1 f2，转step5. Step5 令x2 a 0.618(b a),f2 f(x2)转step3. ②程序：

1）meng.m文件（原函数程序）：

function f=meng(x) f=2*x^2-x-1; 2）grad.m文件（主程序）：

global s t x f1 f2; a=-1;b=1;c=0.16; s=a+0.382*(b-a); t=a+0.618*(b-a); f1=meng(s); f2=meng(t);

while abs(b-a)>c if f1<f2

b=t;t=s;f2=f1;s=a+0.382*(b-a);f1=meng(s); elseif f1==f2

a=s;b=t;t=a+0.618*(b-a);f2=meng(t); else

a=s;s=t;f1=f2;t=a+0.618*(b-a);f2=meng(t); end

end

x=(a+b)/2

f=meng(x)

③Matlab运行结果：x=0.229,f=-1.1235.

④分析：用黄金分割法所得结果为x=0.229,f=-1.1235，与真实结果x=0.25，f=-1.125

比较接近；该方法是线性收敛的，收敛比为

1

;并且精度对精确2

度有一定的影响，但并不是精度越高精确度越好，太高的精度会增加计算量，但对结果精度的提高没有影响或影响甚微。 2、编程求解无约束最优化问题

minz (x12 x2 11)2 (x1 x22 7)2

①算法（牛顿法）： 给定控制误差 0。 Step1 取初始点x0，令k 0； Step2 计算gk g(xk)；

Step3

若gk ，则x* xk，停；否则计算Gk，并由Gkpk gk解出pk； Step4 令xk 1 xk pk,k k 1，转Step2. ②程序：

1）fun.m文件（原函数程序）：

function z=fun(x)

z=(x(1)^2+x(2)-11)^2+(x(1)+x(2)^2-7)^2; 2）tfun.m文件（求梯度函数程序）：

function t=tfun(x)

t=[4*x(1)*(x(1)^2+x(2)-11)+2*(x(1)+x(2)^2-7),2*(x(1)^2+x(2)-11)+4*x(2)*(x(1)+x(2)^2-7)]'; 3）ttfun.m文件（求原函数Hesse矩阵程序）：

function tt=ttfun(x)

tt=[12*x(1)^2+4*x(2)-42,4*x(1)+4*x(2);4*x(1)+4*x(2),4*x(1)+12*x(2)^2-26]; 4）grad.m文件（主程序）：

global p; e=1e-3;k=1; t=tfun(x0);

while norm(t)>e tt=ttfun(x0); p=-inv(tt)*t; x0=x0+p'; k=k+1;

t=tfun(x0);

end k x=x0

z=fun(x0)

③Matlab运行结果：在命令窗口中输入初始点，运行程序。

取初始点为x0=[1 1],则k=8,x=[-3.7793 -3.2832],z=3.0880e-017, 取初始点为x0=[3.1 1.9],则k=3,x=[3.0000 2.0000],z=7.1523e-009 ④分析：当初始点为x0 [11]时，x [ 3.7793 3.2832],z 3.0880e 017，

当初始点为x0 [3.11.9]时，x [3.00002.0000],z 7.1523e-009,而真实结果为x [32],z 0。这表明牛顿法有很快的收敛速度，但它只是局部收敛的，即只有当初始点x0充分接近真实值x*时，才能保证收敛。但若x0离x*较远时，则不能保证迭代点列{xk}收敛，甚至不能保证{f(xk)}单调下降。因此选择初始点时应充分接近真实值。

3、编程求解无约束最优化问题：

minz (x1x2)2(1 x1)2[1 x1 x2(1 x1)5]2

取不同的初始值求解，尽可能求出所有局部极小点，进而找出全局极小点。 ①算法（最速下降法）： 给定控制误差 0。 Step1 取初始点x0，令k 0； Step2 计算gk g(xk)；

Step3 若gk ，则x* xk，停；否则，令pk gk，由一维搜索求步长 k，

使得

f(xk kpk) minf(xk pk)；

0

Step4 令xk 1 xk kpk,k k 1，转Step2.

②程序：

1） meng.m文件（原函数程序）：

function z=meng(x)

z=(x(1)*x(2))^2*(1-x(1))^2*(1-x(1)-x(2)*(1-x(1))^5)^2; 2） tmeng.m文件（求梯度函数程序）：

function t=tmeng(x)

t=[2*x(1)*x(2)*(1-x(1))^2*(1-x(1)-x(2)*(1-x(1))^5)^2-2*(x(1)*x(2))^2*(1-x(1))*(1-x(1)-x(2)*(1-x(1))^5)^2+2*(x(1)*x(2))^2*(1-x(1))^2*(1-x(1)-x(2)*(1-x(1))^5)*(5*x(2)*(1-x(1))^4-1),2*x(1)*x(2)*(1-x(1))^2*(1-x(1)-x(2)*(1-x(1))^5)^2-2*(x(1)*x(2))^2*(1-x(1))^7*(1-x(1)-x(2)*(1-x(1))^5)]';

3） grad.m文件（主程序）：

k=0;a=0.5;b=0.6;c=0.01;maxk=6000; while k<maxk

t=feval('tmeng',x0); d=-t;

if norm(d)<=c,break; end

m=0;mk=0; while m<30 if

(feval('meng',x0+a^m*d')<feval('meng',x0)+b*a^m*t'*d) mk=m;break; end m=m+1; end

x0=x0+a^mk*d'; k=k+1; end k x=x0

z=feval('meng',x)

③Matlab运行结果：在命令窗口中输入不同的初始值，运行程序。 当x0=[3 3]时，k=28,x=[1.7687 2.8641],z= 6.5031e-008; 当x0=[3 2]时，k=25,x=[1.8629 1.8041],z=6.6470e-008; 当x0=[3 1]时，k=20,x=[2.1173 0.6419],z=1.0022e-007;

当x0=[3 1.5]时，k=21,x=[1.9465 1.2465],z=2.0851e-007; 当x0=[3 0]时，k=0,x=[3 0],z=0;

当x0=[2 3]时，k=15,x=[1.7610 2.9821],z=8.8086e-008; 当x0=[2 0]时，k=0,x=[2 0],z=0; 当x0=[1 3]时，k=0,x=[1 3],z=0; 当x0=[1 2]时，k=0,x=[1 2],z=0; 当x0=[0 2]时，k=0,x=[0 2],z=0; 当x0=[0 1]时，k=0,x=[0 1],z=0.

故局部极小点为:x=[2 3],x=[2 2],x=[2 1],x=[3 0],x=[2 0], x=[0 1],x=[0 2],x=[1 3],x=[1 2]；z=0.全局极小点中z=0. ④分析：用最速下降法取不同初始值，结果也不相同。但题中当x1或x2有一个为0时，局部极小点为x=[0 *]或x=[* 0],此时z=0;x1为1时，局部极小点为x=[1 *],此时z=0。都有z=0。最速下降法具有很好的整体收敛性，且收敛速度是线性的，但收敛度慢，是基本算法之一，但不是有效的实用算法。 4、用外罚函数法编程求解约束优化问题：

minf(x) x1x2x3s..72t x1 2x2 2x3 0

①算法（外罚函数法）：

取控制误差 0和罚因子的放大系数c 1(可取 =10 4，c 10)。

Step1 给定初始点x0（可以不是可行点）和初试罚因子 1（可取 1 1），令k 1。 Step2 以xk 1为初始点求无约束问题：

(x) minp(x, k) f(x) kp

(x) ci(x) 其中p

i 1l

j l 1

m

min(0,cj(x)), 1, 1，得最优解xk x( k)。

(x) ，则以xk为近似最优解，停止。否则令 k 1 c k,k k 1，Step3 若 kp

转Step2.

②程序：

1） meng.m文件（原函数程序）： function P=hmeng(x) global a;

P=meng(x)+a*cmeng(x); 2） cmeng.m文件（求惩罚函数程序）： function p=cmeng(x)

p=(abs(72-x(1)-2*x(2)-2*x(3)))^3; 3） hmeng.m文件（无约束问题求解函数）： function P=hmeng(x) global a;

P=meng(x)+a*cmeng(x); 4） grad.m文件（主函数）：

global a;

c=1e-5; %c为控制误差

a=1;b=10; %a为初始罚因子，b为罚因子的放大系数 x0=[5 5 5]; k=1;

while a*cmeng(x0)>=c

x0=fminsearch('hmeng',x0); …… 此处隐藏：2604字，全部文档内容请下载后查看。喜欢就下载吧 ……