大连理工大学优化方法上机大作业.docx
- 文档编号:25272684
- 上传时间:2023-06-06
- 格式:DOCX
- 页数:27
- 大小:4.60MB
大连理工大学优化方法上机大作业.docx
《大连理工大学优化方法上机大作业.docx》由会员分享,可在线阅读,更多相关《大连理工大学优化方法上机大作业.docx(27页珍藏版)》请在冰豆网上搜索。
大连理工大学优化方法上机大作业
2016年理工大学优化方法上机大作业
学院:
专业:
班级:
学号:
:
上机大作业1:
1.最速下降法:
functionf=fun(x)
f=(1-x
(1))^2+100*(x
(2)-x
(1)^2)^2;
end
functiong=grad(x)
g=zeros(2,1);
g
(1)=2*(x
(1)-1)+400*x
(1)*(x
(1)^2-x
(2));
g
(2)=200*(x
(2)-x
(1)^2);
end
functionx_star=steepest(x0,eps)
gk=grad(x0);
res=norm(gk);
k=0;
whileres>eps&&k<=1000
dk=-gk;
ak=1;f0=fun(x0);
f1=fun(x0+ak*dk);
slope=dot(gk,dk);
whilef1>f0+0.1*ak*slope
ak=ak/4;
xk=x0+ak*dk;
f1=fun(xk);
end
k=k+1;
x0=xk;
gk=grad(xk);
res=norm(gk);
fprintf('--The%d-thiter,theresidualis%f\n',k,res);
end
x_star=xk;
end
>>clear
>>x0=[0,0]';
>>eps=1e-4;
>>x=steepest(x0,eps)
2.牛顿法:
functionf=fun(x)
f=(1-x
(1))^2+100*(x
(2)-x
(1)^2)^2;
end
functiong=grad2(x)
g=zeros(2,2);
g(1,1)=2+400*(3*x
(1)^2-x
(2));
g(1,2)=-400*x
(1);
g(2,1)=-400*x
(1);
g(2,2)=200;
end
functiong=grad(x)
g=zeros(2,1);
g
(1)=2*(x
(1)-1)+400*x
(1)*(x
(1)^2-x
(2));
g
(2)=200*(x
(2)-x
(1)^2);
end
functionx_star=newton(x0,eps)
gk=grad(x0);
bk=[grad2(x0)]^(-1);
res=norm(gk);
k=0;
whileres>eps&&k<=1000
dk=-bk*gk;
xk=x0+dk;
k=k+1;
x0=xk;
gk=grad(xk);
bk=[grad2(xk)]^(-1);
res=norm(gk);
fprintf('--The%d-thiter,theresidualis%f\n',k,res);
end
x_star=xk;
end
>>clear
>>x0=[0,0]';
>>eps=1e-4;
>>x1=newton(x0,eps)
--The1-thiter,theresidualis447.213595
--The2-thiter,theresidualis0.000000
x1=
1.0000
1.0000
3.BFGS法:
functionf=fun(x)
f=(1-x
(1))^2+100*(x
(2)-x
(1)^2)^2;
end
functiong=grad(x)
g=zeros(2,1);
g
(1)=2*(x
(1)-1)+400*x
(1)*(x
(1)^2-x
(2));
g
(2)=200*(x
(2)-x
(1)^2);
end
functionx_star=bfgs(x0,eps)
g0=grad(x0);
gk=g0;
res=norm(gk);
Hk=eye
(2);
k=0;
whileres>eps&&k<=1000
dk=-Hk*gk;
ak=1;f0=fun(x0);
f1=fun(x0+ak*dk);
slope=dot(gk,dk);
whilef1>f0+0.1*ak*slope
ak=ak/4;
xk=x0+ak*dk;
f1=fun(xk);
end
k=k+1;
fa0=xk-x0;
x0=xk;
go=gk;
gk=grad(xk);
y0=gk-g0;
Hk=((eye
(2)-fa0*(y0)')/((fa0)'*(y0)))*((eye
(2)-(y0)*(fa0)')/((fa0)'*(y0)))+(fa0*(fa0)')/((fa0)'*(y0));
res=norm(gk);
fprintf('--The%d-thiter,theresidualis%f\n',k,res);
end
x_star=xk;
End
>>clear
>>x0=[0,0]';
>>eps=1e-4;
>>x=bfgs(x0,eps)
4.共轭梯度法:
functionf=fun(x)
f=(1-x
(1))^2+100*(x
(2)-x
(1)^2)^2;
end
functiong=grad(x)
g=zeros(2,1);
g
(1)=2*(x
(1)-1)+400*x
(1)*(x
(1)^2-x
(2));
g
(2)=200*(x
(2)-x
(1)^2);
end
functionx_star=CG(x0,eps)
gk=grad(x0);
res=norm(gk);
k=0;
dk=-gk;
whileres>eps&&k<=1000
ak=1;f0=fun(x0);
f1=fun(x0+ak*dk);
slope=dot(gk,dk);
whilef1>f0+0.1*ak*slope
ak=ak/4;
xk=x0+ak*dk;
f1=fun(xk);
end
k=k+1;
x0=xk;
g0=gk;
gk=grad(xk);
res=norm(gk);
p=(gk/g0)^2;
dk1=dk;
dk=-gk+p*dk1;
fprintf('--The%d-thiter,theresidualis%f\n',k,res);
end
x_star=xk;
end
>>clear
>>x0=[0,0]';
>>eps=1e-4;
>>x=CG(x0,eps)
上机大作业2:
functionf=obj(x)
f=4*x
(1)-x
(2)^2-12;
end
function[h,g]=constrains(x)
h=x
(1)^2+x
(2)^2-25;
g=zeros(3,1);
g
(1)=-10*x
(1)+x
(1)^2-10*x
(2)+x
(2)^2+34;
g
(2)=-x
(1);
g(3)=-x
(2);
end
functionf=alobj(x)%拉格朗日增广函数
%N_equ等式约束个数?
%N_inequ不等式约束个数
N_equ=1;
N_inequ=3;
globalr_alpena;%全局变量
h_equ=0;
h_inequ=0;
[h,g]=constrains(x);
%等式约束部分?
fori=1:
N_equ
h_equ=h_equ+h(i)*r_al(i)+(pena/2)*h(i).^2;
end
%不等式约束部分
fori=1:
N_inequ
h_inequ=h_inequ+(0.5/pena)*(max(0,(r_al(i)+pena*g(i))).^2-r_al(i).^2);
end
%拉格朗日增广函数值
f=obj(x)+h_equ+h_inequ;
functionf=compare(x)
globalr_alpenaN_equN_inequ;
N_equ=1;
N_inequ=3;
h_inequ=zeros(3,1);
[h,g]=constrains(x);
%等式部分
fori=1:
1
h_equ=abs(h(i));
end
%不等式部分
fori=1:
3
h_inequ=abs(max(g(i),-r_al(i+1)/pena));
end
h1=max(h_inequ);
f=max(abs(h_equ),h1);%sqrt(h_equ+h_inequ);
function[x,fmin,k]=almain(x_al)
%本程序为拉格朗日乘子算法示例算法%函数输入:
%x_al:
初始迭代点
%r_al:
初始拉格朗日乘子N-equ:
等式约束个数N_inequ:
不等式约束个数?
%函数输出
%X:
最优函数点FVAL:
最优函数值
%============================程序开始================================
globalr_alpena;%参数(全局变量)
pena=10;%惩罚系数
r_al=[1,1,1,1];
c_scale=2;%乘法系数乘数
cta=0.5;%下降标准系数
e_al=1e-4;%误差控制围
max_itera=25;
out_itera=1;%迭代次数
%===========================算法迭代开始=============================
whileout_itera x_al0=x_al; r_al0=r_al; %判断函数? compareFlag=compare(x_al0); %无约束的拟牛顿法BFGS [X,fmin]=fminunc(alobj,x_al0); x_al=X;%得到新迭代点 %判断停止条件? ifcompare(x_al) disp('wegettheoptpoint'); break end %c判断函数下降度? ifcompare(x_al) pena=1*pena;%可以根据需要修改惩罚系数变量 else pena=min(1000,c_scale*pena);%%乘法系数最大1000 disp('pena=2*pena'); end %%? 更新拉格朗日乘子 [h,g]=constrains(x_al); fori=1: 1 %%等式约束部分 r_al(i)=r_al0(i)+pena*h(i); end fori=1: 3 %%不等式约束部分 r_al(i+1)=max(0,(r_al0(i+1)+pena*g(i))); end out_itera=out_itera+1; end %+++++++++++++++++++++++++++迭代结束+++++++++++++++++++++++++++++++++ disp('theiterationnumber'); k=out_itera; disp('thevalueofconstrains'); compare(x_al) disp('theoptpoint'); x=x_al; fmin=obj(X); >>clear >>x_al=[0,0]; >>[x,fmin,k]=almain(x_al) 上机大作业3: 1、 >>clearall n=3;c=[-3,-1,-3]';A=[2,1,1;1,2,3;2,2,1;-1,0,0;0,-1,0;0,0,-1];b=[2,5,6,0,0,0]'; cvx_begin variablex(n) minimize(c'*x) subjectto A*x<=b cvx_end CallingSDPT34.0: 6variables,3equalityconstraints ------------------------------------------------------------ num.ofconstraints=3 dim.oflinearvar=6 ******************************************************************* SDPT3: Infeasiblepath-followingalgorithms ******************************************************************* versionpredcorrgamexponscale_data NT10.00010 itpstepdsteppinfeasdinfeasgapprim-objdual-objcputime ------------------------------------------------------------------- 0|0.000|0.000|1.1e+01|5.1e+00|6.0e+02|-7.000000e+010.000000e+00|0: 0: 00|chol11 1|0.912|1.000|9.4e-01|4.6e-02|6.5e+01|-5.606627e+00-2.967567e+01|0: 0: 01|chol11 2|1.000|1.000|1.3e-07|4.6e-03|8.5e+00|-2.723981e+00-1.113509e+01|0: 0: 01|chol11 3|1.000|0.961|2.3e-08|6.2e-04|1.8e+00|-4.348354e+00-6.122853e+00|0: 0: 01|chol11 4|0.881|1.000|2.2e-08|4.6e-05|3.7e-01|-5.255152e+00-5.622375e+00|0: 0: 01|chol11 5|0.995|0.962|1.6e-09|6.2e-06|1.5e-02|-5.394782e+00-5.409213e+00|0: 0: 01|chol11 6|0.989|0.989|2.7e-10|5.2e-07|1.7e-04|-5.399940e+00-5.400100e+00|0: 0: 01|chol11 7|0.989|0.989|5.3e-11|5.8e-09|1.8e-06|-5.399999e+00-5.400001e+00|0: 0: 01|chol11 8|1.000|0.994|2.8e-13|4.3e-11|2.7e-08|-5.400000e+00-5.400000e+00|0: 0: 01| stop: max(relativegap,infeasibilities)<1.49e-08 ------------------------------------------------------------------- numberofiterations=8 primalobjectivevalue=-5.39999999e+00 dualobjectivevalue=-5.40000002e+00 gap: =trace(XZ)=2.66e-08 relativegap=2.26e-09 actualrelativegap=2.21e-09 rel.primalinfeas(scaledproblem)=2.77e-13 rel.dual"""=4.31e-11 rel.primalinfeas(unscaledproblem)=0.00e+00 rel.dual"""=0.00e+00 norm(X),norm(y),norm(Z)=4.3e+00,1.3e+00,1.9e+00 norm(A),norm(b),norm(C)=6.7e+00,9.1e+00,5.4e+00 TotalCPUtime(secs)=0.71 CPUtimeperiteration=0.09 terminationcode=0 DIMACS: 3.6e-130.0e+005.8e-110.0e+002.2e-092.3e-09 ------------------------------------------------------------------- ------------------------------------------------------------ Status: Solved Optimalvalue(cvx_optval): -5.4 2、 >>clearall n=2;c=[-2,-4]';G=[0.5,0;0,1]; A=[1,1;-1,0;0,-1];b=[1,0,0]'; cvx_begin variablex(n) minimize(x'*G*x+c'*x) subjectto A*x<=b cvx_end CallingSDPT34.0: 7variables,3equalityconstraints Forimprovedefficiency,SDPT3issolvingthedualproblem. ------------------------------------------------------------ num.ofconstraints=3 dim.ofsocpvar=4,num.ofsocpblk=1 dim.oflinearvar=3 ******************************************************************* SDPT3: Infeasiblepath-followingalgorithms ******************************************************************* versionpredcorrgamexponscale_data NT10.00010 itpstepdsteppinfeasdinfeasgapprim-objdual-objcputime ------------------------------------------------------------------- 0|0.000|0.000|8.0e-01|6.5e+00|3.1e+02|1.000000e+010.000000e+00|0: 0: 00|chol11 1|1.000|0.987|4.3e-07|1.5e-01|1.6e+01|9.043148e+00-2.714056e-01|0: 0: 00|chol11 2|1.000|1.000|2.6e-07|7.6e-03|1.4e+00|1.234938e+00-5.011630e-02|0: 0: 00|chol11 3|1.000|1.000|2.4e-07|7.6e-04|3.0e-01|4.166959e-011.181563e-01|0: 0: 00|chol11 4|0.892|0.877|6.4e-08|1.6e-04|5.2e-02|2.773022e-012.265122e-01|0: 0: 00|chol11 5|1.000|1.000|1.0e-08|7.6e-06|1.5e-02|2.579468e-012.427203e-01|0: 0: 00|chol11 6|0.905|0.904|3.1e-09|1.4e-06|2.3e-03|2.511936e-012.488619e-01|0: 0: 00|chol11 7|1.000|1.000|6.1e-09|7.7e-08|6.6e-04|2.503336e-012.496718e-01|0: 0: 00|chol11 8|0.903|0.903|1.8e-09|1.5e-08|1.0e-04|2.500507e-012.499497e-01|0: 0: 00|chol11 9|1.000|1.000|4.9e-10|3.5e-10|2.9e-05|2.500143e-012.499857e-01|0: 0: 00|chol11 10|0.904|0.904|4.7e-11|1.3e-10|4.4e-06|2.500022e-012.499978e-01|0: 0: 00|chol22 11|1.000|1.000|2.3e-12|9.4e-12|1.2e-06|2.500006e-012.499994e-01|0: 0: 00|chol22 12|1.000|1.000|4.7e-13|1.0e-12|1.8e-07|2.500001e-012.499999e-01|0: 0: 00|chol22 13|1.000|1.000|2.0e-12|1.0e-12|4.2e-08|2.500000e-012.500000e-01|0: 0: 00|chol22 14|1.000|1.000|2.6e-12|1.0e-12|7.3e-09|2.500000e-012.500000e-01|0: 0: 00| stop: max(relativegap,infeasibilities)<1.49e-08 ----------------------------
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 大连理工大学 优化 方法 上机 作业
![提示](https://static.bdocx.com/images/bang_tan.gif)