优化课程实验

实验一一维搜索方法实验实验学时：2h实验类型：常规实验要求：必修一、 实验目的：掌握确定搜索区间的进退法，一维搜索的黄金分割法、牛顿法(切线法)、二次插值法的基本思想及迭代过程。掌握上述常用一维搜索方法的优缺点及程序编制的方法。掌握采用一维优化方法求优时如何提高搜索效率，减少盲目性。二、 实验内容1、 任选实验以下内容之一，自编程序；用0.618法求函数f(x)=x(x+2)在区间［-3,5］中的极小点，要求精度为£=0.01。.用0.618法求函数f(x)=x2+x+2的近似极小点及近似极小值，以致［a，b］=［-1，3］，缩短的相对精度为£=0.08。用牛顿法求函数f(x)=x4-5x3-4x2-6x+60的极小点及极小值，已知［a，b］=［3,4］，要求精度为£=0.01o用二次插值法求f(x)=x2+e-x的极小点及极小值，已知［a，b］=［0，1］，要求精度为£=0.001。用二次插值法求f(x)=3x4-16x3+30x2-24x+8的极小点及极小值，已知［a，b］=［0，3］，要求精度为£=0.001。2、 改变初始搜索区间或精度要求，观察程序运行结果及迭代次数的变化情况。三、 实验原理、方法和手段一维搜索方法的内容包括搜索区间的确定和区间缩小两个内容，首先给定初始点，应用外推法确定一个高一低一高的单谷区间，在采用黄金分割法或插值法进行区间缩小，当区间缩小到收敛精度。四、实验组织运行要求本实验为上机实验，采用单人单机，集中授课形式。五、实验条件人均计算机一台，安装有C语言软件。六、 实验步骤实验课前，明确实验目的和要求，预习相关内容，应用选择的算法语言编写程序；在C语言环境下输入和调试程序；运行程序，分析结果的正确性；改变搜索区间，运行程序，观察结果和迭代次数的变化；改变收敛精度，观察结果和迭代次数的变化。七、 思考题对黄金分割法和二次插值法的特点进行总结，具体使用时应如何选择？八、 实验报告1、 写清上机实验的名称及要求。2、 所选优化方法的基本原理简述。3、 绘制程序框图。4、 程序中变量及参数说明。5、 运行结果分析。(要求输出迭代过程中搜索区间上下限、迭代次数、极小点和极小值)#include<math.h>#include<stdio.h>doubleobfunc(doublex) /*目标函数*/(doubleff;ff=pow(x,2)+9*x;return(ff);}voidjts(doublex0,doubleh0,doubles[],intn,doublea[],doubleb[])/*进退法*/inti;doublex[3],h,fl,f2,f3:h=hO;for(i=0;i<n;i++)x[0]=x0;fl=obfunc(x[0]);for(i=0:i<n;i++)x[l]=x[0]+h^s[i];f2=obfunc(x[1]);if(f2>=fl){h=-hO;for(i=0;i<n;i++)x[2]=x[0]:f3=fl;for(i=0:i<n;i++){x[0]=x[l]:x[l]=x[2]:}fl=f2;f2=f3:}for(;;){h=2.0*h;for(i=0;i<n;i++)x[2]=x[l]+h^s[i]:f3=obfunc(x[2]);if(f2<f3)break;else(for(i=0;i<n;i++)(x[0]=x[1];x[1]=x[2];}f1=f2;f2=f3;}}if(h<0)for(i=0;i<n;i++)(a[i]=x[2];b[i]=x[0];}elsefor(i=0;i<n;i++)(a[i]=x[0];b[i]=x[2];}printf("%4d",n);}doublegold(doublea[],doubleb[],doubleeps,intn,doublexx)/*黄金分割法*/(inti;doublef1,f2,ff,q,w;doublex[3];for(i=0;i<n;i++)(x[0]=a[i]+0.618*(b[i]-a[i]);x[1]=a[i]+0.382*(b[i]-a[i]);}f1=obfunc(x[0]);f2=obfunc(x[1]);do(if(f1>f2)(for(i=0;i<n;i++)(b[i]=x[0];x[0]=x[1];}f1=f2;for(i=0;i<n;i++)x[1]=a[i]+0.382*(b[i]-a[i]);f2=obfunc(x[1]);}else(for(i=0;i<n;i++)(a[i]=x[1];x[1]=x[0];}f2=f1;for(i=0;i<n;i++)x[0]=a[i]+0.618*(b[i]-a[i]);f1=obfunc(x[0]);}q=0;for(i=0;i<n;i++)q=q+(b[i]-a[i])*(b[i]-a[i]);w=sqrt(q);}while(w>eps);for(i=0;i<n;i++)xx=0.5*(a[i]+b[i]);ff=obfunc(xx);printf("xx=ff=%5.2f,,,,%5.2f”,xx,ff);return(ff);}voidmain()(intn=1;doublea[1],b[1],xx;doubles[]={1},x0=0;doubleeps1=0.001,h0=0.1;jts(x0,h0,s,n,a,b);gold(a,b,eps1,n,xx);实验二 多维无约束优化一、 实验目的：掌握常用无约束优化方法的基本思想及迭代过程。掌握最速下降法、共轭梯度法、牛顿法、变尺度法、Powell法、单形替换法等的优缺点及程序编制方法。掌握如何合理选择无约束优化方法，保证迭代效率及解题效果。掌握如何保证所得结果为目标函数的全局最优解的一般方法。二、 实验内容求函数f(x)=2x12+3x22-8x1+1。的近似无约束最小点X*及最小值f(X*)。已知初始点X(o)=[1，2]t，£=0.001。求函数f(x)=1.5x12+0.5x22-x1x2-2x1的近似最小点X*及最小值f(X*)。已知初始点X(0)=[-2，4]t，£=0.02。求函数f(x)=x14-2x12x2-2x22-2x1X2+4.5x「4x2+4的近似最小点X*及最小值f(X*)。已知初始点X(0)=[-0.2，0.2]t，£=0.01。求f(x)=x12+2x22-4X1-2x1X2的近似最小点X*及最小值f(X*)。已知初始点X(0)=[1，1]t，£=0.01。求f(x)=4(x1-5)2+(x2-6)2的最优解。已知初始点X(o)=[8，9]t，£=0.01。三、 实验要求：任选实验内容中之一，自选实验目标2中的一种方法，自编程序；改变初始初始点及精度要求，观察最优解及迭代次数的变化情况。四、 实验报告：写清上机实验的名称及要求。所选优化方法的基本原理简述。绘制程序框图。程序中变量及参数说明。运行结果分析。对所采用的方法的搜索效率及搜索结果进行评价，并提出改进意见及建议。附：Powell法的参考程序：#include"math.h"#include"stdio.h"#include"malloc.h"doubleobfunc(doublex[])(doubleff;ff=x[0]*x[0]+2*x[1]*x[1]-4*x[0]-2*x[0]*x[1];return(ff);}voidjts(doubleX0[],doubleh0,doubles[],intn,doublea[],doubleb[])(inti;double*x[3],h,f1,f2,f3;

