Matlab_2数值2011.ppt

第二讲MATLAB的数值计算,MatLab具有出色的数值计算能力，占据世界上数值计算软件的主导地位,数值运算的功能,创建矩阵矩阵运算多项式运算线性方程组数值统计线性插值函数优化微分方程的数值解,MATLAB矩阵与数组,一、命令行的基本操作,直接输入法规则：1)矩阵元素必须用括住2)矩阵元素必须用逗号或空格分隔3)在内矩阵的行与行之间必须用分号分隔,x=2pi/2;sqrt(3)3+5i矩阵元素可以是任何MatLab表达式，可以是实数，也可以是复数，复数可用特殊函数i，j输入,1.创建矩阵的方法,a=123456,例:a=123;456,符号的作用,逗号和分号的作用逗号和分号可作为指令间的分隔符，MatLab允许多条语句在同一行出现。分号如果出现在指令后，屏幕上将不显示结果。,冒号的作用用于生成等间隔的向量，默认间隔为1。用于选出矩阵指定行、列及元素。循环语句,例：,a=1：5a=1,2,3,4,5,b=6：-2：2b=642,x=10,y=20,z=30;,向量,步长1,步长-2,矩阵,M=1；2；3M=123,N=1，2，3N=123,行向量(1x3阶阵),列向量(3x1阶阵),又例：,A=1234；5678；9101112,A=123456789101112,10,56910,56789101112,1256910,A(3,2),A(2：3，1：2),A(2：3，：),A(：，1：2),forj=1：3A(j，1)=10;end,102341067810101112,：所有行(列),1.只要是赋过值的变量，不管是否在屏幕上显示过，都存储在工作空间中，以后可随时显示或调用。变量名尽可能不要重复，否则会覆盖。,4.当一个指令或矩阵太长时，可用续行,注意：,2.MatLab严格区分大小写字母，因此a与A是两个不同的变量。,3.MatLab函数名必须小写。,2.用MatLab函数创建矩阵,空阵允许输入空阵，当一项操作无结果时，返回空阵。rand随机矩阵eye单位矩阵zeros全部元素都为0的矩阵ones全部元素都为1的矩阵,详见书上P.28.表2-2能生成矩阵的函数,3.矩阵的修改,直接修改可用键找到所要修改的矩阵，用键移动到要修改的矩阵元素上即可修改。指令修改可以用A(,)=来修改。,例如a=120;305;789a=120305789a(3,3)=0a=120305780,还可以用函数subs,find函数修改。,把MatLab工作空间中一些有用的数据长久保存下来的方法是生成mat数据文件。save将工作空间中所有的变量存到默认文件名.mat文件中。,二、数据的保存与获取,savedata将工作空间中所有的变量存到data.mat文件中。savedataab将工作空间中a和b变量存到data.mat文件中。,loadloaddataloaddataabmat文件是标准的二进制文件，还可以ASCII码形式保存。,即可恢复保存过的所有变量,可用load指令调用已生成的mat文件。,%-ssfile.m-a=123456;789101112;131415161718;b=123456;,%-ssfile2.m-c=a*b;c,ssfilessfile2,%-ssfile.m-a=123456;789101112;131415161718;b=123456;savedataab,%-ssfile2.m-loaddataabc=a*b;c,ssfileclearssfile2,data.mat外存,内存,MatLab的数据传递,%-ssfile3.m-n=21；x=zeros(1,n);y=zeros(1,n);fori=1:nx(i)=(i-1)*pi/10;y(i)=sin(x(i);endsavedataxy,%-ssfile4.m-clearloaddataxyplot(x,y);,实现MatLab模块之间的数据传递方式之一,%-ssfile5.m-clearloaddataxz1=sin(x);z2=cos(x);plot3(x,z1,z2,m:p)gridon,MatLab变量的查询与清除,变量查询函数who与whos作用都是列出在matlab工作空间中已经驻留的变量名清单不同的是whos在给出驻留变量的同时，还给出他们的维数及性质清除matlab工作空间中所有变量，使用clear,1.矩阵加、减（,）运算规则：相加、减的两矩阵必须有相同的行和列两矩阵对应元素相加减。允许参与运算的两矩阵之一是标量。标量与矩阵的所有元素分别进行加减操作。,三、矩阵运算,2.矩阵乘（）运算规则：A矩阵的列数必须等于B矩阵的行数标量可与任何矩阵相乘。,例:a=123;456;780;b=123;c=a*bc=143223,d=-1;0;2;f=pi*df=-3.141606.2832,正确否？,b=1;2;3;c=a*bb=123;c=a*bb=1,2,3;c=a*bb=1,2,3;c=a*b,错非列向量,对,注意行向量和列向量的表示,inv矩阵求逆det行列式的值eig矩阵的特征值diag对角矩阵矩阵转置sqrt矩阵开方,4.矩阵的其它运算,详见书上P.29.表2-3常用的矩阵函数,A=17241815235714164613202210121921311182529,d=det(A),t=trace(A),rk=rank(A),c=cond(A),d=5070000t=65rk=5c=5.4618,A=magic(5),例2.8,魔方矩阵,行列式值迹秩条件数,5.矩阵的一些特殊操作,矩阵的变维a=1:12;b=reshape(a,3,4)c=zeros(3,4);c(:)=a(:)矩阵的变向rot90:旋转;fliplr:上翻;flipud:下翻矩阵的抽取diag:抽取主对角线;tril:抽取主下三角;triu:抽取主上三角矩阵的扩展,a=123456789101112,b=123456789101112,关系运算,1).数组加减(.+,.-)可不用a.+ba.-b,6.矩阵的数组运算(点运算),对应元素相加减（与矩阵加减等效）,03-11-13,数组运算指元素对元素的算术运算，与线性代数中矩阵运算不同,aba，b两数组必须有相同的行和列,两数组对应元素相乘。,a./b=b.a-都是b的元素被a的对应元素除a.b=b./a-都是a的元素被b的对应元素除,2).数组乘除(，./，.）,./右除1./2,例,矩阵乘和数组乘的不同,A=1112131415;1617181920;,B=12345;678910;,C=205530280730,D=112439567596119144171200,矩阵规则相乘,对应元素相乘,若a=123;b=456;问:,例,数组除和矩阵除的不同,c1=4/15/26/3=4.00002.50002.0000,c1=a.b=?,c3=b./a=?,c2=b.a=?,c2=1/42/53/6=0.25000.40000.5000,c3=c1,c4=a./b=?,c4=c2,c5=ba=?,c6=b/a=?,要求b是非奇异方阵,要求a是非奇异方阵,a=123;456;789;b=246;135;7910;a.*b,.*,*,例:,例:,ans=253746558510985133172,ans=281841530497290,a./b,a/b,./,/,ans=0.50000.50000.50004.00001.66671.20001.00000.88890.9000,ans=0.5000003.5000-3.00000.00006.5000-6.00000.0000,a*b,a=123456789;,b=2461357910;,1/22/43/64/15/36/57/78/99/10,a*b-1,对应元素,行列元素,请思考：,2.X=AB是方程A*X=B的解？,3.X=A/B是方程X*A=B的解？,1.A/B=AB?,No,No,Yes,因A*B-1不等于A-1*B,因X=A-1*B代入方程,因X=A*B-1代入方程,不成立,A*A-1*B=B成立,A*B-1*A=B,？,改正,X=/是方程X*A=B的解,X=A/B是方程的解,B/A,X*B=A,3.数组乘方(.)元素对元素的幂例:a=123;b=456;z=a.2z=1.004.009.00,z=a.bz=1.0032.00729.00,142536,z=a2错！要求a是方阵,MatLab语言把多项式表达成一个行向量，该向量中的元素是按多项式降幂排列的。f(x)=anxn+an-1xn-1+a0可用行向量p=anan-1a1a0表示.,四、多项式运算,1.求多项式的根,2.构造多项式3.多项式乘运算,多项式除运算,多项式微分,1.roots由多项式求它的根,p(x)=x4-5x3-17x2+129x-180,%-lxmtbp37.m-p=1-5-17129-180;x=roots(p);x,lxmtbp37x=-5.00004.00003.00003.0000,x=-5433;pp=poly(x);pp,2.poly由向量(根)构造多项式,pp=1.0000-5.0000-17.0000129.0000-180.0000,px=x4-5x3-17x2+129x-180,x=-5433;,px=poly2sym(pp,x)px=poly2str(pp,x),例：求特征多项式的根,p=1.00-6.00-72.00-27.00r=roots(p)r=12.12-5.73r是矩阵A的特征值-0.39,p2=poly(r)p2=1.00-6.00-72.00-27.00,规定多项式系数向量用行向量表示，一组根用列向量表示。,A=123;456;780;p=poly(A),3.conv,convs多项式乘运算,例:a(x)=x2+2x+3;b(x)=4x2+5x+6;c=(x2+2x+3)(4x2+5x+6)a=123;b=456;c=conv(a,b)=conv(123,456)c=4.0013.0028.0027.0018.00p=poly2str(c,x)p=4x4+13x3+28x2+27x+18,4.deconv多项式除运算,a=123;c=4.0013.0028.0027.0018.00d=deconv(c,a)d=4.005.006.00,5.polyder多项式微分,提供了polyder函数多项式的微分。命令格式：polyder(p):求p的微分polyder(a,b):求多项式a,b乘积的微分p,q=polyder(a,b):求多项式a,b商的微分例：a=12345;poly2str(a,x)ans=x4+2x3+3x2+4x+5b=polyder(a)b=4664poly2str(b,x)ans=4x3+6x2+6x+4,五、代数方程组求解,有两种除运算左除和右除。对于方程Ax=b，A为Amn矩阵，有三种情况：当rank(A)=n时，此方程组有唯一解当rank(A)n时，此方程组有无穷多解当时rank(A)rank(A|b)时，此方程组无解,MatLab定义的除运算可以很方便地解上述三种方程,恰定方程组的解,方程Ax=b(A为非奇异n=m存在唯一解)x=A-1b矩阵逆两种解法:x=inv(A)b采用求逆运算解方程x=Ab采用左除运算解方程,x=inv(A)*bx=Abx=x=2.002.003.003.00,例:x1+2x2=82x1+3x2=13,恰定方程Ax=b(n=m)A=12;23;b=8;13;,六、微分方程求解,微分方程求解的仿真算法有多种，常用的有Euler（欧拉法）、RungeKutta(龙格-库塔法),当给定仿真步长时：所以yn+1=yn+hf(xn,yn)n=0,1,2y(x0)=y0,Euler法称一步法，用于一阶微分方程,RungeKutta法龙格-库塔法：实际上取两点斜率的平均斜率来计算的，其精度高于欧拉算法。龙格-库塔法：ode23ode45,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文件(wf.m)functionxdot=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=00.25;t,x=ode23(wf,t0,tf,x0)plot(t,x),figure(2),plot(x(:,1),x(:,2),给定区间、初始值;求解微分方程,建立m文件functiondxdt=wf(t,x)dxdt=x(2);x(2)*(1-x(1)2)-x(1);求解微分方程t,x=ode23(wf,030,00.25);plot(t,x);figure(2)plot(x(:,1),x(:,2),命令格式:T,Y=ODE23(ODEFUN,TSPAN,Y0),fmin单变量寻优函数fmins多变量寻优函数constr有约束条件寻优函数,无约束条件,七、函数优化,例1：f(x)=x2+3x+2在-55区间的最小值f=fmin(x2+3*x+2,-5,5)例2：f(x)=100(x2-x12)2+(a-x1)2在x1=a,x2=a2处有最小值functionf=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；y0=-.4471.9783.115.255.024.664.014.583.455.359.22;p=polyfit(x0,y0,3)%p=56.6915-87.117440.0070-0.9043xx=0:0.01:1;yy=polyval(p,xx);plot(xx,yy,-b,x0,y0,or),九、拟合与插值,多项式拟合,2.插值Matlab提供了一维、二维、三次样条等许多插值选择.,t