1、MATLAB程序设计,目录,1 矩阵及其运算 2 语句与函数 3 命令与窗口环境 4 图形功能 5 程序设计,1 矩阵及其运算,1.1 矩阵的建立 1.2 矩阵的裁减与拼接 1.3 矩阵的基本运算 1.4 矩阵的特殊运算,1.1 矩阵的建立,1直接输入法 规则： 矩阵元素必须用 括住 矩阵元素必须用逗号或空格分隔 在 内矩阵的行与行之间必须 用分号分隔,例1-1：直接输入 a=1,2,3;4,5,6;7,8,9;得到3*3矩阵a：,例如 a=1 2 0;3 0 5;7 8 9 a =1 2 03 0 57 8 9 a(3,3)=0 a =1 2 03 0 57 8 0,2 利用矩阵函数生成特殊
2、矩阵,常用的产生特殊矩阵的函数有： zeros：产生全0矩阵(零矩阵)。 ones：产生全1矩阵(幺矩阵)。 eye：产生单位矩阵。 rand：产生01间均匀分布的随机矩阵。 randn：产生均值为0，方差为1的标准正态分布随机矩阵。,3. 利用M文件建立矩阵,对于比较大且比较复杂的矩阵，可以为它专门建立一个M文件。让后利用load 命令导入。过程如下：(1) 启动有关编辑程序或MATLAB文本编辑器，并输入待建矩阵： (2) 把输入的内容以纯文本方式存盘(设文件名为mymatrix.txt)。 (3) 在MATLAB命令窗口中输入 load mymatrix.txt,4行向量的特殊建立方式,
3、用冒号生成，一般格式是： e1:e2:e3 其中e1为初始值，e2为步长，e3为终止值。注： e2=1时可以省略。 用linspace或logspace函数产生行向量。 调用格式为：linspace(a,b,n) 作用：生成a到b共n个数值的等差数值。 用logspace函数产生行向量。 调用格式为：logspace(a,b,n) 作用：生成10a到10b共n个数值的等比数值。,1.2 矩阵的裁减与拼接,从已知矩阵中取出若干行若干列可构成新矩阵。若干个矩阵组合在一起可组成新矩阵。 1 矩阵元素的引用： 通过下标引用，如a(i,j)表示引用矩阵a的第i行第j列元素。 通过元素序号引用，如a(i)
4、表示引用矩阵a的第i个元素。在MATLAB中，矩阵元素按列存储，先第一列，再第二列，依次类推。 若A=1,2,3;4,5,6; 则A(3)=2,2矩阵拆分 (1) 利用冒号表达式获得子矩阵 A(:,j)表示取A矩阵的第j列全部元素；A(i,:)表示A矩阵第i行的全部元素；A(i,j)表示取A矩阵第i行、第j列的元素。 A(i:i+m,:)表示取A矩阵第ii+m行的全部元素；A(:,k:k+m)表示取A矩阵第kk+m列的全部元素，A(i:i+m,k:k+m)表示取A矩阵第ii+m行内，并在第kk+m列中的所有元素。,(2)通过输入行号或列号直接提取矩阵的若干行或列。如A(1,3,2,3)表示提取
5、矩阵A的第1，3行和2，3列组成新矩阵。(3) 利用空矩阵删除矩阵的元素 在MATLAB中，定义为空矩阵。给变量X赋空矩阵的语句为X=。注意，X=与clear X不同，clear是将X从工作空间中删除，而空矩阵则存在于工作空间中，只是维数为0。,例2-3 分别建立33、32的零矩阵和与矩阵A同样大小的零矩阵。 (1) 建立一个33零矩阵。 zeros(3) (2) 建立一个32零矩阵。 zeros(3,2) (3) 设A为23矩阵，则可以用zeros(size(A)建立一个与矩阵A同样大小零矩阵。 A=1 2 3;4 5 6; %产生一个23阶矩阵A zeros(size(A) %产生一个与矩
6、阵A同样大小的零矩阵,例2-4 建立随机矩阵： (1) 在区间20,50内均匀分布的5阶随机矩阵。 (2) 均值为0.6、方差为0.1的5阶正态分布随机矩阵。 命令如下： x=20+(50-20)*rand(5) y=0.6+sqrt(0.1)*randn(5) 此外，常用的函数还有reshape(A,m,n)，它在矩阵总元素保持不变的前提下，将矩阵A重新排成mn的二维矩阵。,1.3 矩阵的基本运算,1基本算术运算 MATLAB的基本算术运算有：(加)、(减)、*(乘)、/(右除)、(左除)、(乘方) 、(转置) 。 注意，运算是在矩阵意义下进行的，单个数据的算术运算只是一种特例。 2点运算
7、在MATLAB中，有一种特殊的运算，因为其运算符是在有关算术运算符前面加点，所以叫点运算。点运算符有.*点乘、./点右除、.点左除和.点乘方。两矩阵进行点运算是指它们的对应元素进行相关运算，要求两矩阵的维数相同。,1 矩阵加减运算,假定有两个矩阵A和B，则可以由A+B和A-B实现矩阵的加减运算。运算规则是：若A和B矩阵的维数相同，则可以执行矩阵的加减运算，A和B矩阵的相应元素相加减。如果A与B的维数不相同，则MATLAB将给出错误信息，提示用户两个矩阵的维数不匹配。,2. 矩阵乘（）运算 规则： A矩阵的列数必须等于B矩阵的行数 标量可与任何矩阵相乘。 假定有两个矩阵A和B，若A为mn矩阵，B
8、为np矩阵，则C=A*B为mp矩阵。 a=1 2 3;4 5 6;7 8 0;b=1;2;3;c=a*b c =143223,3. 矩阵除运算,左除( ): x=ab 相当于解方程：a*x=b 右除(/ ): x=b/a相当于解方程： x*a=b 条件：对左除，a的行数必与b的行数相同对右除，a的列数必与b的列数相同 举例：a=1 2 3;4 2 6;7 4 9;b=4;1;2; c=5,7,9;x=ab ,y=c/a,注：对于含有标量的运算，两种除法运算的结果相同，如3/4和43有相同的值，都等于0.75。又如，设a=10.5,25，则a/5=5a=2.1000 5.0000。,a p a
9、自乘p次幂,方阵,1的整数,4. 矩阵乘方ap,对于p的其它值,计算将涉及特征值 和特征向量，如果p是矩阵，a是标量 ap使用特征值和特征向量自乘到p次 幂；如a,p都是矩阵，ap则无意义。,a=1,2,3;4,5,6;7,8,9;a2ans =30 36 4266 81 96102 126 150,2 语句和函数,赋值语句一般形式：变量名=表达式与C语言不同，此处变量不需先定义。注意：matlab严格区分大小写字母，因此a与A是两个不同的变量 命令语句、函数语句等由命令或函数组成语句。完成命令或函数约定的功能。,标量函数三角函数：sin、cos、tan、cot、sec等其他基本函数：sqrt
10、、pow2、exp、log、log10、abs等 向量函数max、min、sum、length、mean、median、sort、prod等 矩阵函数 zeros、ones、eye、rand、randnsize、det、rank、inv、eig、trace、poly等,3 命令与窗口环境,主窗口 菜单栏、工具栏 命令窗口 工作区窗口、当前目录窗口 命令历史窗口,3.1在线帮助系统 利用help菜单 直接在命令窗口输入help 命令名（或函数名）将显示对应命令或函数的帮助信息。,3.2 数据的显示格式 可以通过菜单file下的子菜单preference来选择显示格式或直接通过format命令来修
11、改。 常用的显示格式及意义如下： format (short)：短格式（5位定点数）99.1253 format long：长格式（15位定点数 99.12345678900000 format short e：短格式e方式 9.9123e+001 format long e：长格式e方式 9.912345678900000e+001 format bank：2位十进制 99.12 format hex：十六进制格式,3.3命令行编辑器类似于一般的文本编辑器 3.4工作区命令 who：显示当前工作空间中所有变量的一个简单列表 whos：则列出变量的大小、数据格式等详细信息 clear ：清除工
12、作空间中所有的变量 clear 变量名：清除指定的变量 save filename variables将变量列表variables所列出的变量保存到磁盘文件filename中 load filename variables将以前用save命令保存的变量variables从磁盘文件中调入MATLAB工作空间 其他：diary、path、cd 、ls、delete等,4 图形功能,4.1 二维绘图 基本绘图函数：plot 对数坐标系绘图函数：semilogx、 semilogy、loglog 极坐标系绘图函数：polar 双纵坐标绘图函数：plotyy 条形图绘图函数：bar 饼状图绘图函数：pi
13、e 4.2 三维绘图 Mesh、meshgrid、contour、quiver等,基本绘图函数：plot,格式：1、PLOT(X,Y)2、PLOT(X,Y,S)3、PLOT(X1,Y1,S1,X2,Y2,S2,X3,Y3,S3,.) 说明：1、第一个变量X作为横坐标，第二个变量Y作为纵坐标,以缺省的方式绘图。2、以S指定的方式绘图，其他同13、同时绘制多条曲线，Xi,Yi,Si的意义同上,PLOT(X,Y),情形1：若X、Y全为向量，则以X中值作为横坐标，Y 中值作为纵坐标。要求 X、Y必有相同的结构。 例：x=0:0.01*pi:2*pi;y=sin(x);plot(x,y) 将在0,2*p
14、i区间绘制sin函数曲线,PLOT(X,Y),情形2：当X和Y为同阶矩阵时，按照矩阵的列操作，在同一幅图中绘出多条曲线（每列对应一组数据）。 例： x=0:0.01*pi:2*pi; y=sin(x),cos(x); plot(x,x,y); %将绘同时绘制sin与cos函数曲线,plot(X,Y,s),说明，依据s指定的格式绘图s参数：- 实线 . 点 * 星号 上三角： 点线 o 圆 s 方形 右三角- 虚线 + +号 v 下三角 p 正五边形y 黄色 m 紫红色 c 蓝绿色 r 红色 g 绿色 b 蓝色 w 白色 k 黑色上述参数的不同组合可以为图形设置不同的线形、颜色和标识，调用时以单
15、引号来引用，各选项直接相连，不需要分隔符。 举例：plot(x,y ,*r) %以红色的*来绘制图像,PLOT(X1,Y1,S1,X2,Y2,S2,X3,Y3,S3,.),同时绘制多条曲线， Xi,Yi,Si的意义与前面相同 举例：x=0:0.01*pi:2*pi;y1=sin(x);y2=cos(x);y3=sin(x).*cos(x);plot(x,y1,-.r,x,y2,-*g,x,y3,b)%分别以红色、绿色、蓝色绘制sin、cos、sin*cos曲线,对数坐标系绘图函数：,Semilogx：横坐标为对数形式semilogy ：纵坐标为对数形式loglog ：横纵坐标均为为对数形式 这
16、些函数的调用格式与plot相同，只是坐标系为对数形式。 例1：x=1:100;y=exp(x);plot(x,y) %以缺省坐标系绘图semilogy(x,y) %纵坐标为对数形式，曲线为直线,对数坐标系绘图函数：,例2：若y=2*x3 ，则有log(y)=log(2)+b*log(x)x=1:50;y=2.*x.3;plot(x,y) %直接绘图，图像不容易看出变量间关系，若用双对数坐标系绘图，关系明显，曲线为直线loglog(x,y) %双对数坐标系形式,极坐标系绘图函数：polar,格式：POLAR(THETA, RHO) 说明：THETA为以弧度表示的角度向量，RHO为对应的幅向量（相
17、当于半径） 例：若x2+y2=9,用极坐标绘图 theta=0:0.1*pi:2*pi; r=3*ones(1,21); polar(theta,r,r) %以红色绘制极坐标图像,为一圆形曲线 若r=1+0.2*theta ,则是一发散的曲线,条形图绘图函数：bar,格式：BAR(X,Y) 例：x=0:0.1*pi:2*pi;y=sin(x);bar(x,y) 将在0,2*pi区间以条形图形式绘制sin函数图像,饼状图绘图函数：pie,格式： PIE(X) 或PIE(X,label) 或PIE(X,explode,label) 说明：X为一非负向量，通过计算X/SUM(X)来确定X向量每一部分
18、所占的比重，然后按照每一部分的比重绘制饼状图。explode为与X结构相同的非负向量，若某值为一非零值，则该部分从图中分离。 例：如某班级，优、良、中、及格、不及格人数分别为5、10、20、15、3，可用以下方式绘制饼状图 pie(5 10 20 15 3,优,良,中,及格 , 不及格 )e=0 1 0 0 0; pie(5 10 20 15 3,优,良,中,及格 , 不及格 ),双纵坐标绘图函数：plotyy,格式： 1、PLOTYY(X1,Y1,X2,Y2) 2、PLOTYY(X1,Y1,X2,Y2,FUN) 3、PLOTYY(X1,Y1,X2,Y2,FUN1,FUN2) 说明： 1、坐标
19、轴左边标识X1,Y1绘图的纵坐标，右边标识X2,Y2绘图的纵坐标 2、用函数FUN来代替plot绘图，其他同1 3、分别用函数FUN1、FUN2来绘制X1,Y1及X2,Y2对应的图。,PLOTYY(X1,Y1,X2,Y2),例1：x=0:0.1*pi:2*pi;plotyy(x,0.1*sin(x),x,cos(x),PLOTYY(X1,Y1,X2,Y2,FUN),例2：若y=a*ebx ，则有log(y)=log(a)+bxx=1:0.1:8; y=2.*exp(3*x); z= 100.*exp(-2*x); plotyy(x,y,x,z,semilogy),PLOTYY(X1,Y1,X2
20、,Y2,FUN1,FUN2),例3：x =0:0.1*pi:2*pi; y=sin(x); z=exp(x); plotyy(x,y,x,z,plot,semilogy),5 Matlab 程序设计,5.1 关系运算符、=、=、=、= 5.2逻辑运算符&与、|或、非、xor异或 5.3条件和循环语句 If语句 switch语句 for 语句 while语句 break和continue语句,MATLAB的两种工作方式,交互式的命令行工作方式：在命令窗口输入matlab命令，直接执行并观察其执 行结果 M文件的程序工作方式先把命令编辑成M文件，然后执行文件。,5.3 M文件,用MATLAB语言编
21、写的程序，称为M文件。M文件有两类：脚本文件和函数文件。 脚本文件：是由matlab命令组成的集合，没有输入参数，也不返回输出参数，可视为批处理文件。 函数文件：可以输入参数，也可返回输出参数。,1 M文件的建立与编辑,建立新的M文件：从MATLAB命令窗口的File菜单中选择New菜单项，再选择M-file命令。,编辑已有的M文件：从MATLAB命令窗口的Flie菜单中选择Open M-file命令。,2 脚本文件,将需要运行的命令编辑到一个脚本文件中，然后在MATLAB命令窗口输入该脚本文件的名字，就会顺序执行命令文件中的命令。 【例1】 建立一个命令文件将变量a,b的值互换。,e1m文件
22、： a=1:9; b=11,12,13;14,15,16;17,18,19; c=a;a=b;b=c;ab 在MATLAB的命令窗口中输入e1，将会执行该命令文件。点击M文件编辑窗口相应的按钮或命令也可执行文件。,函数文件是另一种形式的M文件，每一个函数文件都定义一个函数。事实上，MATLAB提供的标准函数大部分都是由函数文件定义的。,3 函数文件,3.1 函数文件格式,函数文件由function语句引导，其格式为：function 输出形参表=函数名(输入形参表)注释说明部分函数体 注：其中函数名的命名规则与变量名相同。输入形参为函数的输入参数，输出形参为函数的输出参数。当输出形参多于1个时
23、，则应该用方括号括起来。,【例2】 编写函数文件求小于任意自然数n的Fibonacci数列各项。,function f=ffib(n)%用于求Fibonacci数列的函数文件%f=ffib(n) f=1,1;i=1;while f(i)+f(i+1)nf(i+2)=f(i)+f(i+1);i=i+1;end,将以上函数文件以文件名ffib.m存盘，然后在MATLAB命令窗口输入以下命令，可求小于2 000的Fibonacci数。ffib(2000),函数文件举例,3.2 函数调用,函数文件编制好后，就可调用函数进行计算了。如上面定义ffib函数后，调用它求小于2000的Fibonacci数。
24、函数调用的一般格式是：输出实参表=函数名(输入实参表) 直接以这种格式：函数名(输入实参表)此时函数值赋给缺省变量ans。,多项式表示：把多项式表达成一个行向量， 该向量中的元素是按多项式降幂排列的。f(x)=anxn+an-1xn-1+a1x+a0可用行向量 p=an an-1 a1 a0表示,五、 多项式表示及运算,Matlab中常用的多项式函数 roots：多项式求根 poly：由根来创建多项式 polyval：多项式求值 polyfit：多项式拟合,5.1 roots 求多项式的根,p = 1.00 -3 2;%对应多项式p(x)=x2-3x2 r=roots(p) 则有： r =21
25、 p2=poly(r) %可发现p2等于p 注：多项式系数向量用行向量表示，一组根用列向量表示。,格式：r=roots(p),其中p为一行向量，其作用是求系数为行向量p的多项式p(x)=0的根,5.2 poly 生成多项式系数向量,格式：p=poly(A) 若A为nn矩阵，则生成矩阵A的特征多项式 若A为向量，则生成根为A中所列元素的多项式。 如：a=1 2 3;4 5 6;7 8 0; p=poly(a) %求矩阵a的特征多项式 则：p = 1.00 -6.00 -72.00 -27.00 对应多项式：p(x)=x3-6x2-72x-27 r1=root(p) %求方程p(x)=0的根。 r
26、2=eig(a) %求矩阵a的特征根,计算结果r1,r2相同,5.3 polyval：多项式求值,格式：Y = POLYVAL(P,X) 说明：计算自变量取值为X时，多项式P对应的值，若P为N+1维向量，将计算 Y = P(1)*X.N + P(2)*X.(N-1) + . + P(N)*X + P(N+1) 如：p=3,2,5,7;x=0,1,2;y=polyval(p,x) 则y=7 17 49,5.4 polyfit：多项式拟合,格式：P= POLYFIT(X,Y,N) 说明：X、Y是相同结构的向量，通过最小二乘法，寻求拟合程度最好的N次多项式（使得Y与X代入多项式的计算值间的误差平方和
27、最小）,并将多项式系数作为函数返回值赋给P。 如：X=1:11;Y=1,5,7,10,15,20,30,40,60,90,150;P=POLYFIT(X,Y,3) 有：P 0.368 -4.45 19.1 -6.76 把X代入到多项式P求对应多项式值Y1 Y1=POLYVAL(P,X) 有：Y1= -1.72,6.62,10.48,12, 13.56 ,17.2,25.2, 39.7, 63, 97.3, 144.7 比较Y与Y1，有很大偏差，适当调整多项式的次数，可使偏差减小，如用6次以上多项式拟合，仅有很小偏差。通过绘图函数polt（X,Y,r,X,Y1)可更直观的考察其偏差。,七、代数方
28、程组求解,matlab中有两种除运算左除和右除。 对于方程ax=b，a 为anm矩阵，有三种情 况： 当n=m时，此方程成为“恰定”方程 当nm时，此方程成为“超定”方程 当nm时，此方程成为“欠定”方程matlab定义的除运算可以很方便地解上 述三种方程,1.恰定方程组的解,方程ax+b(a为非奇异)x=a-1 b 矩阵逆 两种解: x=inv(a)b 采用求逆运算解方程 x=ab 采用左除运算解方程,方程ax=b a=1 2;2 3;b=8;13; x=inv(a)*b x=abx = x = 2.00 2.003.00 3.00,=,a x = b,例: x1+2x2=8 2x1+3x2
29、=13,2.超定方程组的解,方程 ax=b ,mn时此时不存在唯一解。 方程解 (a a)x=a b x=(a a)-1 a b 求逆法 x=ab matlab用最小二乘法找一个准确地基本解。,例: x1+2x2=1 2x1+3x2=23x1+4x2=3a=1 2;2 3;3 4;b=1;2;3;解1 x=ab 解2 x=inv(aa) a b x = x =1.00 1.000 0.00,=,a x = b,3.欠定方程组的解,当方程数少于未知量个数时,即不定 情况,有无穷多个解存在。 matlab可求出两个解： 用除法求的解x是具有最多零元素的解 是具有最小长度或范数的解，这个解是基于伪逆
30、pinv求得的。,x1+2x2+3x3=12x1+3x2+4x3=2a=1 2 3;2 3 4;b=1;2;x=ab x=pinv(a)b x = x =1.00 0.830 0.330 -0.17,a x = b,六、微分方程求解,微分方程求解的仿真算法有多种，常用的有Euler（欧拉法）、Runge Kutta(龙格-库塔法。 Euler法称一步法，用于一阶微分方程,当给定仿真步长时：所以 yn+1 = yn + hf (xn,yn) n=0,1,2y(x0)=y0,Runge Kutta法龙格-库塔法：实际上取两点斜率 的平均 斜率来计算的，其精度高 于欧拉算法 。 龙格-库塔法：ode
31、23 ode45,k1=hf(xn,yn) k2=hf(xn+h,yn+k),例：x+(x2-1)x+x=0 为方便令x1=x，x2=x分别对x1,x2求一 阶导数，整理后写成一阶微分方程组 形式x1=x2x2=x2(1-x12)-x1 建立m文件 解微分方程,建立m文件 function xdot=wf(t,x) xdot=zeros(2,1) xdot(1)=x(2) xdot(2)=x(2)*(1-x(1)2)-x(1) 给定区间、初始值;求解微分方程 t0=0; tf=20; x0=0 0.25; t,x=ode23(wf, t0, tf, x0) plot(t,x), figure(
32、2),plot(x(:,1),x(:,2),例1：f(x)=x2+3x+2在-5 5区间的最小值 f=fmin(x2+3*x+2,-5,5) 例2：f(x)=100(x2-x12)2+(a-x1)2在x1=a, x2=a2处有最小值 function f=xun(x,a) f=100*(x(2)-x(1).2).2+(a-x(1).2; x=fmins(xun,0,0, , ,sqrt(2),八、数据分析与插值函数,max 各列最大值 mean 各列平均值 sum 各列求和 std 各列标准差 var 各列方差 sort 各列递增排序,九、拟合与插值,1. 多项式拟合 x0=0:0.1:1；
33、y0=-.447 1.978 3.11 5.25 5.02 4.66 4.01 4.58 3.45 5.35 9.22; p=polyfit(x0,y0,3) p = 56.6915 -87.1174 40.0070 -0.9043 xx=0:0.01:1;yy=polyval(p,xx); plot(xx,yy,-b,x0,y0,or),2.插值 插值的定义是对某些集合给定的数据点之间函数的估值方法。 当不能很快地求出所需中间点的函数时，插值是一个非常有价值的工具。 Matlab提供了一维、二维、 三次样条等许多插值选择,table1 table2 intep1 interp2 spline 利用已知点确定未知点 粗糙 精确 集合大的 简化的,插值函数,