for(i=0;i<3;i++)h=h0;for(i=0;i<n;i++)f1=obfunc(x0);for(i=0;i<3;i++)h=h0;for(i=0;i<n;i++)f1=obfunc(x0);for(i=0;i<n;i++)f2=obfunc(x[1]);if(f2>=f1){h=-h0;for(i=0;i<n;i++)f3=f1;*(x[0]+i)=x0[i];*(x[1]+i)=*(x[0]+i)+h*s[i];*(x⑵+i)=*(x[0]+i);for(i=0;i<n;i++){*(x[0]+i)=*(x[1]+i);*(x[1]+i)=*(x[2]+i);}f1=f2; f2=f3;}for(;;){h=2.0*h;for(i=0;i<n;i++) *(x⑵+i)=*(x[1]+i)+h*s[i];f3=obfunc(x[2]);if(f2<f3)elsebreak;if(f2<f3)elsebreak;{for(i=0;i<n;i++){*(x[0]+i)=*(x[1]+i); *(x[1]+i)=*(x[2]+i);}f1=f2;f2=f3; }}if(h<0)for(i=0;i<n;i++) {a[i]=*(x⑵+i); b[i]=*(x[0]+i);}elsefor(i=0;i<n;i++) {a[i]=*(x[0]+i); b[i]=*(x[2]+i); }for(i=0;i<3;i++) free(x[i]);}doublegold(doublea[],doubleb[],doubleeps,intn,doublexx[]){inti;doublef1,f2,ff,q,w,*x[2];for(i=0;i<2;i++)x[i]=(double*)malloc(n*sizeof(double));for(i=0;i<n;i++){*(x[0]+i)=a[i]+0.618*(b[i]-a[i]); *(x[1]+i)=a[i]+0.382*(b[i]-a[i]);}f1=obfunc(x[0]);f2=obfunc(x[1]);do{if(f1>f2){ for(i=0;i<n;i++) {b[i]=*(x[0]+i);*(x[0]+i)=*(x[1]+i);}f1=f2;*(x[1]+i)=a[i]+0.382*(b[i]-a[i]);for(i=0;i<n;i++)*(x[1]+i)=a[i]+0.382*(b[i]-a[i]);f2=obfunc(x[1]);}else{for(i=0;i<n;i++)f2=f1;for(i=0;i<n;i++)f1=obfunc(x[0]);}}else{for(i=0;i<n;i++)f2=f1;for(i=0;i<n;i++)f1=obfunc(x[0]);}q=0;for(i=0;i<n;i++)w=sqrt(q);}while(w>eps);for(i=0;i<n;i++)ff=obfunc(xx);for(i=0;i<2;i++)return(ff);}{a[i]=*(x[1]+i); *(x[1]+i)=*(x[0]+i);}*(x[0]+i)=a[i]+0.618*(b[i]-a[i]);q=q+(b[i]-a[i])*(b[i]-a[i]);xx[i]=0.5*(a[i]+b[i]);free(x[i]);doubleoneoptim(doublex0[],doubles[],doubleh0,doubleepsg,intn,doublex[]){double*a,*b,ff;a=(double*)malloc(n*sizeof(double));b=(double*)malloc(n*sizeof(double));jts(x0,h0,s,n,a,b);ff=gold(a,b,epsg,n,x);free(a);free(b);return(ff);}doublepowell(doublep[],doubleh0,doubleeps,doubleepsg,intn,doublex[]){inti,j,m;double*xx[4],*ss,*s;doublef,f0,f1,f2,f3,fx,dlt,df,sdx,q,d;ss=(double*)malloc(n*(n+1)*sizeof(double));s=(double*)malloc(n*sizeof(double));for(i=0;i<n;i++){for(j=0;j<n;j++) *(ss+i*(n+1)+j)=0;*(ss+i*(n+1)+i)=1;}for(i=0;i<4;i++)xx[i]=(double*)malloc(n*sizeof(double));for(i=0;i<n;i++)*(xx[O]+i)=p[i];for(;;)(for(i=0;i<n;i++)(*(xx[l]+i)=*(xx[O]+i);x[i]=*(xx[l]+i);}fO=f1=obfunc(x);dlt=-l;for(j=0;j<n;j++)(for(i=0;i<n;i++)(*(xx[O]+i)=x[i];*(s+i)=*(ss+i*(n+l)+j);}f=oneoptim(xx[0],s,hO,epsg,n,x);df=fO-f;if(df>dlt)(dlt=df;m=j;}}sdx=O.O;for(i=0;i<n;i++)sdx=sdx+fabs(x[i]-(*(xx[l]+i)));if(sdx<eps)(free(ss);free(s);for(i=0;i<4;i++)free(xx[i]);return(f);}for(i=0;i<n;i++)*(xx⑵+i)=x[i];f2=f;for(i=0;i<n;i++)*(xx[3]+i)=2.0*(*(xx⑵+i)-(*(xx[i]+i)));x[i]=*(xx[3]+i);}fx=obfunc(x);f3=fx;q=(f1-2*f2+f3)*(f1-f2-dlt)*(f1-f2-dlt);d=0.5*dlt*(f1-f3)*(f1-f3);if((f3<f1)ll(q<d)){if(f2<=f3)for(i=0;i<n;i++)*(xx[0]+i)=*(xx⑵+i);elsefor(i=0;i<n;i++)*(xx[0]+i)=*(xx[3]+i);}else{for(i=0;i<n;i++){*(ss+(i+1)*(n+1))=x[i]-(*(xx[1]+i));*(s+1)=*(ss+(i+1)*(n+1));}f=oneoptim(xx[0],s,h0,epsg,n,x);for(i=0;i<n;i++)*(xx[0]+i)=x[i];for(j=m+1;j<=n;j++)for(i=0;i<n;i++)*(ss+i*(n+1)+j-1)=*(ss+i*(n+1)+j);}}}voidmain(){inti;intdimen=2;doublespoint[]={1,1};doubleff,x[20];doubleh0=0.3;doublefepsl=0.001;doublesepsl=0.0001;ff=powell(spoint,h0,fepsl,sepsl,dimen,x);for(i=0;i<dimen;i++)printf("x%d%5.2f\n”,i,x[i]);printf("optimum=%5.2f”,ff);}实验三约束优化方法实验实验学时：2h实验类型：常规实验要求：必修一、 实验目的1、 掌握常用多维约束优化方法的基本思想及特点；2、 掌握约束优化问题的一般求解方法；3、 掌握常用优化方法程序的使用方法；4、 通过上机练习，了解约束优化方法的搜索效率，搜索效果及解的可靠性。二、 实验内容1、任选一种约束优化方法，根据程序框图编写程序，上机调试并运行，对运行结果进行分析。用随机方向法求下面具有不等式约束的优化问题minf(x)=(x1-8)2+(x2-8)2S.T.g1(x)=-x1顷g2(x)=1-x2^0g3(x)=x1+x2-11^0的最优解。已知：a=0.4，£=0.01，x0=[2，3]T，[a，b]=[0.0，11.0]。用复合形法求约束最优解，目标函数及约束函数同1题。用惩罚函数法求下面具有不等式及等式约束的优化问题：minf(x)=4x1-x22-12S.T.g1(x)=34-10x1-10x2+x12+x22^0g2(x)=-x1^0g3(x)=-x2^0h(x)=x12+x22-25=0已知：x0=[3.0，3.0]T，r0=1.0，c=0.2，£=0.00005,[a1，b1]=[0.0，6.0]，[a2,b2]=[0.0，8.0]。2、改变化初始值，分析最优解的变化。三、实验原理、方法和手段选择间接法(惩罚函数法、拉格朗日乘子法)和直接法(随机方向发、复合形法、可行方向法)中的一种，理解基本原理和程序框图，考虑用选择的语言如何实现。MATLAB是“矩阵实验室(Matrixlaboratory)”的缩写，是由美国Mathworks公司推出的一种以矩阵运算为基础，集通用数学运算、图形交互、程序设计和系统建模为一体的著名软件，它具有功能强、使用简单、容易扩展等优点。与其他计算机语言相比，MATLAB表达方式与人们在数学、工程计算中常用的书写格式十分相似，它以解释方式工作，输入程序后就可得结果.人机交互性好，易于调试，用户学习和使用它都极为方便。MATLAB的基本部分有：矩阵运算和各种变换，代数和超越方程的求解，数据处理和傅里叶变换，数值积分等：除此之外，为了支持不同专业领域的用户，MATLAB还提供了大量的面向专业领域的工具箱，工具箱(toolbox)包含一系列专用的MATLAB函数库，以解决特定领域的问题。工具箱主要有：通讯工具箱(Communicationtoolbox)、控制系统工具箱(ControlSystemtoolbox)、信号处理工具箱(Signalprocessingtoolbox)、图像处理上具箱(Imaginetoolbox)、系统辨识工具箱(SystemIdentificationToolbox)、模糊逻辑工具箱(FuzzyLogictoolbox)、遗传算法最优化工具箱(geneticA1gorithmoptimizationToolbox)、最优化工具箱(optimizationToolbox)、数理统计工具箱(StatisticsToolbox)、小波分桥工具箱(waveletToo1box)等等。这些工具箱通常表现为M文件和高级MATLAB语言的集合形式，允许用户修改函数的源代码或增加新的函数来适应自己的应用，允许用户方便地综合使用不同工具箱中的技术来设计针对某个问题的用户解决方案。使用MATLAB语言和MATLAB工具箱，用户可以专注于算法研究，编程只需要几行就可以完成，而且可以很快的画出图形，从而迅速地进行多种算法的比较，从中找出最好的方案。MATLAB被广泛用于研究和解决各种具体的工程问题，它不光已成为欧美等发达国家各设计单位和科研部门的基本工具，也成为了各高等院校大学生、硕士生和博士生必须掌握的工具。优化上具箱(optimizationToolbox)涉及函数的最小化和最大化问题，也就

