共轭梯度法c++程序10页.docx
- 文档编号:5298321
- 上传时间:2023-05-08
- 格式:DOCX
- 页数:7
- 大小:16.68KB
共轭梯度法c++程序10页.docx
《共轭梯度法c++程序10页.docx》由会员分享,可在线阅读,更多相关《共轭梯度法c++程序10页.docx(7页珍藏版)》请在冰点文库上搜索。
共轭梯度法c++程序10页
最优化课程设计
题目:
共轭梯度法
姓名:
田鑫
指导老师:
智红英
学号:
201118030216
班级:
信息与计算科学111802班
共轭梯度法(ConjugateGradient)
是介于最速下降法与牛顿法之间的一个方法,它仅需利用一阶导数信息,但克服了最速下降法收敛慢的缺点,又避免了牛顿法需要存储和计算Hesse矩阵并求逆的缺点,共轭梯度法不仅是解决大型线性方程组最有用的方法之一,也是解大型非线性最优化最有效的算法之一。
在各种优化算法中,共轭梯度法是非常重要的一种。
其优点是所需存储量小,具有步收敛性,稳定性高,而且不需要任何外来参数。
设我们要求解下列线性系统
其中n-×-n矩阵A是对称的(也即,AT=A),正定的(也即,xTAx>0对于所有非0向量x属于Rn),并且是实系数的。
将系统的唯一解记作x*。
最后算法
经过一些简化,可以得到下列求解Ax=b的算法,其中A是实对称正定矩阵。
x0:
=0
k:
=0
r0:
=b
repeatuntilrkis"sufficientlysmall":
k:
=k+1
ifk=1
p1:
=r0
else
endif
xk:
=xk-1+αkpk
rk:
=rk-1-αkApk
endrepeat
结果为xk
共轭梯度法程序源代码
#include
#include
#defineN10
#defineepspow(10,-6)
doublef(doublex[],doublep[],doublet)
{
doubles;
s=pow(x[0]+t*p[0],2)+25*pow(x[1]+t*p[1],2);
returns;
}
/*以下是进退法搜索区间源程序*/
voidsb(double*a,double*b,doublex[],doublep[])
{
doublet0,t1,t,h,alpha,f0,f1;
intk=0;
t0=2.5; /*初始值*/
h=1; /*初始步长*/
alpha=2; /*加步系数*/
f0=f(x,p,t0);
t1=t0+h;
f1=f(x,p,t1);
while
(1)
{
if(f1 { h=alpha*h;t=t0; t0=t1; f0=f1; k++; } else { if(k==0) {h=-h;t=t1;} else { *a=t t: t1; *b=t>t1? t: t1; break; } } t1=t0+h; f1=f(x,p,t1); } }/*以下是黄金分割法程序源代码*/ doublehjfg(doublex[],doublep[]) { doublebeta,t1,t2,t; doublef1,f2; doublea=0,b=0; double*c,*d; c=&a,d=&b; sb(c,d,x,p);/*调用进退法搜索区间*/ printf("\nx1=%f,x2=%f,p1=%f,p2=%f",x[0],x[1],p[0],p[1]); printf("\n[a,b]=[%f,%f]",a,b); beta=(sqrt(5)-1.0)/2; t2=a+beta*(b-a);f2=f(x,p,t2); t1=a+b-t2; f1=f(x,p,t1); while (1) { if(fabs(t1-t2)<0) break; else { if(f1 { t=(t1+t2)/2; b=t2; t2=t1; f2=f1; t1=a+b-t2; f1=f(x,p,t1); } else { a=t1; t1=t2; f1=f2; t2=a+beta*(b-a); f2=f(x,p,t2); } } } t=(t1+t2)/2; returnt; } /*以下是共轭梯度法程序源代码*/ voidgtd() { doublex[N],g[N],p[N],t=0,f0,mod1=0,mod2=0,nanda=0; inti,k,n; printf("请输入函数的元数值n=2"); scanf("%d",&n); printf("\n请输入初始值: 2,2"); for(i=0;i scanf("%lf",&x[i]); f0=f(x,g,t); g[0]=2*x[0];g[1]=50*x[1]; mod1=sqrt(pow(g[0],2)+pow(g[1],2));/*求梯度的长度*/ if(mod1>0) { p[0]=-g[0];p[1]=-g[1];k=0; while (1) { t=hjfg(x,p);/*调用黄金分割法求t的值*/ printf("\np1=%f,p2=%f,t=%f",p[0],p[1],t); x[0]=x[0]+t*p[0];x[1]=x[1]+t*p[1]; g[0]=2*x[0]; g[1]=50*x[1]; /*printf("\nx1=%f,x2=%f,g1=%f,g2=%f",x[0],x[1],g[0],g[1]);*/ mod2=sqrt(pow(g[0],2)+pow(g[1],2));/*求梯度的长度*/ if(mod2<=0) break; else { if(k+1==n) { g[0]=2*x[0];g[1]=50*x[1]; p[0]=-g[0];p[1]=-g[1];k=0; } else { nanda=pow(mod2,2)/pow(mod1,2); printf("\nnanda=%f,mod=%f",nanda,mod2); p[0]=-g[0]+nanda*p[0]; p[1]=-g[1]+nanda*p[1]; mod1=mod2; k++; } } printf("\n--------------------------"); } } printf("\n最优解为x1=%f,x2=%f",x[0],x[1]); printf("\n最终的函数值为%f",f(x,g,t)); } main() { gtd(); } 运行结果如下: 请输入函数的元数值n=2 请输入初始值: 22 x1=2.000000,x2=2.000000,p1=-4.000000,p2=-100.000000 [a,b]=[-4.500000,1.500000] p1=-4.000000,p2=-100.000000,t=0.020030 nanda=0.001474,mod=3.842730 -------------------------- x1=1.919879,x2=-0.003022,p1=-3.845655,p2=0.003665 [a,b]=[-4.500000,1.500000] p1=-3.845655,p2=0.003665,t=0.499240 -------------------------- x1=-0.000026,x2=-0.001192,p1=0.000052,p2=0.059610 [a,b]=[-4.500000,1.500000] p1=0.000052,p2=0.059610,t=0.020000 nanda=0.000000,mod=0.000050 -------------------------- x1=-0.000025,x2=-0.000000,p1=0.000050,p2=0.000001 [a,b]=[-4.500000,1.500000] p1=0.000050,p2=0.000001,t=0.495505 -------------------------- x1=-0.000000,x2=0.000000,p1=0.000000,p2=-0.000023 [a,b]=[-4.500000,1.500000] p1=0.000000,p2=-0.000023,t=0.020007 最优解为x1=-0.000000,x2=-0.000000 最终的函数值为0.000000 希望以上资料对你有所帮助,附励志名言3条: 1、生气,就是拿别人的过错来惩罚自己。 原谅别人,就是善待自己。 2、未必钱多乐便多,财多累己招烦恼。 清贫乐道真自在,无牵无挂乐逍遥。 3、处事不必求功,无过便是功。 为人不必感德,无怨便是德。
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 共轭 梯度 c+ 程序 10