第一节 MATLAB 中的矩阵的输入 - 湖北民族学院(共31页)

1、第一节 MATLAB 中的矩阵(j zhn)的输入1 直接(zhji)输入一、直接(zhji)在工作窗中输入：A=2, 4, 6, 8;1 3 5 7; 0 0 0 0;1,0,1,0 其意义是定义了矩阵 二、如果矩阵中的元素是等步长的，可以用下面的方法 A=1:0.2:2;1:6;2:2:12A=1:5 “”号在这里表示为转置，而 1:5 中间少了一个循环步长，此时将步长自动取为 1 。2 增删改设已经定义 A=1 2 3 4 5；10 8 6 4 2; B=0 1;1 0; C=1 2;2 4，即已定义A= B= C= 1 2 3 4 5 0 1 1 28 6 4 2 1 0 2 4则命令

2、：A=A(:,1:4);C,B,0 2 0 4 将 A 定义成：A= 而 A(:,3)=； 将删除 A 的第三列 ，得 2 3 4 0 A= 1 2 4 010 8 6 4 2 10 8 4 22 0 1 0 1 2 1 04 1 0 4 2 4 0 43 命令生成使用(shyng) MATLAB 命令生成矩阵(j zhn)一般使用下面的命令 命令(mng lng) linspace，它有两个格式：a1=linspace(1,100) %生成一个从1到100的有100 个元素的向量a2=linspace(0,1)%仍然是有 100 个元素但是是从 0 到 1 的向量 a3=linspace(0

3、,-1) %请与上一个向量进行比较上面是第一种格式 linspace(a,b)，它是将 a 到 b 等分成 100份形成的向量。第二种格式 linspace(a,b,n) 中的 n 为一个正整数，表示是从 a 到 b 等分成 n 份后形成的向量。例如a4=linspace(1,100,11) %从 1 到 100 但只形成11 个元素的向量 a5=linspace(1,100,10) %自己体会这个命令作用 a6=linspace(0,1,11) %加上了“”表示转置a7=linspace(0,-1,10) %自己体会这个命令作用2 命令 ones,zeros 分别形成元素全为 1或全为零的矩

4、阵它也有两种格式。请观察它们的作用： ones(6,3) %生成 63 阶元素全为 1 的矩阵 ones(5) %生成 5 阶元素全为 1 的方阵 zeros(3,6) %生成 36 阶元素全为零的矩阵 zeros(4) %生成四阶元素全为零的方阵 3 命令 diag 生成对角阵及从矩阵的主对角线生成向量，例如：diag(1 3 5 7) %生成了以 1 3 5 7 为主对角线的方阵： ans= 1 0 0 0 0 3 0 0 0 0 5 00 0 7相反如果先定义了一个三阶方阵：A=1 2 3;4 5 6;7 8 9显示(xinsh)： A= 1 2 3 4 5 6 7 8 9则命令(mng

5、 lng) a8=diag(A) 将用 A 的主对角线生成(shn chn)新的列向量：a8= 1 5 9命令 eye(n) 生成 n 阶单位方阵，即主对角线上元素为 1，其余元素为零的方阵。例如键入：A=eye(5) 将得到：A= 1 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 1 00 0 0 0 1 第二节 MATLAB 文件处理 1 文件编辑如果要在 MATLAB 的工作窗定义矩阵，则用鼠标点击屏幕左上方的 File 选择项，再从中选择 New 中的M-file 项并且用鼠标点击它，就打开了 MATLAB 文件编辑窗并且可以在此窗中定义 MATLAB 矩阵了（注

6、意对于已有的文件，可以选择 open 来打开它，然后对其进行修改）。在 MATLAB 文件编辑窗中定义的矩阵与工作窗中定义的方法是完全一样。并且可以在MATLAB 文件编辑窗的菜单中使用菜单命令直接运行。 可以(ky)在MATLAB中使用菜单中的“File”中的“Set path”将当前工作文件夹定义在你正在(zhngzi)工作的文件夹。2 MATLAB工作窗中变量值的保存(bocn)与调用MATLAB 工作窗中的变量在退出 MATLAB 工作状态后值不能保存，如果需要保存，可以使用命令 save 将其存储到磁盘上，命令格式有两种：第一种是用二进制格式来存储。例如先定义三个矩阵：A1=0:3;