是函数的极值问题oMATLAB的优化工具箱内一些对普通非线性函数求解最小化成最大化极值的函数和解决诸如线性规划等标准矩阵问题的函数组成。下面简单地对优化上具箱中的基本函数作一介绍。优化工具箱中求非线性函数极小值的函数如下表所示。表一求非线性函数极小值的函数列表类型含义调用格式无条件标量问题minf(x)，其中x为标量xx=fmin(‘f’,x)无限定条件矩阵问题minf(x)，其中x为矩阵xx=fminu(‘f’,x)有限定条件问题minf(x)，条件为G(x)<0xx=constr(‘fg’,x)目标条件问题minf(x)，条件为加-foalxx=attgoal(‘f’，x,goal,w)最大最小极值min{maxf(x)},条件为G(x)-0xx=minimax(‘fg’,x)非线性二次平方极值minZ(f(x)xf(x))x=leastsq(‘f’,x)非线性方程f(x)=0x=fsolve(‘f’,x)半无穷条件minf(x) ,条件为x©(x,o)-0,VoX=seminf(‘ft’,n,x)问题的原形算法是Nelder-Mead单纯形搜索方法和BFGS拟牛顿(qusi-Newton)方法。限定条件下的最小、最大最小、目标法和半无穷优化等问题，所用的原理算法是二次规划法，非线件二次平方问题的原理算法是Gausss-Newton法和Levenberg—Marquardt法。非线性最小和非线性二次平方问题，可进行线性搜索策略的选择，线件搜索策略使用的是二次或三次内插和外插方法。优化上具箱处能解决几类求矩阵的极小值问题，此时仅需要将相应的系数矩陈和向量传递到函数中。MATLAB系统能解决的矩阵问题如表所示。类型含义语法非负二次平方问题minAx-b2，条件为x>0xx=nnls(A,b)二次问题.「1 K条件为min-xtHx-ctx，志斤如2xL匕 _1Ax<bx=qp(H,c,A,b)线性规划问题min(fTx)，条件为Ax<bxx=lp(f,A,b)函数功能举例1无限定条件的极值问题考虑如下的问题：求合适的集合[%,%],使问题min/(x)=exl(4x^4-2x|+4xrx2+2x24-1]q工X ' /成立。为求解此问题，可先编写一个能返回函数值的M文件，将函数表达式与写入MATLAB系统中，然后调用非限定最小程序fminu。具体的步骤如下；利用文件编辑器编写M文件functionf=fun(x).f=exp(x(1))*(4*x(1)^2+2 (2)a2-f4*x(1)Xx(2)+2*x(2)+l);程序中functionf=fun(x)表示将此M文件定义为名为fun(x)的文件。2)在命令窗口中调用优化程序fminux=[-l,1]； %估计初值x=fminu(‘fun’,X)经过36次函数调用后，得到极值点为x=0.5000 -1.0000在命令窗口中接着输入fun(x)得到极值点处的函数值为：ans1.3029e-0l0注意：当函数存在一个以上的局部最小点时，的初值会影响选代的次数和解的值，本例中的初值为[-1，1]。为了改变优化解的特征，可在函数fminu的说明参数中加入可选的参数说明(options)，语句如下：X=fminu['fun',x,options)options是一个包含允许误差值和算法选择等选项的向量。它的第一个元素所起的作用是控制在优化周期中显示输出的数量，将它设置为1时会以表格形式显示函数值和收敛信息。options的第二和第三个元素建立终止判据。options中的其他元素用于选择算法和设置迭代的最大次数等。2限定条件极值问题使constr函数，可在上面例子中加入不等式限定条件。如对问题2minf(x)=exl(4x^+2对+4x,x2+2x2+1)限定条件为：1.5+x{x2一巧—x2<0—X]x2—10<0将上面编写的M文件进行修改，使其能返回目标函数值和限定函数值下的极值函数constr求出极值。具体过程如下：function■[f,g]=fun(X)f=exp(x(l))*(4*x(1)A2+2*x(2)A2+4*x(1) (2)+2*x(2)+1);g(l)=1.5+x(l)*x(2)-x(l)-x(2);g(2)—x(l)*x(2)-10;在命令窗口中调用优化程序x=[-11]； %估计初值x=constr('fun'，x)得到的结果为：X=-9.5474 1.0474在命令窗口中接着输入[f,g]=fun(x)得到的结果是在极值点处的函数值和限定条件的值f=0.0236g=-0.8882 03下界条件与上界条件使用函数的有界语法可将变量x限定在某一区域之内。例如对于函数constr,下面的语句X=constr('fun',x,optlons,vlb,vub)可以将x限定在vlbVxVvub范同之内。如在例2中，要将x限定为大于0,可在命令窗口中，使用如下的命令

