《实验报告阻尼牛顿法.docx》由会员分享，可在线阅读，更多相关《实验报告阻尼牛顿法.docx（13页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、太原理工高校机械学院机测系课程上机试验报告课程名称：机械优化设计班 级日 期成果评 定姓 名试验 室老师签 名试 验名 称用阻尼牛顿法程序解题所 用 软 件C+ DEV实 验 目 的 及 内 容试验目的：L把握并能够建立最优化基本类型问题的数学模型。2 .把握最优化方法的基本概念、基本理论和基本方法，奠定最优化 的理论基础。3 .能够娴熟编制和调试最优化方法的程序，奠定解决实际中的优化 问题的基础试验内容：理解阻尼牛顿法并编写相关程序求其最优解。牛顿法程序考核题于(X) = 4(玉 +1)2 + 2(x2 -1)2 + . +(2 +10初始点：X =0,0T,梯度精度 = 0.01ebs=s
2、qrt(ebs)/(N*N* 100);for(i=0;iN;i+)bii=1.0;for(i=0;iN;i+)(if(fabs(aii)ebs)(for(k=i;kebs)(for(m=0;mN;m+)(fact=aim;aim=akm;akm=fact;fact=bim;bim=bkm; bkm=fact;) break;)fact=aii;for(k=0;kN;k+)(aik/=fact;bik/=fact;)for(k=0;ki;k+)(fact=aki;for(m=0;mN ;m+)(aim-=aifm*fact;bi m-=bi m *fact;)for(k=i+l；kN;k+)(
3、fact=aki;for(m=0;mN ;m+)(ai m-=ai m *fact;return;)/*阻尼牛顿法*/void newtow_method(double x,double h,double ebsin,int yw) (double dN,t,f,sN,heiNN,hei_lNN;int i,j,k=l;fprintf(fp,”* 阻尼牛顿法计算结果 *nn);ywddf(yw);/* 输出一维迭代方法*/fprintf(通”初始坐标:n);f=xfout(x,0);fprintf(fp,nn);gread(x,d);do(heisen(x,hei);atoa_l(hei,he
4、i_l);for(i=0;iN;i+)(si=0.0;for(j=0;jN;j+)si+=hei_l皿*d;)switch(yw)(goldcut(x,s,h,ebsin);break;case 1: j=rccz(x,s,h,ebsin);break;)gread(x,d)；fprintf(fp,”迭代轮数 k=%3d nM,k);f=xfout(x,k);t=0;for(i=0;iebsin);fprintf(fp n*n)fprintf（fp,”阻尼牛顿法法优化最优点及目标函数值为：n）f=xfout(x,-l);fprintfCfp/迭代精度：”)；fprintf(fp,15.9mn”
5、,t);return;main()(double xON,h,ebsin;int yw;csd_x(xO); /* 初始坐标*/h=H_QJ;ebsin=EPSIN;yw=Y_F;fp=fopen(outname,wn);newtow_method(xO,h,ebsin,yw); /* 调用阻尼牛顿法*/ fclose(fp);return 0;备注：不交此报告者，本次试验为“不合格”。试验原理:实 验 原 理 步 骤试验步骤：1,画流程图，编写程序；2,将目标函数代入；3,编译运行，将结果保存*阻尼牛 顿法 计 算 结果*+一维搜寻方法：黄金分割法+Tf( 0)=X( 0)=16.00000
6、000.0000000,迭代轮数k= 1x( 1)=-1.1249988,f( 1)二9.8125000迭代精度：0.000009986迭代轮数k= 2x( 2)=-1.1250000,f( 2)二9.8125000迭代精度：0.000000114实*1* 7，*1*1*1* *1* *1*rTw rTw eiw卜卜卜.卜.卜卜4、rTw rTw *Tw一卜验阻尼牛顿法法优化最优点及目标函数值为：x( * )=-1.1250000,结f( * )=9.8125000迭代精度：0.000000114果初始坐标:0.0000000,0.7499992,0.7500000,0.7500000,算法程
7、序实现/*csssqj.cpp */#include #include #include #include #include #define N 2#define EPSIN 0.000001#define H_QJ 1.0#define Y_F 1/*优化设计维数*/*迭代精度*/*初始区间搜寻步长*/* 一维搜寻方法选择：1 黄金分割法*/*2二次插值法*/FILE *fp;char outname50=阻尼牛顿法计算结果/*计算结果输出文件*/东给出初始点坐标*/void csd_x(double x0)(int i;for(i=0;iN;i+)/*初始点为坐标原点的状况*/x0i=0.
8、0;return;)/*目标函数号double hanshu(double x |) (double f;/*阻尼/*阻尼f=4.0*(x0+ 1.0)*(x0+1 .0)+2,0*(xl -1.0)*(xl -1.0) +x0+xl+10.0;牛顿法*/return f;)/*函数的梯度*/void gread(double x,double gff)(/*阻尼牛顿法*/gf0=8.0*(x0+1.0)+1.0;gfl=4.0*(xl-1.0)+1.0;return;)/*Heisen 矩阵*/void heisen(double x,double heiflfN)(hei00=8.0;he
9、iOJl=O.O;heil|0=0.0;heill=4.0;return;)/*计算 f(xk+as)*/double xkadd(double x,double d,double a)(int i;double xl|N;for(i=0;iN;i+)xli=xi4-a*di;return hanshu(xl);)/*输出选定的一维迭代方法*/void ywddf(int yw)(switch(yw)(黄金分割法+nn);二次插值法+nn”);fprintf(fp, +一维搜寻方法: break;fprintf(fp, ”+一维搜寻方法: break;)return;)/*输出当前迭代点坐标及
10、目标函数值*/double xfout(double x,int m)(intj;double f;f=hanshu(x);if(-l=m)fprintf(fp；x( * )=)；elsefprintf(fp,x(%3d)=n,m);for(j=0;jf2)(a2=a3;a3=0.0;fl=f2;f2=f3;f3=fl;h=h;)do(al=a2;a2=a3;fl=f2;f2=f3;a3=a2+h;f3=xkadd(x,d,a3);h=2*h;while(f30.0)abO=al;abl=a3;)else(ab0=a3;ab|l=al;)return;)/*黄金分割法*/void goldcu
11、t(double x,double df,double h,double ebsin) (double a 1 ,a2,f 1 ,f2,a,b,ab2;int i;csssqj(x,d,h,ab);a=ab0;b=abl;al=b-0.618*(b-a);fl=xkadd(x,d,al);a2=a+0.618*(b-a);f2=xkadd(x,d,a2);do(if(flf2)(a=al;al=a2;fl=f2;a2=a+0.618*(b-a);f2=xkadd(x,d,a2);)else(b=a2;a2=al;f2=fl;al=b-0.618*(b-a);fl=xkadd(x,d,al);w
12、hile(b-aebsin);for(i=0;iN;i+)xi+=(b+a)*di/2;return;)/*二次插值法*/int rccz(double x,double d,double h,double ebsin) (double a 1 ,a2,a3,a4,fl ,f2,f3,f4,c l,c2,ab2;int i,p=0,k=0;csssqj(x,d,h,ab);al=ab0;a3=abl;a2=(al+a3)/2;fl=xkadd(x,d,al);f2=xkadd(x,d,a2);f3=xkadd(x,d,a3);while(l)(cl=(f3-fl)/(a3-al);c2=(f2
13、-fl)/(a2-a l)-c l)/(a2-a3);if(0=c2)(P=l；break;a4=0.5*(a 1 +a3-c 1 /c2);if(a4-al )*(a3-a4)=0.0)(P=2；break;)f4=xkadd(x,d,a4);if(l =k&fabs(a4-a2)=ebsin)break;if(a4f4)(fl=f2;al=a2;k=l;a2=a4;f2=f4;)else(a3=a4;f3=f4;)else(if(f2f4)(f3=f2;a3=a2;k=l;a2=a4;f2=f4;)else(al=a4;fl=f4;if(O=p)fl=f4-f2?a4:a2;elsefl=a2;for(i=0;iN;i+) xi+=return p;)/*矩阵求逆*/void atoa_l(double aN,double bN) (int i,j,k,m;double fact,ebs=0.0;for(i=0;iN;i+)for(j=0;jN;j4-+)biUl=O.O;ebs+=aij*aij;