7、2*ones(4);4:1:1 ; A2=1 3 2 4;A3=zeros(3,1);生成下列矩阵与向量：键入：save file1 A1 A2 A3 %用二进制格式以文件名 file1.mat 存储 A1,A2,A3save file2.m A1 A3 ascii%用 ascii 码以文件名 file2.m 存储 A1,A3我们还要注意：用二进制格式存储的文件连变量名一起存储并可再重新调入时恢复变量的值，而用 ASCII码存储的文件只存储了变量的值，而变量名是没有的。用二进制格式存储的变量，可用命令 load 调用，调用格式为： load 例如，前面用 save file1 存储了所有变量

8、A1,A2,A3，调用时只要键入 load file1 即可。 第三节 MATLAB 中的矩阵运算 1 矩阵运算命令与通常线性代数(xin xn di sh)命令运算的异同一、MATLAB 在运行时是以矩阵为单位进行运算的。它通常有两种运算，第一种是矩阵运算，运算时满足线性代数中矩阵运算所规定(gudng)一切运算法则，如加、减、乘，乘方即幂运算(yn sun)(当然运算要符合规定的条件，例如矩阵 A 与矩阵 B 相乘，必须 A 的列数等于 B的行数)，运算符号：A+B， A-B， A B (注意“*”不能少) An二、不同之处：1、与通常矩阵运算不同之处在：在线性代数中矩阵不能与数相加减，而

9、在MATLAB的矩阵运算中允许矩阵与数相加减。2、函数命令可以直接作用到矩阵的每一个元素。3、线性代数中矩阵没有除法，而MATLAB中有矩阵除法，例如：输入 A=1:3;4:6;7:9;b=14,32,50;c=Ab显示： c=204 4、函数作用到矩阵的每一个元素，例如如果令 A=1 1/2 1/3; 1/2 1/4 1/8*pi，即定义 则 三、MATLAB中除法运算的规定与意义：1、运算定义：设已经定义好矩阵 A 与矩阵 B，如果矩阵 A 与矩阵B 的行的维数相同，则 MATLAB 中可以用矩阵 A 左除矩阵 B，即可以令： X=AB如果(rgu)矩阵 A 与矩阵(j zhn)B 的列的

10、维数相同(xin tn)，则 MATLAB 中可以用矩阵 A 右除矩阵 B，即可以令： X=B/A 2、矩阵除法的意义 给出线性方程组 AX=B，则 X=A 给出线性方程组的一个解。而 X=B/A 给出了线性方程组 XA=B 的一个解。 目前我们用 MATLAB 求线性方程组的解只有三个方法:当 A 是可逆方阵时，X=inv(A)*B 给出线性方程组AX=B 的唯一解。但是 A 不可逆时方法失效。可用命令 rref 化线性方程组 AX=B 的增广矩阵为行最简阶梯矩阵方法来求解，但线性方程组可能出现无解(称为“超定”)、唯一解(称为“恰定”)及无穷多解(称为“欠定”)的情形。无论 A 是否可逆都

11、可用MATLAB 除法求解，并且无论何种情形都是唯一解。当方程存在唯一解时，三种方法求出的解是一样的。但是用除法作的解一般精度更高。方程为“超定”或者“欠定”时解意义则不同。线性方程组 AX=b 为欠定时有无穷个解，X=Ab得到其中解分量中零元素为最多的一个特解。线性方程组 AX=b 为超定时是无解的, 用 X=Ab得到的是使范数 |AX-b| 为最小的解。我们不详细说明这个范数的意义，可理解为使 AX-b 最接近于零的解。例如方程组 通解为：输入(shr) A=1:3;4:6;7:9;b=14,32,50;c=Ab显示(xinsh)： c=204是其中(qzhng)有一个零分量的特解。输入

12、d=0 32 50;g=Ad显示 g= 1.0e+017 * 0.6305 -1.26100.6305再输入 h=A*g-d显示 h=329614因此 g 不满足 A*g=d，只是使 A*g-d 尽可能接近于零。2 常用的数学运算命令 以下是一些在MATLAB中最常用的数学运算命令(均用小写字母，命令的作用可在MATLAB中用 help 查询其作用与格式)：一、基本函数运算命令 1、三角函数： sin cos tan cot sec csc 2、反三角函数：asin acos atan acot asec acsc atan2(四象限反正切)3、双曲函数(hnsh)： sinh cosh ta

