数学实验8月13日作业

1、1.取不同的初值计算下列平方和形式的非线性规划，尽可能求出所有局部极小点，进而找出全局极小点，并对不同算法（搜索方向、搜索步长、数值梯度与分析梯度等）的结果进行分析、比较。(2). ,(4).，其中.解：(2)由题意设先求（文件名为a.m），再求（文件名为b.m）,解得近似解为或或或或或或或或,取初值（exam0705grad_run.m）利用matlab计算可得以下表中数据：搜索方向步长搜索最优解最优解最优值目标函数调用次数Case 1bfgs混合二、三次插值-2.8996e+0002.1584e-0014.2939e-0065.0000e+000Case 2dfp-2.8996e+0002

2、.1584e-0014.2939e-0065.0000e+000Case 3bfgs-2.8996e+0002.1584e-0014.2939e-0065.0000e+000Case 4steepdesc-2.8996e+0002.1584e-0014.2939e-0065.0000e+000Case 5bfgs三次插值-2.8996e+0002.1584e-0014.2939e-0065.0000e+000Case 6dfp-2.8996e+0002.1584e-0014.2939e-0065.0000e+000Case 7bfgs-2.8996e+0002.1584e-0014.2939e

3、-0065.0000e+000Case 8steepdesc-2.8996e+0002.1584e-0014.2939e-0065.0000e+000当取初值（exam07051grad_run.m）利用matlab计算可得以下表中数据：搜索方向步长搜索最优解最优解最优值目标函数调用次数Case 1bfgs混合二、三次插值3.4578e+000-7.9722e-0028.3883e-0121.4000e+001Case 2dfp3.4578e+000-7.9722e-0027.2387e-0101.6000e+001Case 3bfgs3.4578e+000-7.9722e-0028.3883

4、e-0121.4000e+001Case 4steepdesc3.4578e+000-7.9722e-0025.0577e-009 1.7000e+001Case 5bfgs三次插值3.4578e+000-7.9722e-0028.3883e-0121.4000e+001Case 6dfp3.4578e+000-7.9722e-0027.2387e-0101.6000e+001Case 7bfgs3.4578e+000-7.9722e-0028.3883e-0121.4000e+001Case 8steepdesc3.4578e+000-7.9722e-0025.0577e-009 1.700

5、0e+001因为不同步长搜索到相同的最优解，但最优值不同，所以此时初值离最优解很远。继续取不同初值求解（只需在程序中不断改变初值即可）。当初值取时可得到以下数据搜索方向步长搜索最优解最优解最优值目标函数调用次数Case 1bfgs混合二、三次插值-2.1181e+001-3.6468e+0011.3691e-0061.7000e+001Case 2dfp-2.1158e+001-3.6504e+0011.9079e+0001.9500e+002Case 3bfgs-2.1181e+001-3.6468e+0011.3691e-0061.7000e+001Case 4steepdesc-2.11

6、74e+001-3.6478e+0011.5885e-001 7.7100e+002Case 5bfgs三次插值-2.1181e+001-3.6468e+0011.3691e-0061.7000e+001Case 6dfp-2.1158e+001-3.6504e+0011.9079e+0001.9500e+002Case 7bfgs-2.1181e+001-3.6468e+0011.3691e-0061.7000e+001Case 8steepdesc-2.1174e+001-3.6478e+0011.5885e-001 7.7100e+002综上全局极小点为，局部极小点，。（2）8 取不同的

