MATLAB实用技巧精讲之二ppt课件.ppt
Page 1,第二章 MATLAB的数值计算功能多项式的运算。MATLAB矢量化。MATLAB的符号计算功能。,Page 2,MATLAB的数值计算功能：多项式运算,多项式的表示方法,例如:输入多项式 p(x)=x3-5x2+6x-33p=1 -5 6 -33;poly2sym(p)ans= x3-5*x2+6*x-33,Page 3,MATLAB的数值计算功能：多项式运算,求多项式的根用函数roots求解多项式的根例:求解方程 p(x)=2x4-5x3+6x2-x+9 的所有根p=2 -5 6 -1 9;roots(p)ans= 1.6024 + 1.2709i 1.6024 - 1.2709i -0.3524 + 0.9755i -0.3524 - 0.9755i,Page 4,MATLAB的数值计算功能：多项式运算,多项式的四则运算（1）加法和减法 对于多项式的加法和减法，MATLAB 不提供专用的函数。（2）乘法和除法乘法由函数conv来实现，其使用格式为 c=conv(a,b)除法由deconv实现,其使用格式为 d=deconv(a,b),Page 5,MATLAB的数值计算功能：多项式运算,多项式的乘除运算乘法由函数conv来实现 除法由deconv实现例：计算两多项式的乘除法pd= 6 -195 432 -453 9 -792 -162;poly2sym(pd)ans= 6*x6-195*x5+432*x4-453*x3+9*x2-792*x-162d=3 -90 -18;poly2sym(d)ans=3*x2-90*x-18p1=deconv(pd,d) p1 =2 -5 6 -1 9,Page 6,MATLAB的数值计算功能：多项式运算,对多项式求导数采用函数polyder对多项式求导数，polyder(p)例： p= 6 -195 432 -453 9 -792 -162;poly2sym(p)ans= 6*x6-195*x5+432*x4-453*x3+9*x2-792*x-162Dp=polyder(p)Dp = 36 -975 1728 -1359 18 -792 poly2sym(Dp)ans=36*x5-975*x4+1728*x3-1359*x2+18*x-792,Page 7,MATLAB的数值计算功能：多项式运算,多项式拟合多项式拟合用拟合函数polyfit实现，调用格式polyfit(X, Y, n) X和Y为拟合数据，n为拟合多项式的阶数例：用5阶多项式对0, /2上的正弦函数进行拟合，并图示拟合情况。x=0:pi/20:pi/2; y=sin(x); a=polyfit(x,y,5);x1=0:pi/30:pi*2; y1=sin(x1);y2=a(1)*x1.5+a(2)*x1.4+a(3)*x1.3+a(4)*x1.2+a(5)*x1+a(6);plot(x1,y1,b-,x1,y2,r*);legend(原曲线,拟合曲线)axis(0,7,-1.2,4),Page 8,MATLAB的数值计算功能：保存与再用,把Matlab工作空间中一些有用的数据长久保存下来的方法是生成mat数据文件 save 将工作空间中所有的变量存到matlab.mat文件中 save data 将工作空间中所有的变量存到data.mat文件中 save data a b 将工作空间中a和b变量存到data.mat文件中 load data 将data.mat记录的变量调入工作空间,Page 9,MATLAB的数值计算功能： Matlab矢量化,何为矢量化问题尽量使用矩阵表述避免出现太多（两重或以上）的循环嵌套矢量化技术的目的在于改善程序的性能例1：有一200400矩阵data,for i = 1:200, for j = 1:400 , if data(i,j) 0 data(i,j) = 0; end end end,cputime=0.0470,data(data 0) = 0;,cputime=0.0160,Page 10,MATLAB的数值计算功能： Matlab矢量化,例2：上万个点的计算采用矢量化技术可以大幅度提高运算速度,%一般循环编程i=0;for t=0:.01:10000, i=i+1; y(i)=sin(t);end,%矢量化编程t=0:.01:10000;y=sin(t);,cputime20秒,cputime=0.2340秒,Page 11,MATLAB的数值计算功能： Matlab矢量化,使用向量作索引如果 X 和 V 都为向量X(V) 就是 X(V(1), X(V(2), ., X(V(n)例：X=2 5 8 11 14 17 20 23 26 29 V=4 2 6X(V) 就是X(4) X(2) X(6)，即 11 5 17,Page 12,MATLAB的数值计算功能： Matlab矢量化,创建和操作矩阵从向量构建矩阵例1：A = ones(5,5)*10;A = 10; A = A(ones(5,5) ;例2：M = repmat(1:5, 1,3);,M = 1 1 1 2 2 2 3 3 3 4 4 4 5 5 5,怎样用repmat生成55全10矩阵？,A = 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10,Page 13,MATLAB的数值计算功能： Matlab矢量化,排序、设置和计数max 最大元素 min 最小元素 sort 递增排序 unique 递增排序并且去掉相同元素例1：找出向量中最大元，返回其下标a=5 2 7 89 12 7 525 78;Y,I=max(a)Y = 525I = 7,Page 14,MATLAB的数值计算功能： Matlab矢量化,排序、设置和计数例2：向量元素排序a=5 2 7 89 12 7 525 78;sort(a)ans = 2 5 7 7 12 78 89 525 unique(a)ans = 2 5 7 12 78 89 525,Page 15,MATLAB的符号计算功能,1993年之前的符号计算语言：Maple、Mathematic和MathCAD等1993年Mathworks公司为了解决MATLAB用户对符号计算的需求，从加拿大滑铁卢大学购买了Maple的使用权，并在此基础上开发了符号计算工具箱(Symbolic Toolbox)至此Matlab 不仅具有数值运算功能，具有符号计算功能,Page 16,MATLAB的符号计算功能,本部分包含内容符号表达式、符号矩阵的创建符号线性代数因式分解、展开和简化符号代数方程求解符号微积分符号微分方程,Page 17,MATLAB的符号计算功能,符号运算与数值运算的区别数值运算中必须先对变量赋值，然后才能参与运算。符号运算无须事先对独立变量赋值，运算结果以标准的符号形式表达。符号运算的特点：运算对象可以是没赋值的符号变量可以获得任意精度的解尽管MATLAB 7的符号运算功能很强大，但是一门优秀的语言应该能够博采众家之长。因此MATLAB7语言提供了和MAPLE语言的良好接口，通过maple.m和map.m两个专用的M文件来实现。,Page 18,MATLAB符号变量和符号表达式的生成和使用,定义基本符号对象的指令有两个：sym, syms。它们的常用使用格式如下：f =sym(arg) %把数字、字符串或表达式arg 转换为符号对象f。syms(argv1, argv2, argvk) %把字符arg1, arg2, argk定义为 基本符号对象。syms argv1 argv2 argvk %上述格式的简洁形式。,Page 19,MATLAB符号变量和符号表达式的生成和使用,【例】符号常数形成中的差异。 a1=1/3,pi/7,sqrt(5),pi+sqrt(5)a1 = 0.3333 0.4488 2.2361 5.3777 a2=sym(1/3,pi/7,sqrt(5),pi+sqrt(5)a2 = 1/3, pi/7, sqrt(5), pi+sqrt(5) a3=sym(1/3,pi/7,sqrt(5),pi+sqrt(5)a3 = 1/3, pi/7, sqrt(5), 6054707603575008*2(-50) a23=a2-a3a23 = 0, 0, 0, pi+5(1/2)-189209612611719/35184372088832,Page 20,MATLAB符号变量和符号表达式的生成和使用,【例】使用sym函数定义符号表达式, a=sym(a);b=sym(b);c=sym(c); x=sym(x); f=a*x2+b*x+cf =a*x2+b*x+c,也可以采取整体定义法 f=sym(a*x2+b*x+c) f = a*x2+b*x+c,Page 21,MATLAB符号变量的生成和使用,【例题】用符号计算验证三角等式, syms fai1 fai2 y=simple(sin(fai1)*cos(fai2)-cos(fai1)*sin(fai2)y =sin(fai1-fai2), clear y=simple(sin(fai1)*cos(fai2)-cos(fai1)*sin(fai2)? Undefined function or variable fai1.,Page 22,符号表达式（符号函数）的操作,符号表达式的四则运算 syms x y a b fun1=sin(x)+cos(y); fun2=a+b; fun1+fun2ans =sin(x)+cos(y)+a+b fun1*fun2ans =(sin(x)+cos(y)*(a+b) fun1/fun2ans =(sin(x)+cos(y)/(a+b),Page 23,MATLAB符号矩阵的创建,用Matlab函数sym创建矩阵命令格式：A=sym( ) 需用sym函数 需用 标识注意：数值矩阵A=1,2;3,4有效，A=a,b;c,d 则不识别！例：B = sym(a , 2*b ; 3*a , 0) B = a, 2*b 3*a, 0这就完成了一个符号矩阵的创建。,符号矩阵的每一行的两端都有方括号，这是与 Matlab数值矩阵的一个重要区别。,Page 24,MATLAB符号矩阵的创建,用生成子矩阵的方法生成符号矩阵 命令格式：A= ； ,【例题】用生成子矩阵的方法生成符号矩阵 h1=55,ttt;44,1? Error using = vertcatAll rows in the bracketed expression must have the same number of columns. h2=55,ttt;44,1 h2 =55,ttt44,1 ,该方法不需要调用sym函数，但是要保证同一列的元素具有相同的长度。当要输入的字符串长度不一致时，可以在较短的字符串前边或者后边加上空格符补充,Page 25,MATLAB符号矩阵的创建,由数值矩阵转换成符号矩阵将数值矩阵转化为符号矩阵，函数调用格式：B=sym(A)例：A=1/3,2.5;1/0.7,2/5 A = 0.3333 2.5000 1.4286 0.4000 B=sym(A) B = 1/3, 5/2 10/7, 2/5,A+B=?ans = 2/3, 5 20/7, 4/5,Page 26,MATLAB符号矩阵及符号数组的运算,符号矩阵的四则运算A+B和A-B命令可以实现符号阵列的加法和减法。A*B命令可以实现符号矩阵的乘法。AB命令可以实现矩阵的左除法。 AB为符号线性方程组A*X=B的解。B/A命令可以实现矩阵的右除法。 B/A为符号线性方程组X*A=B的解。,Page 27,MATLAB符号矩阵及符号数组的运算,【例】 m=sym(a,a2,a*2,1/a)m = a, a2, a*2, 1/a n=sym(2*a,b,a,a2)n = 2*a, b, a, a2 m+n ans = 3*a, a2+b, 3*a, 1/a+a2,Page 28,MATLAB符号矩阵及符号数组的运算, x=mnx = 0, 0, 0, 0 0, 0, 0, 0 0, 0, 0, 0 2*a2, b*a, a2, a3y =m/nInfWarning: System is inconsistent. Solution does not exist.y =Inf,Page 29,MATLAB符号矩阵及符号数组的运算,符号矩阵的逆和行列式运算 inv函数可以用来求方正的逆矩阵，命令为：inv(X)； det函数可以求方正的行列式，命令为：det(X)。,Page 30,MATLAB符号矩阵及符号数组的运算,符号数组的四则运算A.*B命令用于符号数组的乘法。A.*B为按参量A 与B对应的分量进行相乘。A.B命令用于数组的左除法运算。A.B为按对应的分量进行相除。A./B命令用于数组的右除法运算。A./B为按对应的分量进行相除。,Page 31,MATLAB符号矩阵及符号数组的运算,【例】符号数组的四则运算 q=sym(3,4,9,6;x,y,z,w;a,b,c,d) q = 3, 4, 9, 6 x, y, z, w a, b, c, d p=sym(x,1/x,x2,x3;a,b,c,d;5,2,3,6)p = x, 1/x, x2, x3 a, b, c, d 5, 2, 3, 6 r=q*p? Error using = sym.mtimesInner matrix dimensions must agree.,Page 32,MATLAB符号矩阵及符号数组的运算, q.*pans = 3*x, 4/x, 9*x2, 6*x3 x*a, y*b, z*c, w*d 5*a, 2*b, 3*c, 6*d q.p ans = 1/3*x, 1/4/x, 1/9*x2, 1/6*x3 a/x, b/y, c/z, d/w 5/a, 2/b, 3/c, 6/d q./pans = 3/x, 4*x, 9/x2, 6/x3 x/a, y/b, z/c, w/d 1/5*a, 1/2*b, 1/3*c, 1/6*d,Page 33,MATLAB符号极限,符号极限limit (F,a)命令中设定变量为x，计算当 时F的极限值。limit(F)命令中设定变量为x，计算当 时F的极限值。limit(F,x,a,right)或limit (F,x,a,left)来计算符号函数F的单侧极限：左极限 或右极限 。,Page 34,MATLAB符号极限,【例】使用limit函数来求符号极限。 syms x a t h; limit(sin(x)/x)ans =1 limit(x-2)/(x2-4),2)ans =1/4 limit(1/x,x,0,right)ans =Inf limit(1/x,x,0,left)ans =-Inf,Page 35,MATLAB符号微积分,符号微积分diff(f) 对缺省变量求一阶导数diff(f,v) 对指定变量v求一阶导数diff(f,v,n) 对指定变量v求n阶导数int(f) 对f表达式的缺省变量求积分int(f,v) 对f表达式的v变量求积分int(f,v,a,b) 对f表达式的v变量在(a,b) 区间求定积分int(被积表达式，积分变量，积分上限，积分下限) 定积分缺省时为不定积分,Page 36,符号微分和求导,【例】求单个自变量时的微分 syms x diff(x3+2*x2+4*x+6) ans = 3*x2+4*x+4 diff(sin(x3),4) ans = 81*sin(x3)*x8-324*cos(x3)*x5-180*sin(x3)*x2,Page 37,符号微分和求导,【例】对多个自变量的函数中的某个自变量求导。 syms x y diff(x*y+y2+sin(x)+cos(y),y)ans = x+2*y-sin(y) diff(x*y+y2+sin(x)+cos(y),y,3) ans =sin(y),Page 38,MATLAB的符号积分,例：计算二重不定积分,Page 39,MATLAB的符号积分,例：计算二重不定积分 syms x yF=int(int(x*exp(-x*y),x),y) F = 1/y*exp(-x*y),Page 40,MATLAB符号积分变换,Fourier变换及其逆变换 时域中的f(t)与它在频域中的Fourier变换F(w)之间存在如下关系：,Fw=fourier(ft,t,w) %求时域函数ft的Fourier变换Fwft=ifourier(Fw,w,t) %求频域函数Fw的Fourier反变换ft,Page 41,MATLAB符号积分变换,Fourier变换及其逆变换【例】求 的Fourier变换。, syms t w ut=sym(Heaviside(t) UT=fourier(ut)UT =pi*dirac(w)-i/w vt=ifourier(UT,w,t) vt =heaviside(t),Page 42,MATLAB符号积分变换,Laplace变换及其逆变换的定义为：,Fs=laplace(ft,t,s) %求时域函数ft的Laplace变换Fsft=ilaplace(Fs,s,t) %求频域函数Fs的Laplace反变换ft,Page 43,MATLAB符号积分变换,【例】求 的Laplace变换。, syms t s; syms a b positive; Dt=sym(Dirac(t-a); Ut=sym(Heaviside(t-b); Mt=Dt,Ut;exp(-a*t)*sin(b*t),t2*exp(-t); Ms=laplace(Mt,t,s)Ms = exp(-s*a), exp(-s*b)/s 1/b/(s+a)2/b2+1), 2/(s+1)3,Page 44,MATLAB符号代数方程的求解,符号代数方程求解Matlab符号运算能够解一般的线性方程、非线性方程及一般的代数方程、代数方程组。命令格式：solve(f) 求一个方程的解solve(f1,f2, fn) 求n个方程的解solve(eq1,eq2, eqn, v1, v2, vn) 求方程组关于指定变量的解（推荐格式）。solve(exp1,exp2,.,expn,v1,v2,vn) 求方程组关于指定变量的解（推荐格式)。,Page 45,MATLAB符号代数方程的求解,符号代数方程求解例1：f = ax2+bx+c 求解f=a*x2+b*x+c;solve(f) %对缺省变量x求解ans = 1/2/a*(-b+(b2-4*a*c)(1/2) 1/2/a*(-b-(b2-4*a*c)(1/2),计算机格式,一般格式,Page 46,MATLAB符号代数方程的求解,符号代数方程求解solve(f , b) %对指定变量b求解ans = -(a*x2+c)/x例2：符号方程cos(x)=sin(x) 求解f1=solve(cos(x)=sin(x)f1 = 1/4*pi,Page 47,MATLAB符号代数方程的求解,符号代数方程求解例：解方程组g1=x+y+z=1; g2=x-y+z=2; g3=2*x-y-z=1;f=solve(g1,g2,g3)f = x: 1x1 sym y: 1x1 sym z: 1x1 symf.x, f.y, f.zans= 2/3, -1/2, 5/6,f=solve(x+y+z=1,x-y+z=2,2*x-y-z=1),x,y,z=solve(x+y+z=1,x-y+z=2,2*x-y-z=1),x =2/3y =-1/2z =5/6,Page 48,MATLAB符号微分方程求解,符号微分方程求解用一个函数可以方便地得到微分方程的符号解，求解指令：dsolve调用格式：dsolve(f, g)f微分方程，可多至12个微分方程的求解；g为初始条件默认自变量为 x,可任意指定自变量t, u等微分方程的各阶导数项以大写字母D表示y的一阶导数 Dyy的二阶导数 D2yy的 n 阶导数 Dny,Page 49,MATLAB符号微分方程求解,符号微分方程求解一阶微分方程例1：求 x,y,z=dsolve(Dx=y,Dy=-z,Dz=x,y(0)=1,x(0)=0,z(0)=0)x =-1/3*exp(-t)+1/3*3(1/2)*exp(1/2*t)*sin(1/2*3(1/2)*t)+1/3*exp(1/2*t)*cos(1/2*3(1/2)*t) y =1/3*exp(-t)+2/3*exp(1/2*t)*cos(1/2*3(1/2)*t) z =1/3*exp(-t)+1/3*3(1/2)*exp(1/2*t)*sin(1/2*3(1/2)*t)-1/3*exp(1/2*t)*cos(1/2*3(1/2)*t),Page 50,MATLAB符号微分方程求解,符号微分方程求解 二阶微分方程求解例3：求该方程的解y=dsolve(D2y+2*Dy+2*y=0,y(0)=1,Dy(0)=0)ans = exp(-x)*cos(x)+exp(-x)*sin(x)ezplot(y) %方程解y(t)的时间曲线图,Page 51,补充说明,1=conv(1,0,2,conv(1,4,1,1);p2=1,0,1,1;q,r=deconv(p1,p2)cq=商多项式为;cr=余多项式为disp(cq,poly2str(q,s),disp(cr,poly2str(r,s)商多项式为 s + 5余多项式为 5 s2 + 4 s + 3,q,r=deconv(p1,p2) %多项式p1被p2除的商项式为 q，而余多项式为r,Page 52,作业,1.求d+n/2+p/2=q, n+d+q-p=0, q+d-n/4=p, q+p-n-8d=1的线性方程组的解。 2.验证Laplace时移性质：,3.求 的Fourier变换，在此x是参量， t是时间变量。,Page 53,作业,4.对函数sin(x)+cos(x)取任意多个数据点，求拟合6阶多项式，并图示拟合情况。5.求符号表达式 的5次微分。,
