阻尼牛顿法

阻尼牛顿法求解二维函数极小值的程序说明一题目利用阻尼牛顿法求函数2212112(,)22fxxxxx（）（）的极小值点（迭代两次，一维搜索任选一种方法）。二阻尼牛顿法基本思想:1)给定初始点0x，收敛精度，置0k。2)计算()kfx、2()kfx、21(())kfx和21(())()kkkdfxfx3)求1kkkkxxd，其中k为沿kd进行一维搜索的最佳步长。4)检查收敛精度。若1kknxx，则*1kxx，停机；否则置1kk，返回步骤2，继续进行进行搜索。改进后的阻尼牛顿法程序框图如下：开始0,x给定0k1?kknxx*1kxx1kk是否结束12()()kkkfxfxd1:()kkkkkkkdxxminfdx三用阻尼牛顿法求函数程序如下：//阻尼牛顿法.cpp:Definestheentrypointfortheconsoleapplication.//#includestdio.h#includemath.h#includeconio.h#includeiostreamdoublefun1(doubleq1,doubleq2){return(pow((q1-2),4)+pow((q1-2*q2),2));//修改函数f(x1,x2)=(x1-2)*(x1-2)*(x1-2)*(x1-2)+(x1-2*x2)*(x1-2*x2)}doublefun2(doubleg,doublex,doubley,doubler1,doubler2){return(pow((x+g*y-2),4)+pow((x+g*y-2*(r1+g*r2)),2));//关于阻尼因子的函数}voidmain(){doubleA[2][1],B[2][2],C[2][1],D[2][1],X[2][1];doubleE[2][1]={4,3};//迭代的初始点x0intt=0,i=0,j=0;doubleE0,x1,x2,x3,h(0.1);doubley1,y2,y3,m;doublea,b,k=0.618,a1,a2,f1,f2;printf(输入收敛精度：);//输入标准收敛精度std::cinE0;do{D[0][0]=E[0][0];D[1][0]=E[1][0];A[0][0]=4*(D[0][0]-2)*(D[0][0]-2)*(D[0][0]-2)+2*D[0][0]-4*D[1][0];A[1][0]=-4*(D[0][0]-2*D[1][0]);//A[0][0],A[1][0]为原函数梯度的各项B[0][0]=1.0/(12.0*(D[0][0]-2)*(D[0][0]-2));B[0][1]=1.0/(24.0*(D[0][0]-2)*(D[0][0]-2));B[1][0]=1.0/(24.0*(D[0][0]-2)*(D[0][0]-2));B[1][1]=(6.0*(D[0][0]-2)*(D[0][0]-2)+1)/(48.0*(D[0][0]-2)*(D[0][0]-2));//B[0][0],B[0][1],B[1][0],B[1][1]分别代表原函数的海赛矩阵的逆阵的各项C[0][0]=-(B[0][0]*A[0][0]+B[0][1]*A[1][0]);C[1][0]=-(B[1][0]*A[0][0]+B[1][1]*A[1][0]);//C[0][0],C[1][0]为搜索方向dk的各项//下面利用外推法寻找函数2的区间，找单谷区间x1=0;x2=x1+h;y1=fun2(x1,D[0][0],C[0][0],D[1][0],C[1][0]);y2=fun2(x2,D[0][0],C[0][0],D[1][0],C[1][0]);if(y2y1){h=-h;x3=x1,y3=y1;x1=x2,y1=y2;x2=x3,y2=y3;}x3=x2+h;y3=fun2(x3,D[0][0],C[0][0],D[1][0],C[1][0]);while(y3y2){h=2*h;x1=x2,y1=y2;x2=x3,y2=y3;x3=x2+h;y3=fun2(x3,D[0][0],C[0][0],D[1][0],C[1][0]);i++;}//下面利用黄金分割法寻找函数2极值a=x1;b=x3;a1=b-k*(b-a);a2=a+k*(b-a);f1=fun2(a1,D[0][0],C[0][0],D[1][0],C[1][0]);f2=fun2(a2,D[0][0],C[0][0],D[1][0],C[1][0]);do{if(f1=f2){a=a1;a1=a2;f1=f2;a2=a+k*(b-a);f2=fun2(a2,D[0][0],C[0][0],D[1][0],C[1][0]);}else{b=a2;a2=a1;f2=f1;a1=b-k*(b-a);f1=fun2(a1,D[0][0],C[0][0],D[1][0],C[1][0]);}j++;}while(fabs((b-a)/b)=E0&&fabs((f2-f1)/f2)=E0);m=0.5*(a+b);//m为阻尼因子E[0][0]=D[0][0]+m*C[0][0];E[1][0]=D[1][0]+m*C[1][0];printf(%d%15f10%15f10\n,t,E[0][0],E[1][0],fun1(E[0][0],E[1][0]));t++;}while(fabs(E[0][0]-D[0][0])=E0&&fabs(E[1][0]-D[1][0])=E0);X[0][0]=E[0][0];X[1][0]=E[1][0];printf(迭代了%d次\n,t);printf(极小点(x1,x2)=(%f10,%f10)\n,X[0][0],X[1][0]);printf(极小值f(x1,x2)=%f10\n,fun1(X[0][0],X[1][0]));}程序运行结果：四结论由该程序的运行结果可知，要求迭代两次后函数的极小值点在（2.52939210,1.46320410）处，经验证，运算结果完全正确。验证了该程序的可行性。