7、初值计算下列非线性规划，尽可能求出所有局部极小点，进而找出全局极小点，并对不同算法（搜索方向、搜索步长、数值梯度与分析梯度等）的结果进行分析、比较。 (2),其中，,。解：(2)1. 设计程序如examfun.m.建立功能函数，便于调用。functionf,g=examfun(x)a1=4,4;a2=2.5,3.8;c=0.7,0.73;f=1/(x-a1)*(x-a1)+c(1)-1/(x-a2)*(x-a2)+c(2);if nargout1 g(1)=2*(x(1)-a1(1)/(x(1)-a1(1)2+(x(2)-a1(2)2+c(1)2+2*(x(1)-a2(1)/(x(1)-a2(

8、1)2+(x(2)-a2(2)2+c(2)2 %该为梯度函数的两个方向g(2)=2*(x(2)-a1(2)/(x(1)-a1(1)2+(x(2)-a1(2)2+c(1)2+2*(x(2)-a2(2)/(x(1)-a2(1)2+(x(2)-a2(2)2+c(2)2end设计脚本程序如exam.m，opt1=optimset(LargeScale,off,MaxfunEvals,1000,Tolfun,1e-8,TolX,1e-8,Gradobj,on); opt2=optimset(opt1,HessUpdate,dfp); opt3=optimset(opt1,HessUpdate,steep

9、desc); opt4=optimset(opt1,lineSearchType,cubicpoly); opt5=optimset(opt5,HessUpdate,dfp);opt6=optimset(opt5,HessUpdate,steepdesc); %x0=3,3; %x0=4,4;%x0=5,5;%x0=2,7;%x0=10,10;%x0=50,20;%x0=1,10;%x0=50,100;%x0=43,2;x0=20,43; x1,v1,exit1,out1=fminunc(examfun,x0,opt1) x2,v2,exit2,out2=fminunc(examfun,x0,

10、opt2)x3,v3,exit3,out3=fminunc(examfun,x0,opt3)x4,v4,exit4,out4=fminunc(examfun,x0,opt4)x5,v5,exit5,out5=fminunc(examfun,x0,opt5)x6,v6,exit6,out6=fminunc(examfun,x0,opt6)运行结果及分析。得出答案；制成表格初值设定为（20，43）搜索方向 步长搜索最优解最优值 目标函数调用次数BFGS混合二、三次插值0.2024e-004 0.4775e-0042.0423e-00512DFP混合二、三次插值0.2024e-004 0.4775e

11、-0042.0423e-00512最速下降混合二、三次插值0.2024e-004 0.4775e-0042.0423e-00512BFGS三次插值0.2024e-004 0.4775e-0042.0423e-00512DFP三次插值0.2024e-004 0.4775e-0042.0423e-00512最速下降三次插值0.2024e-004 0.4775e-0042.0423e-00512 可以看到采用分析方法(给定梯度函数) 比起默认的数值方法得到结果所需要的调用次数明显减少。DFP算法的函数调用次数较少。 改变初值，采用默认搜索方向BFGS，默认搜索步长，给定梯度函数，得到的结果如下:X(

12、1)X(2)最优值迭代次数(3,3)2.49333.7345-1.03292(4,4)3.67113.95610.76401(5,5)5 50.25161(2,7)2.98994.1982-0.31781可以看到，改变初值对结果的影响很大.结论： 由上面计算结果得出，无论采用那种计算方法最终都可以得到函数的极小值的，而采用DFP算法函数的调用次数较少；分析法计算梯度函数效果较好；最终结果取哪一个局部极值同初值密切相关。 5. 某分子由25 个原子组成，并且已经通过实验测量得到了其中某些原子对之间的距离（假设在平面结构上讨论），如下表所示。请你确定每个原子的位置关系。原子对距离原子对距离原子对距

13、离原子对距离（4,1）0.9607(5,4)0.4758(18,8)0.8363(15,13)0.5725(12,1)0.4399(12,41.3402(13,9)0.3208(19,13)0.7660(13,1)0.8143(24,4)0.7006(15,9)0.1574(15,14)0.4394(17,1)1.3765(8,6)0.4945(22,9)1.2736(16,14)1.0952(21,1)1.2722(13,6)1.0559(11,10)0.5781(20,16)1.0422(5,2)0.5294(19,6)0.6810(13,10)0.9254(23,16)1.8255(16

14、,2)0.6144(25,6)0.3587(19,10)0.6401(18,17)1.4325(17,2)0.3766(8,7)0.3351(20,10)0.2467(19,17)1.0851(25,2)0.6893(14,7)0.2878(22,10)0.4727(20,19)0.4995(5,3)0.9488(16,7)1.1346(18,11)1.3840(23,19)1.2277(20,3)0.8000(20,7)0.3870(25,11)0.4366(24,19)1.1271(21,3)1.1090(21,7)0.7511(15,12)1.0307(23,21)0.7060(24,3

15、)1.1432(14,8)0.4439(17,12)1.3904(23,22)0.8052解:问题分析每个原子的位置都是未知的，在坐标系中只有相对的位置参数，不妨固定原子1的坐标为。并且分子可以在平面内任意旋转，很难确定每个原子的绝对位置，表中给出了52组数据，可以得到52个方程，而未知数个数为48个（每个原子的，坐标）。属于超定方程组，没有确定的解，只能求得最优解。故采用最小方差标准来求最小二乘解得到最优解。模型建立设第个点所在的位置为，因为所求的是各原子间的位置关系，可以设定第一个点坐标为，然后再计算其他原子的位置，使它在最大程度上满足上表中提供的数据，即让达到最小，求得的坐标值就是最优的

16、解。其中表示第个原子和第个原子之间的距离，数据如上表所示，问题转化为无约束优化；算法设计根据题目中所给的数据，利用matlab编写函数文件（文件名为diatance.m）,再主程序为exam5.m,输出结果每个原子的位置如下（即第二个原子的位置为（1.2378,0.0746），第三个原子为（1,6142,1,。1211），）p = 0.3339 -1.2619 1.5092 -1.3714 0.8438 -0.4442 0.6956 -0.8857 1.3388 -0.6804 0.6158 -0.6627 0.9014 -0.9378 0.3552 -0.4180 1.2183 -0.652

17、6 0.6709 -0.4580 -0.4305 -0.0207 0.2989 -0.7737 0.5112 -0.7080 0.5654 -0.2765 0.6723 -1.7884 -0.1065 -1.3723 1.2895 -1.6922 0.9565 -1.2228 0.9735 -0.7787 1.2460 -0.2878 1.4318 -1.0946 1.9280 -0.4652 1.5063 -0.2325 0.9215 -0.8439为了形象直观地表示出各点的位置，画出如下的散点图7. 经济学中著名的柯布道格拉斯（CobbDouglas）生产函数的一般形式为，其中,分别表示产

18、值、资金、劳动力，式中,要由经济统计数据确定。现有中国统计年鉴（2003）给出的统计数据如下表，请用非线性最小二乘拟合求出式中的，并解释，的含义。年份总产值（万亿元）资金（万亿元）劳动力（亿人）19840.71710.09104.817919850.89640.25434.987319861.02020.31215.128219871.19620.37925.278319881.49280.47545.433419891.69090.44105.532919901.85480.45176.474919912.16180.55956.549119922.66380.80806.615219933

19、.46341.30726.680819944.67591.70426.745519955.84782.00196.806519966.78852.29146.895019977.44632.49416.982019987.83452.84067.063719998.20682.98547.139420009.94683.29187.208520019.73153.73147.3025200210.47914.35007.3740其中总产值取自“国内生产总值”，资金取自“固定资产投资”，劳动力取自“就业人员”。解：模型分析题目已经确定了此数据的拟合形式为Cobb-Douglas生产函数方程，再以

20、这些数据分别给出线性与非线性的最小二乘拟合，则确定式中的,。 在构造函数时，以向量的分量分别表示,以矩阵的第一列表示资金，第二列表示劳动力。非线性拟合程序代码先用mutlab编写函数（文件名为examfun7.m）,主程序为exam7.m,输出结果如下x = 0.8345 0.7736 0.7312norm = 2.7489res = Columns 1 through 8 -0.3048 0.0404 0.1000 0.1339 0.1254 -0.1437 -0.0866 -0.0577 Columns 9 through 16 0.1528 0.6527 0.4133 -0.0455 -

21、0.2861 -0.4392 -0.0194 -0.0217 Columns 17 through 19 -1.0567 0.1571 0.7350ef = 3out = firstorderopt: 4.7566e-004 iterations: 14 funcCount: 60 cgiterations: 0 algorithm: large-scale: trust-region reflective Newton message: 1x460 char非线性拟合给出的结果，, ,，其误差平方和，稍大。在接下来的分析中将看到它与原数据的拟合程度。 线性拟合程序代码1. 作图观察QKL的趋

22、势用mutlab编写程序（文件名为exam71.m）,运行，得到如下图像2. 对数处理从形式上观察，我们不妨将生产函数适当处理，将所给的Cobb-Douglas方程两边取对数，则可将非线性方程化为线性形式：再令,则方程可以转化为。令，则问题转化为超定方程根据最小二乘准则计算可以得到矩阵中的待定系数。由于要求线性拟合，故采用lsqlin 函数。 3. 线性拟合代码用matlab编写线性拟合代码（文件名为exam72.m），输出结果如下，x = -0. 0.73489 0.64329 resnorm = 0.35783 residual = Columns 1 through 4 0.42453

23、-0.13017 -0.16924 -0.17176 Columns 5 through 8 -0.13504 0. 0. -0. Columns 9 through 12 -0. -0.16972 -0. 0. Columns 13 through 16 0. 0. 0. 0. Columns 17 through 19 0.15767 0. -0. exit = 3 out = iterations: 6 algorithm: 1x43 char firstorderopt: 0.0114 cgiterations: 5 message: 1x123 char a = 0.99338 alpha = 0.73489 beta = 0.64329 4.直接左除程序名为exam73.m,代码如下：k=log(K);l=log(L); f=log(Q); M= ones(19,1),k,l; s=Mf a=exp(s(1) alpha=s(2) beta=s(3)输出结果为：s = -1.1464 0.6605 1.2627a = 0.3178alpha = 0.6605beta =1.2627不符合实际情况，应是左除无法加入参数条件限制导致的。下对此情况不予讨论4. 作图对比程序为exam74.m,代码如下：for t=1:19 Q1(t)=0.8