options=[];与使用默认选顼options=[];与使用默认选顼%没有-上界v'jb=[]f%没有-上界r('funTtx,options,vlb,vub)经过七坎选代后，问题的解为：X=0 1.5000在MATLAB命令行窗口中输入如下内容：[f,g]=fun(x)得到的结果是在极值点处的函数值和限定条件的值：f=8.5000g=0 -10在上面的例子中，x没有上升，因此vub被置为空矩阵。也可以使用如下形式的命令省掉上界，其调用格式为：r(rfunTfxroptionsf )当vlb或vub的元素数目比向量x的元素数目少时，则x中只有前n个元素被限定有界。另外，上界和下界也可以用线性不等式表示。当仅有几个元素有界时，用这样的方法是比较方便的。具体的调用格式为：上界：x<U，表示为：x-U<0下界：x下界：x>L表示为：-x+L<04梯度计算一般说来，最小化程序都要使用由有限微分近似法计算出的梯度，这一过程会自动作用到每个变量以计算函数和限定条件的偏导数。若所编写的程序能返回函数和限定条件的偏导数，则问题的求解将会更精确和有效。用梯度来解决的前面的问题2具体过程如下：(l)利用文件编辑器为目标函数编写M文件。function[f,g]=fun(x)f^exp(x(l))* ^2+2*x(2)A2+4*x(1)*x(2)+2*x(2)-»-l);g(l)-l.5+x(1)*x(2)-x(l)-x(2);g(2)=-x(l)*x(2)-10;(2)编写求梯度的M文件。function[dffdg]=grad(x)t^exp(x(l))*(4*x(l)A2+2*x(2)A2+4*x(1) (2)+2*x(2)+1);df=[t+4*exp(x(D)*(2*x(l)+x(2)),4*exp(x(l))*(x(l)+x(2)+0.5)]; 电定义目标函数的梯度函数dg=[x⑵⑵x(l)-1,-x(1)]; 号定义限定函数的梯度在命令窗口中调用函数求限定条件下的极值问题。1] %初始估值options=[] 若使用默认选项vlb-[]； %没有■下界vub=[]; %没有上界x^constr1.1fun' options,vlb,vubf/