13、nh coth sech csch4、反双曲函数(hnsh)：asinh acosh atanh acoth asech acsch5、复数(fsh)：real imag conj angle6、小数取整：fix(朝零方向)，ceil(朝正无穷方向)，floor(朝负无穷方向)，round(四舍五入)7、对数与指数：log10 log exp8、其它：sqrt abs sign sum prod solve二、线性代数运算命令： 1、det inv rank rref eig cond行列式求值命令 det， 矩阵求秩命令 rank，求方阵的逆方阵命令，inv 求矩阵行最简阶梯阵命令 rref

14、，求特征值与特征向量命令 eig，求矩阵条件数命令 cond三、多项式运算命令 MATLAB中用向量表示多项式，如 a= 1 2 3 0 4 表示多项式：常用的多项式运算命令有：1、多项式加减法：在次低的向量前面加零后使其元素个数相同，再按向量加减法运算命令进行。 2、多项式 a 与多项式 b 相乘：c=conv(a, b)； 3、多项式 a除以多项式 b： q, r=deconv(a, b) q 是商， r 余项) 4、多项式求值： p1=polyval(a, x) (多项式 a 在点 x 处的值) 5、方阵多项式求值： p2=polyvalm(a, A) (矩阵多项式 a 在方阵 A 处的

15、值)6、多项式求根： p3=roots(a) (求多项式 a 的根)7、多项式微分： p4=polyder(a)8、多项式拟合： p5=polyfit(x,y,n)例：x=1 1.2 3 4.2 5; y=x.3-2*x.2+3*x-1; p5=polyfit(x,y,3) 将得出： p5=1.0000 -2.0000 3.0000 -1.0000四、高等数学中的数值运算(yn sun)命令在MATLAB数值计算(j sun)中很少有高等数学中的运算命令，下面只介绍两个：1、求数值积分： trapz2、差分(ch fn)： diffMATLAB中也有很多符号运算命令。但不在本课程中。在那些符号

16、运算命令中能实现许多高等数学中的运算。五、显示格式 有 format rat, format long, format short, format long g, format short g, format long e , format short e , format hex等等。 第四节 MATLAB 的数组操作前面我们谈到了 MATLAB 中有数与矩阵的加减法，这是线性代数中所没有的。MATLAB 中还有许多与高等数学中内容不同的运算。这些运算对于求解实际问题很有作用，它使编程比其它语言要简单方便得多。MATLAB 中的数组运算就是其它计算机语言中所没有的使编程变得十分简单的运算。M

17、ATLAB 中的数组运算只在同阶矩阵中进行。设 A,B 是两个已经定义好的同阶矩阵： 则数组运算 A.*B,A./B,A.B 的值分别是：除了(ch le)前面讲过的数与矩阵的加减法(不算(b sun)数组运算) 外，数与矩阵之间还有下面(xi mian)的数组运算： 从上面叙述可见矩阵与同阶矩阵之间(及数与矩阵之间)的数组运算的定义特点是： 1 定义在两个同阶矩阵或数与矩阵之间进行；2 定义方法是前面一项与后面一项用小数点加运算符连结；3 定义的实质是：当定义的运算在矩阵之间时，是相同位置上的元素进行由运算符规定的运算。当定义的运算是矩阵与数之间，则是矩阵中的每一个元素与数进行运算符规定的运

18、算。对于其他计算机语言(yyn)，这些运算常常要通过双重循环才能完成，而在 MATLAB 中只要(zhyo)简单地用“.”加运算符即可表示(biosh)，这是 MATLAB 特别方便之处。数组运算及函数运算定义容易记忆与使用，它在编程中十分方便，我们看后面的例。如果要画函数 e-xsinx2 +1 在区间 0,3 上的图形，操作过程是：x=0:0.01:3; %定义自变量 x 在区间 0,3 上的值y=exp(x).*sin(x.2)+1 %用数组运算计算函数值 yplot(x,y) %作平面上的曲线图。只用简单三行命令完成了在其他计算机语言中要用很长程序才能完成的程序，这是 MATLAB 语

