最优化 多目标优化 惩罚函数法 梯度法 牛顿法

2008-12-08 12:30

利用梯度法和牛顿法编程求最优解(matlab)

f(x)=x1^2+4*x2^2 x0=[2;2] e=0.002

利用梯度法和牛顿法编程求最优解

方法一.梯度法

function y=fun(x1,x2)

y=x1^2+4*x2^2; %定义fun.m函数

clc

syms x1 x2 d;

f=x1^2+4*x2^2;

fx1=diff(f,'x1');

fx2=diff(f,'x2');

x1=2;

x2=2;

for n=1:100

f0=subs(f);

f1=subs(fx1);

f2=subs(fx2);

if (double(sqrt(f1^2+f2^2)) <= 0.002)

n

vpa(x1)

vpa(x2)

vpa(f0)

break;

else

D=fun(x1-d*f1,x2-d*f2);

Dd=diff(D,'d');

dd=solve(Dd);

x1=x1-dd*f1;

x2=x2-dd*f2;

end

end %结果n=10,x1=0.2223e-3,x2=-0.1390e-4,f0=0.5021e-7. 方法二.牛顿法

clc

syms x1 x2 ;

f=x1^2+4*x2^2;

fx1=diff(f,'x1'); fx2=diff(f,'x2');

fx1x1=diff(fx1,'x1');fx1x2=diff(fx1,'x2');fx2x1=diff(fx2,'x1');fx2x2= diff(fx2,'x2');

x1=2;

x2=2;

for n=1:100

f0=subs(f);

f1=subs(fx1);

f2=subs(fx2);

if (double(sqrt(f1^2+f2^2)) <= 0.002)

n

x1=vpa(x1,4)

x2=vpa(x2,4)

f0=vpa(f0,4)

break;

else

X=[x1 x2]'-inv([fx1x1 fx1x2;fx2x1 fx2x2]) *[f1 f2]';

x1=X[1,1];

x2=X[2,1];

end

end %结果 n=2,x1=0,x2=0,f0=0.

惩罚函数法（内点法、外点法）求解约束优化问题最优值编程 matlab

1 用外点法求下列问题的最优解

方法一：外点牛顿法：

clc

m=zeros(1,50);a=zeros(1,50);b=zeros(1,50);f0=zeros(1,50);%a b为最优点坐标，f0为最优点函数值，f1 f2最优点梯度。

syms x1 x2 e; %e为罚因子。m(1)=1;c=10;a(1)=0;b(1)=0; %c为递增系数。赋初值。

f=x1^2+x2^2+e*(1-x1)^2;f0(1)=1;

fx1=diff(f,'x1');fx2=diff(f,'x2');fx1x1=diff(fx1,'x1');fx1x2=diff(fx1 ,'x2');fx2x1=diff(fx2,'x1');fx2x2=diff(fx2,'x2');%求偏导、海森元素。for

k=1:100 %外点法e迭代循环.

x1=a(k);x2=b(k);e=m(k);

for

n=1:100 %梯度法求最优值。

f1=subs(fx1); %求解梯度值和海森矩阵

f2=subs(fx2);

f11=subs(fx1x1);

f12=subs(fx1x2);

f21=subs(fx2x1);

f22=subs(fx2x2);

if(double(sqrt(f1^2+f2^2))<=0.001) %最优值收敛条件

a(k+1)=double(x1);b(k+1)=double(x2);f0(k+1)=double(subs(f ));

break;

else

X=[x1 x2]'-inv([f11 f12;f21 f22])*[f1 f2]';

x1=X(1,1);x2=X(2,1);

end

end

if(double(sqrt((a(k+1)-a(k))^2+(b(k+1)-b(k))^2))<=0.001)&&(double(abs ((f0(k+1)-f0(k))/f0(k)))<=0.001) %罚因子迭代收敛条件

a(k+1) %输出最优点坐标，罚因子迭代次数，最优值

b(k+1)

k

f0(k+1)

break;

else

m(k+1)=c*m(k);

end

end

方法二：外点梯度法：

clc

m=zeros(1,50);a=zeros(1,50);b=zeros(1,50);f0=zeros(1,50); syms d x1 x2 e;

m(1)=1;c=10;a(1)=0;b(1)=0; f=x1^2+x2^2+e*(1-x1)^2; f0(1)=1;

fx1=diff(f,'x1');

fx2=diff(f,'x2');

for k=1:100

x1=a(k);x2=b(k);e=m(k);

for n=1:100

f1=subs(fx1);

f2=subs(fx2);

if(double(sqrt(f1^2+f2^2))<=0.002)