matlab教程第五讲

第五讲 数值计算(二)一、 线性优化用命令x=lp(C,A,b,vlb,vub)。

[例] 最小值线性优化 f(x)=-5x 1-4x 2-6x 3 x 1-x 2+x 3≦20 3x 1+2x 2+4x 3≦42 3x 1+2x 2≦30(0≦x 1, 0≦x 2,0≦x 3)First, enter the coefficients: f = [-5; -4; -6] A = [1 -1 1 3 2 4 3 2 0]; b = [20; 42; 30];lb = [0,0,0]; % x 的最小值 [0,0,0] ub = [inf,inf,inf];Next, call a linear programming routine: x= lp(f,A,b,lb,ub);Entering x x =0.0000 15.0000 3.0000实际此命令改为：x = linprog(f,A,b,Aeq,beq)x = linprog(f,A,b,Aeq,beq,lb,ub) 对以上的问题可做如下的操作： First, enter the coefficients:f = [-5; -4; -6];vubx vlb ..,min ≤≤≤∈b Ax t s Rx x CnTA = [1 -1 13 2 43 2 0];b = [20; 42; 30];lb = zeros(3,1);Next, call a linear programming routine:[x,fval,exitflag,output,lambda] = linprog(f,A,b,[],[],lb)；x =0.000015.00003.0000fval =-78.0000exitflag =1output =iterations: 6cgiterations: 0algorithm: 'lipsol'lambda =ineqlin: [3x1 double]eqlin: [0x1 double]upper: [3x1 double]lower: [3x1 double][例]线性优化Min -400x1-1000x2-300x3+200x4-2x2 + x3 + x4=02x1 +3x2 <=163x1 +4x2 <=24x1, x2, x3, x4>=0; x3<=5c=[-400,-1000,-300,200]; %目标函数系数A=[0 -2 1 1; 2 3 0 0; 3 4 0 0]; %约束条件系数 b=[0; 16; 24];xLB=[0,0,0,0]; % x 取值范围的最小值 xUB=[inf,inf,5,inf]; % x 取值范围的最大值 x0=[0,0,0,0]; % x 取迭代初始值nEq=1; % 约束条件中只有一个 = 号,其余为<= x=lp(c,A,b,xLB,xUB,x0,nEq)disp(['最优值为: ',num2str (c*x)])结果: x =3.4483 3.0345 5.0000 1.0690最优值为 -5700二、 非线性优化用命令x=constr('f ',x0)。