19、言非常简便好用的原因。 第五节 MATLAB 作图下面(xi mian)介绍 MATLAB 中的作图命令(mng lng)。1 二维作图命令(mng lng) plot一、MATLAB 二维作图命令最常用的是 plot，除作图外，还有如下可控制的操作：1、一张图上画多条曲线(可以使用hold命令)。2、一个屏幕上开多个窗口作图(使用subplot命令)。3、曲线选择线型、颜色(在每条曲线后加,.b 或者 ,color,1,0,0.5 等)。4、标注文字，使用命令：title(字符) 图形标题xlabel(字符) x 轴标注，ylabel(字符) y 轴标注，zlabel(字符) z 轴标注，l

20、egend(字符, 字符,) 可移动标注text(x,y,(z),字符) 在指定坐标标注gtext(字符) 用鼠标在图形中选择地方标注5、坐标网格线显示与关闭命令 grid二、坐标轴控制(也适用于三维作图) axis 1、坐标轴刻度范围控制axis(xmin, xmax, ymin, ymax (, zmin, zmax)。 2、关闭坐标轴显示开关 axis(off) 或者 axis off，axis(on) 或者 axis on 3、坐标轴等刻度命令 axis(equal) 或者 axis equal 4、坐标轴等长命令 axis(square) 或者 axis squal 5、使坐标等刻度

21、与等长命令失效命令 axis(normal) 或者 axis normal三、下面介绍一个在实际应用中非常有用的命令：ginput该命令在三维空间中仍然可用，但效果不如二维空间中那么好。四、背景色控制命令figure(color, r, g, b) 三维空间中此命令仍然(rngrn)适用。五、程序(chngx)标注在程序(chngx)中“%”后的部分用作标注。六、还有许多二维作图命令，请大家参看有关参考书。我们仅就使用二维等值线 contour命令用于求解线性方程组的例子。 MATLAB 二维作图命令 contour可以用于求二维超越方程组的数值解。请看 exmaple1.m:x=0:0.01

22、:2.4995;y=0:0.01:2.4995;A=(0*y+1)*sin(x)+y.2*(0*x+1)+log(5-(0*y+1)*x-y*(0*x+1)-7;contour(x,y,A,0,0)hold onB=3*(0*y+1)*x+2.y*(0*x+1)-(5-(0*y+1)*x-y*(0*x+1).3+1;contour(x,y,B,0,0)hold off此程序用于求方程组：的数值解。此方法的优点是：可以看到图形全貌，可以求出多解情形，可以获极高精度。练习 1：求曲线与曲线在区间0, 上的所有交点。练习 2：求方程组 的数值(shz)解。要求将解代入方程后每个方程左方减去右方后的绝

23、对值小于 0.00001。2 三维作图命令(mng lng) 一、MATLAB 三维作图曲线作图命令(mng lng)最常用的是 plot3，contour3。除作图外，还有如下可控制的操作：三维曲线作图是对参数式为作图。例如螺旋线： 作图命令是 t=0:0.01:10*pi;x=cos(t);y=sin(t);z=t/2; plot3(x,y,z)contour3是三维等值线作图。二、三维曲面作图 常用的三维曲面作图命令有：mesh , surf, meshc, meshz , waterfall, 我们重点介绍 mesh 与 surf。这两个命令都是作曲面图形，但是 mesh 是用网格作图

24、，而 surf 是用曲面片作图。效果有少许不同，但作图编程是相同的(这个作图过程比较复杂)。事实上它们都是用曲面参数式 作的图，因此要画曲面(qmin)图形，首先得会写上面的参数式。然后会用上面的参数式编程。注意在编程时 z 都是二维矩阵(j zhn)，行数与向量 u 的维数相同(xin tn)，列数与 v 的维数相同，x, y 可以是与 z 维数相同的二维矩阵也可以是向量。请看下面马鞍面 的作图过程：y= -1 :0.1: 1; x=y; z=(0*y+1)*x.2- y.2* (0*x+1); mesh(x,y,z);这个命令等价于y= -1 :0.1: 1; x=y; z= ones(s

25、ize(y)*x.2 - y.2*ones(size(x); surf(x,y,z);上面作的图是以 z 的高度作为图形色彩数据的(这是缺省时的数据)。如果加上颜色、光照、纹理、视角等设置，可以得出十分漂亮的图形。请看 example2.m。有时作图的参数式要使用比较复杂的数学运算，例如 2001 年全国数学建模竞赛时要求从血管的 100 张切片中读出血管的中心与半径。竞赛中要求将血管中心线向三个坐标平面投影。但是我们根据血管中心线可画出血管的三维图象。显然这更好看更直观。曲面参数式的推导见广西大学学报2003年第3期，而程序请看 example3.m。练习3：画一个求生圈的图。提示：求生圈的

26、曲面参数方程为： 练习4：画出在二维区域 x2 + y21 上的马鞍面 z= x2- y2。练习5：作曲面 z=x2+y2 与平面 4x+2y+5z=20 相交的图(要求平面只画出被曲面截出的部分，参考右图)。练习6：作圆柱面 x2+z2=1 与圆柱面 y2+z2=1 相交的图(要求平面只画出在每卦限的部分)。3 MATLAB 动画 MATLAB 中的动画功能不是太强，并且需要(xyo)使用 MATLAB 的底层作图命令(mng lng)，即图形的句柄操作。MATLAB 中的作图命令如果使用了赋值语句，即将作图赋给了某个变量，则该变量即为此张图形的句柄。可以对此句柄设置一些图形属性，对于句柄图

27、形我们在前面已经用到了，例如(lr)马鞍面的图形中的赋值语句aa=surf(X,Y,Z);就是将曲面的句柄赋给 aa 了，然后用set(aa, Facecolor, 0,0.7,0, bb, Facecolor, 0.7,0,0.7)设置曲面颜色。MATLAB 动画实际上就是在句柄中设置了图形将是被复盖还是改写等属性，然后在作图时按不同的位置画图(同时擦除前面的图形或者复盖前面的图形以实现动画)。由于作动画一般要涉及编程，我们将此部分内容放到第四节中。在图形的属性控制中要注意区分各种图形对象：figure, axes, line, patch, surface, image, text 等等。

28、4 MATLAB 图象处理 MATLAB 可对图象作处理，常用的简单命令有：imread, image, colormapimread 可以读出图象数据，而image可以将图象数据显示为图象，而colormap可以设置图象色彩。练习7：本文件夹中有一张图“figure0”，是2002年全国数学建模大赛中的一张图，要求从图中读出图形的内切圆圆心与半径。 在MATLAB图象处理中还有光照、色彩、视角等许多处理。由于水平与时间关系，只能简单介绍。第六节 MATLAB 编程一、MATLAB 编程 MATLAB 可以单步(dn b)运行，也可以编程运行。 与所有软件编程一样，MATLAB中也有条件(ti

29、ojin)、分支、循环等语句。条件(tiojin)语句： if end 或者 if else end循环语句： for end 或者 while end需要指出，在其他许多语言中需要循环的地方在MATLAB中是不需要循环的。下面我们看一个动画程序：画单摆的运动程序：pend.m。figure(color,1 1 1)plot(-0.2;0.2,0;0,color,y,linestyle,-,linewidth,10);%画悬线的横梁g=9.8; %重力加速度l=1; %线长theta0=pi/6; %初始角度x0=l*sin(theta0); %初始 x 值y0=-l*cos(theta0);

30、 %初始 y 值axis(-0.75,0.75,-1.25,0); %画坐标范围axis(off); %关闭坐标显示head=line(x0,y0,color,r,linestyle,.,erasemode,xor,markersize,40); %定义小球格式body=line(0;x0,0;y0,color,b,linestyle,-,erasemode,xor); %定义悬线格式t=0; %时间初值dt=0.0005; %时间增量while 1 %死循环(关闭窗口后中止循环) t=t+dt; theta=theta0*cos(sqrt(g/l)*t); %计算时刻 t 时的幅角 x=l*

31、sin(theta); %计算时刻 t 时的 x 坐标 y=-l*cos(theta); %计算时刻 t 时的 y 坐标 set(head,xdata,x,ydata,y); %计算时刻 t 时小球的位置 set(body,xdata,0;x,ydata,0;y); %计算时刻 t 时悬线的位置 drawnow; %重画小球(xio qi)与悬线end二、程序(chngx)与函数将MATLAB运行(ynxng)中的每一句写在扩展名为“m”的文本文件中，就得到MATLAB程序文件。程序中的变量值都在工作窗中。如果是写成function u , v, =fun(x , y, )则是MATLAB函数

32、，函数参数 x , y 可以是变量、向量或者矩阵，返回值 u, v 也可以是变量、向量或者矩阵，存储时函数名必须与存储文件名相同。 例如，运行 redotree(n, B) (设的初值为 B=0 0;0 1，取 n=5 或者 n=6)，可以画出一棵树的图形。而运行leafage1(n)可以画出一张树叶。练习8：作曲面 z=x2+y2 与平面 4x+2y+5z=c，c 从 40 变化到 20 的动画图 第七节 数据拟合MATLAB有数据插值与数据拟合功能，使用的命令是：interpn(注意这里的“n”必须是一个具体正整数，一般为 n=1, 2, 3, ，表示)。例如：运行exmaple51.m，

33、对南半球气旋作可视化显示图。%exmaple51.m-南半球气旋数据，程序实现气旋数据可视化y=5:10:85;x=1:12;z= 2.4, 1.6, 2.4, 3.2, 1.0, 0.5, 0.4, 0.2, 0.5, 0.8, 2.4, 3.6; 18.7, 21.4, 16.2, 9.2, 2.8, 1.7, 1.4, 2.4, 5.8, 9.2, 10.3, 16; 20.8, 18.5, 18.2, 16.6, 12.9, 10.1, 8.3, 11.2, 12.5, 21.1, 23.9 ,25.5;22.1, 20.1, 20.5, 25.1, 29.2, 32.6, 33.0,

34、 31.0, 28.6, 32.0, 28.1, 25.6; 37.3, 28.8, 27.8, 37.2, 40.3, 41.7, 46.2, 39.9, 35.9, 40.3, 38.2, 43.4; 48.2, 36.6, 35.5, 40, 37.6, 35.4, 35, 34.7, 35.7, 39.5, 40, 41.9; 25.6, 24.2, 25.5, 24.6, 21.1, 22.2, 20.2, 21.2, 22.6, 28.5, 25.3, 24.3; 5.3, 5.3, 5.4, 4.9, 4.9, 7.1, 5.3, 7.3, 7, 8.6, 6.3, 6.6; 0

35、.3, 0, 0, 0.3, 0, 0, 0.1, 0.2, 0.3, 0, 0.1, 0.3; xi,yi=meshgrid(1:12,5:85); zi=interp2(x,y,z,xi,yi,cubic); surf(xi,yi,zi)练习9：一天(y tin)24小时每两小时(xiosh)的温度记录如下：9 9 10 18 24 28 27 25 20 18 15 13请用插值方法(fngf)用出一天24小时的温度变化图。用命令 curvefit(fun,x0 , x , y) 可以作曲线拟合，其拟合功能比 MATHMATICA 更强。如某数据近似符合规律y(t)=a+be-kt，已知

36、从100到1000时分别对应于cdata=0.00001*454 499 535 565 590 610 626 639 650 659;求 a, b, k 使 尽可能地小。此问题的解可由运行exmaple52.m得到。附录： 内插及曲线拟合我们常常会有从实验或是对物理现象观察所收集到的数据，分析这些数据需要有辅助的工具。通常(tngchng)这些原始 数据如果会制成图，是成点状分布。我们需要将这些点连成曲线，再进一步将曲线转换成有意义的数学函数 ，才能对数据做比较和分析。上述的过程涉及到内插(interpolation)及曲线拟合(n h)(curve-fitting)，我们(w men)会

37、在以下介绍。1 内插假设一组已知的数据其型态为 f(xk), =1, 2, , n, x1=a, x xn=c, 假设某些点 xi并不属于上述的 xk但是 axic，如果我们要估计这些点的函数值 f(xi) 就须要做内插 (interpolation)。我们可视原数据所描述的函数复杂程度来选择不同的数值内插方法。 1 内插 1.1 HYPERLINK /kejian/matlabcomplex/chapter7/ch7_1_1.htm 一维内插 1.2 HYPERLINK /kejian/matlabcomplex/chapter7/ch7_1_2.htm 二维内插 1.3 HYPERLINK

38、 /kejian/matlabcomplex/chapter7/ch7_1_3.htm Spline 内插 1.1 一维内插线性内插是假设在二个已知数据中的变化为线性关系，因此可由已知二点的座标(a, b)去计算通过这二点的 斜线，公式如下： 其中(qzhng) abc 在上式的 b 点即是代表(dibio)要内插的点，f(b) 则是要计算的内插函数值。下图即是一个(y )以二种内插法的比较 pcxfile12cm,5cmfig9_1.pcx caption线性式与 spline 函数的曲线拟合线性内插是最简单的内插方法，但其适用范围很小；如果原来数据的函数f有极大的变化，假设其数据点之 间为

39、线性变化并不合理。所以我们可以用二次、三次方程组或是另一种称为spline函数来近似原来数据的函 数。MATLAB的一维内插函数是interp1，其语法为interp1(x,y,xi),interp1(x,y,xi,method)；其中的x,y是原已知的 数据的x,y值，而xi则是要内插的数据点，另外method可以设定内插方法有 linear,cubic,spline，分别是一次、三 次方程组和spline函数，其中预设方法是linear。如果数据的变化较大，以 spline函数内插所形成的曲线最平滑 ，所以效果最好。而三次方程组所得到的内插曲线平滑度，则介于线性与spline函数之间。我们

40、以下面的例子说明。假设有一个汽车发动机在定转速下，温度与时间（单位为sec）的三次量测值如下timetemp1 temp2temp3 00 00 120 110176 260 180220 368 240349 477 310450 5110 405503 其中(qzhng)温度的数据从 20oC变化(binhu)到 503oC，如果要估计在t=2.6, 4.9 sec 的温度，可以下列指令(zhlng)计算 x=0 1 2 3 4 5; % 键入时间 y=0 20 60 68 77 110; % 键入第一组时间 y1=interp1(x,y,2.6) % 要内插的数据点为 2.6 y1 =

41、% 对应 2.6 的函数值为 64.8 64.8 y1=interp1(x,y,2.6 4.9) %内插数据点为 2.6, 4.9，注意用 将多个内插点放在其中 y1 = 64.8 106.7 y1=interp1(x,y,2.6,cubic) % 以三次方程组对数据点 2.6 作内插 y1 = % 对应 2.6 的函数值为 66.264 66.264 y1=interp1(x,y,2.6,spline) % 以spline函数对数据点 2.6 作内插 y1 = % 对应 2.6 的函数值为 66.368 66.368以下的例子还配合绘图功能，用以比较不同(b tn)内插方法的差异。 h=1:

42、12; temp=5 8 9 15 25 29 31 30 22 25 27 24; % 这组温度数据(shj)变化较大 plot(h,temp,-,h,temp,+) % 将线性内插结果(ji gu)绘图 h_3=1:0.1:12 % 要每0.1小时估计一次温度值 t_3=interp1(h,temp,h_3,cubic) % 以三次方程组做内插 t_s=interp1(h,temp,h_3,spline) % 以spline函数做内插 hold on subplot(1,2,1) plot(h,temp,-,h,temp,+,h_3,t_3) % 将线性及三次方程组内插绘图 subplot

43、(1,2,2) plot(h,temp,-,h,temp,+,h_3,t_s) % 将线性方程组及spline内插绘图 hold off7.1.2 二维内插二维内插与一维内插的区别(qbi)是二维内插数据为二维,语法结构为interp2(X,Y,Z,XI,YI),其中X,Y,Z为已知数据,Z=Z(X,Y),而XI,YI 为要插值的数据(shj)点;如果语法结构为interp2(X,Y,Z,XI,YI,method),其中method可以为linear,cubic表示线形或三次方插值,我们以下例说明:假设(jish)一汽车的转速(单位为:rpm)、温度（单位为：oC）、时间（单位为：sec）如下

44、表：time speed 0 2000 rpm 3000 rpm 4000 rpm 1 20 110 176 2 60 180 220 3 68 240 349 4 77 310 450 5 110 405 503 其中温度的数据为20oC到503oC，如果要估计t=2.6, sec, rpm=2500的温度，可以利用下面的语句：d2(:,1)=0 1 2 3 4 5; % 将时间输入 d2(:,2)=2000 20 60 68 77 110; % 将 rpm=2000的温度输入 d2(:,3)=3000 110 180 240 310 405; % 将 rpm=3000 的温度输入 d2(:

45、,4)=4000 176 220 349 450 503; % 将 rpm=4000 的温度输入 t=d2(2:6,1); %选择做内插的时间rpm=d2(1,2:4); % 选择做内插的 rpm temp=d2(2:6,2:4); % 选择做内插的温度temp_i=interp2(rpm,t,temp,2500,2.6) % 以线形内插决定 rpm=2500,t=2.6 的温度 temp_i = 140.4000 1.3 Spline 内插关于spline内插我们(w men)在1.1节已介绍过，它可以用interp1指定内插方式为spline来做。另一种方式也可以用 spline(x,y,

46、xi)来做，其中的x,y,xi的用法与interp1中的语法相同。事实上这二种方法采用相同的spline 函数做 运算，也就是当我们执行interp1(x,y,xi,spline)时，MATLAB即呼叫spline(x,y,xi)做运算，再将计算结果传回 interp1。 我们还是以在 1.1 节相同的数据(shj)做 spline 内插示范， x=0 1 2 3 4 5; y=0 20 60 68 77 110; y1=spline(x,y,2.6) y1 = 67.3 y1=spline(x,y,2.6,4.9) y1 = 67.3 105.23 曲线(qxin)拟合曲线拟合 (curve

47、-fitting) 故名思义就是要将一组离散的数据以一个近似的曲线方程组来代表。有了这个解析形 态的方程组，我们就可以很方便去运用它。曲线拟合与前述的内插有许多相似之处，但是这二者最大的区 别在于曲线拟合要找出一个曲线方程组而内插仅是要求出内插数值即可。 要找出近似一组数据（也就是最能代表(dibio)这些数据）的曲线方程组，有许多选择(xunz)，从简单的一阶线性方程组 到高阶的多项式都有。我们以下就介绍(jisho)二种方法：线性回归(linear regression)，多项式回归(polynomial regression)。3 曲线拟合 3.1 HYPERLINK /kejian/m

48、atlabcomplex/chapter7/ch7_3_1.htm 线性回归 3.2 HYPERLINK /kejian/matlabcomplex/chapter7/ch7_3_2.htm 多项式回归 3.3 HYPERLINK /kejian/matlabcomplex/chapter7/ch7_3_3.htm 多项式拟合及函数计算 3.1 线性回归我们以一简单数据组来说明什么是线性回归。假设有一组数据型态为 y=y(x)，其中 x=0, 1, 2, 3, 4, 5, y=0, 20, 60, 68, 77, 110 如果我们要以一个最简单的方程组来近似这组数据，则非一阶的线性方程组莫属。

49、先将这组数据绘图如下 plot(x, y)图中的斜线是我们随意假设一阶线性方程组 y=20 x，用以代表这些数据的一个方程组。以下将上述绘图的 MATLAB 指令列出，并计算这个线性方程组的 y 值与原数据 y 值间误差平方的总合。 x=0 1 2 3 4 5; y=0 20 60 68 77 110; y1=20*x; % 一阶线性方程组的 y1 值 sum_sq = sum(y-y1).2); % 误差平方总合为 573 axis(-1,6,-20,120) plot(x,y1,x,y,o), title(Linear estimate), grid如此(rc)任意的假设一个线性方程组并无

50、根据(gnj)，如果换成其它人来设定就可能采用不同的线性方程组；所以我们 须要有比较精确(jngqu)方式决定理想的线性方程组。我们可以要求误差平方的总合为最小，做为决定理想的线性方 程序的准则，这样的方法就称为最小平方误差(least squares error)或是线性回归。MATLAB的polyfit函数提供了 从一阶到高阶多项式的回归法，其语法为polyfit(x,y,n)，其中x,y为输入数据组n为多项式的阶数，n=1就是一阶 的线性回归法。polyfit函数所建立的多项式可以写成 从polyfit函数得到的输出值就是上述的各项系数，以一阶线性回归为例n=1，所以只有 二个输出值。如

51、果指令为coef=polyfit(x,y,n)，则coef(1)= , coef(2)= ,.,coef(n+1)= 。注意上式对n 阶的多 项式会有 n+1 项的系数。我们来看以下的线性回归的示范： x=0 1 2 3 4 5; y=0 20 60 68 77 110; coef=polyfit(x,y,1); % coef 代表线性回归的二个输出值 a0=coef(1); a1=coef(2); ybest=a1*x+a0; % 由线性回归产生的一阶方程组 sum_sq=sum(y-ybest).2); % 误差(wch)平方总合为 356.82 axis(-1,6,-20,120) plot(x,ybest,x,y,o), ti