[例] 最小值非线性优化Min f(x)=-x1x2x3, -x1-2x2-2x3≤0, x1+2x2+2x3≤72,初值: x = [10; 10; 10]第一步：编写M 文件 myfun.m function [f,g]=myfun(x) f=-x(1)*x(2)*x(3);g(1)=-x(1)-2*x(2)-2*x(3); g(2)=x(1)+2*x(2)+2*x(3)-72;x g ≤∈)(..),(min t s R x x f n第二步：求解在MATLAB工作窗中键入x0=[10,10,10];x=constr('myfun',x0)即可.x =24.0000 12.0000 12.0000[例] 非线性优化Min f(x)=-x1x2(x1+ x2)x3<=0;x1, x2>=0; x3>=2;第一步：编写M文件 fxxgh.mfunction [F,G]=fxxgh(x)F=-x(1)*x(2);G(1)=(x(1)+x(2))*x(3)-120;第二步：求解在MATLAB工作窗中键入x=[1,1,1]; % x取迭代初始值options(13)=0; % 约束条件中有0个 = 号,其余为<= XL=[0,0,2]; % x取值范围的最小值XU=[inf;inf;inf]; % x取值范围的最大值[x,options]=constr('fxxgh',x,options,XL,XU);options(8) %输出最小值xans =-900.0000x =30.000030.00002.0000三、曲线拟合与插值曲线拟合和插值函数polyfit(x, y, n) 对描述n阶多项式y=f(x)的数据进行最小二乘曲线拟合interp1(x, y, xo) 1维线性插值interp1(x, y, xo, ' spline ') 1维3次样条插值interp1(x, y, xo, ' cubic ') 1维3次插值interp2(x, y, Z, xi, yi) 2维线性插值interp2(x, y, Z, xi, yi, ' cubic ') 2维3次插值interp2(x, y, Z, xi, yi, ' nearest ') 2维最近邻插值在大量的应用领域中，人们经常面临用一个解析函数描述数据(通常是测量值)的任务。

对这个问题有两种方法。

在插值法里，数据假定是正确的，要求以某种方法描述数据点之间所发生的情况。

曲线拟合或回归是人们设法找出某条光滑曲线，它最佳地拟合数据，但不必要经过任何数据点。

图1说明了这两种方法。

标有'o'的是数据点；连接数据点的实线描绘了线性内插，虚线是数据的最佳拟合。

曲线拟合曲线拟合涉及回答两个基本问题：最佳拟合意味着什么？应该用什么样的曲线？可用许多不同的方法定义最佳拟合，并存在无穷数目的曲线。

所以，从这里开始，我们走向何方？正如它证实的那样，当最佳拟合被解释为在数据点的最小误差平方和，且所用的曲线限定为多项式时，那么曲线拟合是相当简捷的。

数学上，称为多项式的最小二乘曲线拟合。

如果这种描述使你混淆，再研究图1。

虚线和标志的数据点之间的垂直距离是在该点的误差。

对各数据点距离求平方，并把平方距离全加起来，就是误差平方和。

这条虚线是使误差平方和尽可能小的曲线，即是最佳拟合。

最小二乘这个术语仅仅是使误差平方和最小的省略说法。

在MA TLAB中，函数polyfit求解最小二乘曲线拟合问题。

为了阐述这个函数的用法，让我们以上面图1中的数据开始。

»x=[0 .1 .2 .3 .4 .5 .6 .7 .8 .9 1];»y=[-.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2];为了用polyfit，我们必须给函数赋予上面的数据和我们希望最佳拟合数据的多项式的阶次或度。

如果我们选择n=1作为阶次，得到最简单的线性近似。

通常称为线性回归。

相反，如果我们选择n=2作为阶次，得到一个2阶多项式。

现在，我们选择一个2阶多项式。

» n=2; % polynomial order»p=polyfit(x, y, n)p =-9.8108 20.1293 -0.0317polyfit的输出是一个多项式系数的行向量。

其解是y = －9.8108x2＋20.1293x－0.0317。

为了将曲线拟合解与数据点比较，让我们把二者都绘成图。

ezplot('-9.8108*x*x+20.1293*x-0.0317')» xi=linspace(0, 1, 100); % x-axis data for plotting» z=polyval(p, xi);为了计算在xi数据点的多项式值，调用MA TLAB的函数polyval。

» plot(x, y, ' o ' , x, y, xi, z, ' : ' )画出了原始数据x和y，用'o'标出该数据点，在数据点之间，再用直线重画原始数据，并用点' : '线，画出多项式数据xi和z。

» xlabel(' x '), ylabel(' y=f(x) '), title(' Second Order Curve Fitting ')将图作标志。

这些步骤的结果表示于前面的图1中。

多项式阶次的选择是有点任意的。

两点决定一直线或一阶多项式。

三点决定一个平方或2阶多项式。

按此进行，n+1数据点唯一地确定n阶多项式。

于是，在上面的情况下，有11个数据点，我们可选一个高达10阶的多项式。

然而，高阶多项式给出很差的数值特性，人们不应选择比所需的阶次高的多项式。

此外，随着多项式阶次的提高，近似变得不够光滑，因为较高阶次多项式在变零前，可多次求导。

例如，选一个10阶多项式» pp=polyfit(x, y, 10) ;» format short e % change display format» pp.' % display polynomial coefficients as a columnans =-4.6436e+0052.2965e+006-4.8773e+0065.8233e+006-4.2948e+0062.0211e+006-6.0322e+0051.0896e+005-1.0626e+0044.3599e+002-4.4700e-001要注意在现在情况下，多项式系数的规模与前面的2阶拟合的比较。

还要注意在最小(-4.4700e-001)和最大(5.8233e+006)系数之间有7个数量级的幅度差。

将这个解作图，并把此图与原始数据及2阶曲线拟合相比较，结果如何呢？» zz=polyval(pp, xi); % evaluate 10th order polynomial» plot(x, y, ' o ' , xi, z, ' : ' , xi, zz) % plot data» xlabel(' x '), ylabel(' y=f(x) '), title(' 2nd and 10th Order curve Fitting ')在下面的图11.2中，原始数据标以'o'，2阶曲线拟合是虚线，10阶拟合是实线。