二讲Matlab符号计算

《二讲Matlab符号计算》由会员分享，可在线阅读，更多相关《二讲Matlab符号计算（127页珍藏版）》请在装配图网上搜索。

1、参考文献:MATLAB程序设计教程程序设计教程 李海涛，邓樱 编著高等教育出版社，2002 所谓符号计算是指在运算时，无须事先对变量赋值，而将所得到结果以标准的符号形式来表示。Matlab本身并没有符号计算功能，MathWorks公司于1993年购买了著名的符号数学软件Maple的使用权，并利用Maple已有的函数库，开发了符号数学工具箱。符号运算与数值运算的区别：符号运算与数值运算的区别：数值运算中必须先对变量赋值，然后才能参与运算。符号运算无须事先对独立变量赋值，运算结果以标准的 符号形式表达。符号计算是精确计算，符号计算是精确计算，可以获得任意精度的解可以获得任意精度的解；Symboli

2、c Math Toolbox符号运算工具包通过调用Maple软件实现符号计算的。Maple软件主要功能是符号运算，它占据符号软件的主导地位。符号运算的特点：符号运算的特点：运算对象可以是没赋值的符号变量运算对象可以是没赋值的符号变量。符号计算定义在符号变量的基础上，符号表达式计算前必须定义符号变量符号表达式计算前必须定义符号变量;Matlab的符号数学工具箱的主要功能包括：符号表达式的创建、符号矩阵的运算、符号表达式的化简、符号微积分、符号方程（主要指代数方程与常微分方程）的求解、符号函数绘图等。符号计算的计算速度较慢符号计算的计算速度较慢；符号计算的运算符和基本数学函数与数值计算中的运算符号

3、计算的运算符和基本数学函数与数值计算中的运算符和基本数学函数几乎完全相同符和基本数学函数几乎完全相同。一、符号对象的创建一、符号对象的创建（一）（一）创建符号变量和表达式创建符号变量和表达式参与符号运算的对象可以是符号变量、符号表达式或符号矩阵。符号变量要先定义，后引用。符号变量要先定义，后引用。Matlab提供了两个建立符号变量的函数：sym和syms，两个函数的用法不同。sym函数函数 sym函数的主要功能是创建创建单个单个符号变量符号变量，以便进行符号运算，也可以用于创建符号表达式或符号矩阵。用sym函数创建符号变量的一般格式为：符号量名符号量名=sym(符号字符串符号字符串)符号字符串

4、符号字符串可以是常量、变量、函数、表达式常量、变量、函数、表达式或方程方程。如如:f1=ax2+bx+c 二次三项式 f2=ax2+bx+c=0 方程 f3=Dy+y2=1 微分方程u 符号表达式或符号方程可以赋给符号变量，以后调用方便；也可以不赋给符号变量直接参与运算a=sym(a);b=sym(b);c=sym(c);d=sym(d);%定义4个符号变量w=10;x=5;y=-8;z=11;%定义4个数值变量A=a,b;c,d%建立符号矩阵AB=w,x;y,z%建立数值矩阵Bdet(A)%计算符号矩阵A的行列式det(B)%计算数值矩阵B的行列式 使用符号变（常）量进行代数运算和数值变（常

5、）量进行的运算是不同的。例例1 考察符号变量和数值变量的差别。例例2 定义一个符号函数 fxy=(a*x2+b*y2)/c2，分别求该函数对x、y的导数和对x的积分。a=sym(a);b=sym(b);c=sym(c);%定义符号变量x=sym(x);y=sym(y);fxy=(a*x2+b*y2)/c2;%生成符号函数 diff(fxy,x)%符号函数fxy对x求导数diff(fxy,y)%符号函数fxy对y求导数 int(fxy,x)%符号函数fxy对x求积分ans=2*a*x/c2 ans=2*b*y/c2 ans=1/c2*(1/3*a*x3+b*y2*x)例例3 比较符号常数与数值在

6、代数运算时的差别。pi1=sym(pi);k1=sym(8);%定义符号变量k2=sym(2);k3=sym(3);pi2=pi;r1=8;r2=2;r3=3;%定义数值变量sin(pi1/3)%计算符号表达式值 sin(pi2/3)%计算数值表达式值sqrt(k1)%计算符号表达式值sqrt(r1)%计算数值表达式值sqrt(k3+sqrt(k2)%计算符号表达式值sqrt(r3+sqrt(r2)%计算数值表达式值ans=1/2*3(1/2)ans=ans=2*2(1/2)ans=ans=(3+2(1/2)(1/2)ans=符号表达式由符号变量、函数、算术运算符等组成。符号表达式的书写格式与

7、数值表达式相同。利用单引号来生成；用sym函数建立；使用已定义的符号变量组成。建立符号表达式的方法：建立符号表达式的方法：例例4 已知一复数表达式 z=x+i*y，试求其共轭复数，并求该表达式与其共轭复数乘积的多项式。为了使乘积表达式 x2+y2非负，这里，把变量x和y定义为实数。x=sym(x,real);y=sym(x,real);若要去掉x的属性,可以使用下面语句x=sym(x,unreal)将x创建为纯格式的符号变量。z=x+i*y;%定义复数表达式conj(z);%求共轭复数expand(z*conj(z)%求表达式与其共轭复数乘积的多项式ans=x2+y2 syms函数函数 函数s

8、ym一次只能定义一个符号变量，使用不方便。Matlab 提供了另一个函数syms，一次可以定义多个符号变量。syms函数的一般调用格式为：syms 符号变量名符号变量名1 符号变量名符号变量名2 符号变量名符号变量名n用这种格式定义符号变量时不要在变量名上加字符串分界符()，变量间用空格而不要用逗号分隔变量间用空格而不要用逗号分隔。U=sym(3*x2+5*y+2*x*y+6)%定义符号表达式Usyms x y;%建立符号变量x、yV=3*x2+5*y+2*x*y+6%定义符号表达式V2*U-V+6%求符号表达式的值例例5 用两种方法建立符号表达式。U=3*x2+5*y+2*x*y+6V=3*

9、x2+5*y+2*x*y+6ans=3*x2+5*y+2*x*y+12（二）（二）符号函数的生成符号函数的生成 将表达式中的自变量定义为符号变量后，赋值给符号函数名，即可生成符号函数。例如有一数学表达式：syms a b c x y%定义符号fxy=(a*x2+b*y2)/c2%生成符号函数其用符号表达式生成符号函数fxy的过程为：生成符号函数fxy后，即可用于微积分等符号计算。符号表达式中变量的确定 Matlab 中的符号可以表示符号变量和符号常数。findsym可以帮助用户查找一个符号表达式中的符号变量。该函数的调用格式为：findsym（S,n）函数返回符号表达式S中的n个符号变量，若没

10、有指定n，则返回 S 中的全部符号变量。在求函数的极限、导数和积分时，如果用户没有明确指定自变量，Matlab将按缺省原则确定自变量并对其进行相应微积分运算。可用可用findsym(S,1)查找系统的缺省变查找系统的缺省变量。事实上，量。事实上，Matlab 按离字符按离字符x最近原则确定缺省变量最近原则确定缺省变量例例6 查询符号函数 f=xn g=sin(at+b)中的系统缺省变量。表示函数 f 中查询的 1 个系统缺省变量为 x。syms a b n t x%定义符号变量f=xn;%给定符号函数g=sin(a*t+b);findsym(f,1)%在f函数中查询1个系统缺省变量ans=x（

11、三）（三）符号矩阵的创建符号矩阵的创建 符号矩阵也是一种符号表达式，所以前面介绍的符号表达式运算都可以在矩阵意义下进行。但应注意这些函数作用于符号矩阵时，是分别作用于矩阵的每一个元素。用字符串直接创建符号矩阵用字符串直接创建符号矩阵例例7 A=a,2*b;3*a,0 A=a,2*b3*a,0u 模仿matlab数值矩阵的创建方法；u 每一行都要用方括号括起来；每一行都要用方括号括起来；u 需保证同一列中各行元素字符串有相同的长度需保证同一列中各行元素字符串有相同的长度（较短的字符串要用空格补齐长度）。用用sym函数创建符号矩阵函数创建符号矩阵命令格式：A=sym()u 符号矩阵内容同数值矩阵：

12、矩阵元素可以是任何不带等号的符号表达式，各矩阵元素之间用逗号或空格分隔，各行之间用分号分隔，各元素字符串的长度可以不相等各元素字符串的长度可以不相等；u 需用sym函数定义；u 需用 标识。例例8 A=sym(a,2*b;3*a,0)A=a,2*b3*a,0 注意：符号矩阵的每一行的两端都有方括号符号矩阵的每一行的两端都有方括号，这是与 Matlab数值矩阵的一个重要区别。将数值矩阵转化为符号矩阵将数值矩阵转化为符号矩阵A=1/3,2.5;1/0.7,2/5A=sym(A)ans=1/3,5/210/7,2/5 函数调用格式：sym(A)例例9A=sym(1/3,5/2;10/7,2/5)ev

13、al(A)将符号矩阵转化为数值矩阵函数调用格式：eval(A)A=1/3,5/210/7,2/5ans=例例10 numeric(A)或 符号矩阵的修改u 指令修改指令修改用A(,)=new 来修改:用new替换 A 中的某个元素；用A1subs(A,new)来修改:用new替换 A中的自由变量；用A1=subs(A,new,old)来修改:用new替换 A中的变量oldu 直接修改直接修改 可用、键找到所要修改的矩阵，直接修改；例例11A2=subs(A,sym(c),b)A2=a,2*c 3*a,4*c A=sym(a,2*b;3*a,0)A=a,2*b 3*a,0A(2,2)=4*bA=

14、a,2*b 3*a,4*bsyms a x;f=a*sin(x)+5 f=a*sin(x)+5 f1=subs(f,sin(x),sym(y)f1=a*y+5f2=subs(f,a,x,2,pi/3)f2=3(1/2)+5f3=subs(f,a,x,2,pi/3)f3=f4=subs(subs(f,a,2),x,0:pi/6:pi)f4=f5=subs(f,a,x,0:6,0:pi/6:pi)f5=实际上，Matlab的符号数学工具箱提供了两个符号表达式的替换函数subexpr和subs，可以通过符号替换使表达式的输出形式简化，以得到一个简单的表达式。l 符号表达式的替换符号表达式的替换 将表

15、达式中重复出现的字符串用变量代替将表达式中重复出现的字符串用变量代替的函数为subexpr，其调用形式为：此函数用变量Sigma（字符或字符串）的值代替符号表达式S中重复出现的字符串，Y返回替换后的结果。Y,Sigma=subexpr(S,Sigma)例例12 求解并化简 3 次方程 x3+ax+1=0 的符号解。t=solve(x3+a*x+1=0)t=1/6*(-108+12*(12*a3+81)(1/2)(1/3)-2*a/(-108+12*(12*a3+81)(1/2)(1/3)-1/12*(-108+12*(12*a3+81)(1/2)(1/3)+a/(-108+12*(12*a3+

16、81)(1/2)(1/3)+1/2*i*3(1/2)*(1/6*(-108+12*(12*a3+81)(1/2)(1/3)+2*a/(-108+12*(12*a3+81)(1/2)(1/3)-1/12*(-108+12*(12*a3+81)(1/2)(1/3)+a/(-108+12*(12*a3+81)(1/2)(1/3)-1/2*i*3(1/2)*(1/6*(-108+12*(12*a3+81)(1/2)(1/3)+2*a/(-108+12*(12*a3+81)(1/2)(1/3)r,s=subexpr(t,s)r=1/6*s(1/3)-2*a/s(1/3)-1/12*s(1/3)+a/s(

17、1/3)+1/2*i*3(1/2)*(1/6*s(1/3)+2*a/s(1/3)-1/12*s(1/3)+a/s(1/3)-1/2*i*3(1/2)*(1/6*s(1/3)+2*a/s(1/3)s=-108+12*(12*a3+81)(1/2)而函数 subs 用指定符号替换符号表达式中的某一特定符号。数值运算中，所有矩阵运算操作指令都比较直观、简单。而符号运算就不同了，所有涉及符号矩阵运算的操作都有专用符号矩阵函数来进行，这些函数作用于单个的数据无意义这些函数作用于单个的数据无意义。符号矩阵的运算 symsize 求符号矩阵维数；charploy 特征多项式；inverse 逆矩阵；jord

18、an 约当标准型；simple 符号矩阵简化。u 符号矩阵运算函数：transpose(S)返回S矩阵的转置矩阵；determ(S)返回S矩阵的行列式值；colspace(S)返回S矩阵列空间的基；Q,D=eigensys(S)Q返回S矩阵的特征向量，D返回S矩阵的特征值。例例13 计算3阶范得蒙矩阵行列式的值。设A是一个由符号变量a,b,c确定的范得蒙矩阵。其实，曾介绍过的许多应用于数值矩阵的函数，如diag、triu、tril、inv、det、rank、eig等，也可直接应用于符号矩阵。syms a b c;U=a,b,c;A=1,1,1;U;U.2%建立范得蒙符号矩阵det(A)%计算A

19、的行列式值A=1,1,1 a,b,c a2,b2,c2 ans=b*c2-c*b2-a*c2+a*b2+a2*c-a2*b二、二、符号表达式运算符号表达式运算（一）（一）符号表达式的四则运算符号表达式的四则运算symadd(a,d)符号符号表达式表达式的加；的加；symsub(a,b)符号符号表达式表达式的减；的减；symmul(a,b)符号符号表达式表达式的乘；的乘；symdiv(a,b)符号符号表达式表达式的除；的除；sympow(a,b)符号符号表达式表达式的幂运算；的幂运算；symop(a,b)符号符号表达式表达式的综合运算。的综合运算。符号表达式运算的函数：例例14 符号表达式的四则

20、运算示例。syms x y;f1=2*x+x2*x-5*x+x3%符号表达式的结果为最简形式f2=2*x/(5*x)%符号表达式的结果为最简形式f3=(x+y)*(x-y)%符号表达式的结果不是x2-y2，%而是(x+y)*(x-y)f1=-3*x+2*x3 f2=2/5 f3=(x+y)*(x-y)（1）f=cos(x);g=sin(2*x);symop(f,/,g,+,f,*,g)ans=cos(x)/sin(2*x)+cos(x)*sin(2*x)（2）f=2*x2+3*x-5;g=x2+x-7;h=symadd(f,g)h=3*x2+4*x-12syms xf=2*x2+3*x-5;g

21、=x2+x-7;h=f+gsyms xf=cos(x);g=sin(2*x);f/g+f*g(1）f=2*x2+3*x-5;g=x2+x-7;例例15（2）f=cos(x);g=sin(2*x);h=3*x2+4*x-12ans=cos(x)/sin(x)+cos(x)*sin(x)该函数将符号表达式s转换为分子和分母都是整系数分子和分母都是整系数的最佳多项式，并分别存放在n与d中。（二）（二）符号表达式的分式通分符号表达式的分式通分 如果符号表达式是一个有理分式或可以展开为有理分式，可利用numden函数来提取符号表达式中的分子或分母。其一般调用格式为：n,d=numden(s)例例16 对

22、表达式 进行通分。syms x yf=x/y+y/x;n,d=numden(f)n=x2+y2d=y*xfactor(S)对S分解因式，S是符号表达式或符号矩阵;pretty(S)将符号表达式S按照类似书写习惯的方式显示;expand(S)对S进行展开，S是符号表达式或符号矩阵;collect(S)对S合并同类项，S是符号表达式或符号矩阵;collect(S,v)对S按变量v合并同类项，S是符号表达式或符号矩阵。（三）（三）符号表达式的符号表达式的因式分解、展开与合并因式分解、展开与合并 Matlab 提供了符号表达式的因式分解与展开的函数，函数的调用格式为：例例17 将表达式(x9-1)分解

23、为多个因式。syms xfactor(x9-1)ans=(x-1)*(x2+x+1)*(x6+x3+1)syms xf=(x+1)5;expand(f)syms a b x y;A=2*a2*b3*x2-a*b4*x3+10*a*b6*x4,3*x*y-5*x2;4,a3-b3;factor(A)%对A的每个元素分解因式pretty(A)例例18 对符号矩阵A的每个元素分解因式。例例19 展开表达式f=(x+1)5ans=a*b3*x2*(10*b3*x2-b*x+2*a),-x*(-3*y+5*x)4,(a-b)*(a2+a*b+b2)2 3 2 4 3 6 4 2 2 a b x -a b

24、 x +10 a b x 3 y x-5 x 3 3 4 a -b ans=x5+5*x4+10*x3+10*x2+5*x+1syms x t;f=x*(x*(x-6)+12)*t;collect(f)collect(f,t)例例20 对于表达式 f=x(x(x-6)+12)t，分别将自变量x和t的同类项合并。ans=t*x3-6*t*x2+12*t*xans=x*(x*(x-6)+12)*t（四）（四）符号表达式的化简符号表达式的化简Matlab提供的对符号表达式化简的函数有：syms x y;s=(x2+y2)2+(x2-y2)2;r,how=simple(s)%Matlab自动调用多种函

25、数对s进行化简，并显示每步结果。simple(s)：调用Matlab 的其他函数对表达式进行综合化简，以寻求s的最简形式，并显示化简过程。调用格式为：simplify(s)：应用函数规则对s进行化简;r,how=simple(s)返回s的最简化形式，r为返回的简化形式，how为化简过程中使用的主要方法，包括simplify,radsimp,combine,collect,factor,convert,expand等。例例21 化简r=2*x4+2*y4 how=simplify（五）（五）符号表达式符号表达式嵌套形式重写嵌套形式重写将符号多项式s用嵌套形式表示，即用多层括号的形式表示。Horn

26、er函数可以实现此功能。该函数的调用格式为：syms xhorner(x3-6*x2+11*x-6)horner(s)例例22 将表达式x3-6*x2+11*x-6用嵌套形式表示。ans=-6+(11+(-6+x)*x)*x（六）任意精度的数学运算（六）任意精度的数学运算 在symbolic中有三种不同的算术运算：u数值类型matlab的浮点算术运算；u有理数类型maple的精确符号运算；uvpa类型maple的任意精度算术运算。l浮点算术运算1/2+1/3 (定义输出格式format long)ans=l符号运算sym(1/2)+(1/3)ans=5/6 精确解l任意精度算术运算digits

27、(n)设置可变精度，缺省16位vpa(x,n)显示可变精度计算digits(25)vpa(1/2+1/3)vpa(5/6,40)ans=.8333333333333333333333333333333333333333 a=sym(1/4,exp(1);log(3),3/7)a=1/4,exp(1)log(3),3/7vpa(a,10)ans=.2500000000,2.718281828 1.098612289,.4285714286例例23ans=.8333333333333333333333333 时域中的f(t)与它在频域中的Fourier变换F()之间存在如下关系：例例24 求函数的

28、傅立叶变换及其逆变换。（七）积分变换（七）积分变换傅立叶傅立叶(Fourier)变换变换u fourier(ft,t,x):求“时域”函数f(t)的Fourier变换F(x)u ifourier(Fx,x,t):求“频域”函数F(x)的Fourier反变换f(t)syms x t;y=abs(t);Fx=fourier(y,t,x)%求y的傅立叶变换ft=ifourier(Fx,x,t)%求Fx的傅立叶逆变换Fx=-2/x2 ft=t*(2*heaviside(t)-1)例例25 根据Fourier变换定义，用积分指令求方波脉冲 的Fourier变换。syms A t wsyms tao po

29、sitiveyt=heaviside(t+tao/2)-heaviside(t-tao/2);Yw=fourier(A*yt,t,w)Yt=ifourier(Yw,w,t)（1）求Fourier变换Yw=2*A/w*sin(1/2*tao*w)（2）用反变换验算Yt=A*(heaviside(t+1/2*tao)-heaviside(t-1/2*tao)yt3=subs(yt,tao,3)yt3=heaviside(t+3/2)-heaviside(t-3/2)Yw3=subs(Yw,A,tao,1,3)Yw3=2/w*sin(3/2*w)subplot(2,1,1)Ht=ezplot(yt3

30、,-3,3);set(Ht,Color,r,LineWidth,3)subplot(2,1,2),ezplot(Yw3)（3）曲线绘制（设tao=3,A=1）时域方波及其Fourier变换 u ilaplace(Fs,s,t)：求“频域”函数F(s)的Laplace反变换f(t)例例26 计算y=x2的拉普拉斯变换及其逆变换。拉普拉斯拉普拉斯(Laplace)变换变换Laplace变换和反变换的定义为：u laplace(ft,t,s)：求“时域”函数f(t)的Laplace变换F(s)syms t y s;y=t2;Fs=laplace(y,t,s)%对函数y进行拉普拉斯变换ft=ilapl

31、ace(Fs,s,t)%对函数Ft进行拉普拉斯逆变换Fs=2/s3 ft=t2 u ztrans(fn,n,z)：求“时域”序列fn的Z变换Fzsyms n zfn=exp(-n);Fz=ztrans(fn,n,z)%求fn的Z变换f=iztrans(Fz,z,n)%求Fz的逆Z变换Z变换变换 当函数f(x)呈现为一个离散的数列f(n)时，数列f(n)的Z变换及其反变换的定义为：u iztrans(Fz,z,n)：求“频域”序列Fz的Z反变换fn例例27 求数列 fn=e-n的Z变换及其逆变换。Fz=z/exp(-1)/(z/exp(-1)-1)f=exp(-1)n三、三、符号微积分符号微积分

32、（一）（一）符号极限符号极限limit(f,x,a,left)：求符号函数f在x=a处的左极限。limit(f,x,a)：求符号函数f(x)的极限值。即计算当变量x趋近于常数a时，函数f(x)的极限值；limit(f,a)：求符号函数f(x)的极限值。由于没有指定符号函数f(x)的自变量,则使用该格式时,符号函数f(x)的变量为函数findsym(f)确定的默认自变量,既变量x趋近于a;limit(f)：求符号函数f(x)的极限值。符号函数f(x)的变量为函数findsym(f)确定的默认变量；没有指定变量的目标值时，系统默认变量趋近于0，即a=0的情况；limit(f,x,a,right)：

33、求符号函数f在x=a处的右极限；函数 limit 用于求符号函数 f 的极限。系统可以根据用户要求，计算变量从不同方向趋近于指定值的极限值。例例28 求极限syms xf=(x*(exp(sin(x)+1)-2*(exp(tan(x)-1)/sin(x)3;w=limit(f)w=-1/2例例29 求极限。lsyms a m x;f=(x*(exp(sin(x)+1)-2*(exp(tan(x)-1)/(x+a);limit(f,x,a)ans=(1/2*a*exp(sin(a)+1/2*a-exp(tan(a)+1)/alsyms x t;limit(1+2*t/x)(3*x),x,inf)

34、ans=exp(6*t)lsyms x;f=x*(sqrt(x2+1)-x);limit(f,x,inf,left)ans=1/2lsyms x;f=(sqrt(x)-sqrt(2)-sqrt(x-2)/sqrt(x*x-4);limit(f,x,2,right)ans=-1/2（二）（二）符号导数符号导数diff函数用于对符号表达式求导数。其一般调用格式为：diff(s)：没有指定变量和导数阶数，则系统按findsym函数指示的默认变量对符号表达式s求一阶导数。diff(s,v)：以v为自变量，对符号表达式s求一阶导数。diff(s,n)：按findsym函数指示的默认变量对符号表达式s求n

35、阶导数，n为正整数。diff(s,v,n)：以v为自变量，对符号表达式s求n阶导数。例例30 求导数：x=sym(x);diff(sin(x2)ans=2*cos(x2)*x例例31 求函数的导数。syms a b t x y z;f=sqrt(1+exp(x);diff(f)%未指定求导变量和阶数，按缺省规则处理f=x*cos(x);diff(f,x,2)%求f对x的二阶导数diff(f,x,3)%求f对x的三阶导数ans=1/2/(1+exp(x)(1/2)*exp(x)ans=-2*sin(x)-x*cos(x)ans=-3*cos(x)+x*sin(x)f1=a*cos(t);f2=b

36、*sin(t);diff(f2)/diff(f1)%按参数方程求导公式求y对x的导数(diff(f1)*diff(f2,2)-diff(f1,2)*diff(f2)/(diff(f1)3%求y对x的二阶导数f=x*exp(y)/y2;diff(f,x)%z求f对x的偏导数diff(f,y)%z求f对y的偏导数f=x2+y2+z2-a2;zx=-diff(f,x)/diff(f,z)%按隐函数求导公式求z对x的偏导数zy=-diff(f,y)/diff(f,z)%按隐函数求导公式求z对y的偏导数ans=-b*cos(t)/a/sin(t)ans=-(a*sin(t)2*b+a*cos(t)2*b

37、)/a3/sin(t)3ans=exp(y)/y2 ans=x*exp(y)/y2-2*x*exp(y)/y3zx=-x/z zy=-y/z（三）（三）符号符号积分积分int(s,v,a,b)：求定积分运算。a,b分别表示定积分的下限和上限。该函数求被积函数在区间a,b上的定积分。a和b可以是两个具体的数，也可以是一个符号表达式，还可以是无穷(inf)。当函数f关于变量x在闭区间a,b上可积时，函数返回一个定积分结果。当a,b中有一个是inf时，函数返回一个广义积分。当a,b中有一个符号表达式时，函数返回一个符号函数。符号积分由函数 int 来实现。该函数的一般调用格式为：int(s)：没有指

38、定积分变量和积分阶数时,系统按findsym函数指示的默认变量对被积函数或符号表达式s求不定积分int(s,v)：以v为自变量，对被积函数或符号表达式s求不定积分。例例32 求不定积分。x=sym(x);f=(3-x2)3;int(f)%求不定积分f=sqrt(x3+x4);int(f)%求不定积分g=simple(ans)%调用simple函数对结果化简ans=27*x-1/7*x7+9/5*x5-9*x3 ans=-1/48*(x3+x4)(1/2)*(-16*(x2+x)(3/2)+12*(x2+x)(1/2)*x+6*(x2+x)(1/2)-3*log(1/2+x+(x2+x)(1/2

39、)/x/(x+1)*x)(1/2)g=1/3*(x+1)*x)(1/2)*x2+1/12*(x+1)*x)(1/2)*x-1/8*(x+1)*x)(1/2)+1/16*log(1/2+x+(x+1)*x)(1/2)x=sym(x);t=sym(t);int(abs(1-x),1,2)%求定积分f=1/(1+x2);int(f,-inf,inf)%求定积分int(4*t*x,x,2,sin(t)%求定积分f=x3/(x-1)100;I=int(f,2,3)%用符号积分的方法求定积分double(I)%将上述符号结果转换为数值例例33 求定积分。ans=1/2 ans=pi ans=2*t*(si

40、n(t)2-4)I=ans=例例34 计算二重不定积分syms x yF=int(int(x*exp(-x*y),x),y)F=1/y*exp(-x*y)例例35 轴的长度为10米，若该轴的线性密度计算公式是f(x千克/米(其中x为距轴的端点距离)，求轴的质量syms x;f=6+0.3*x;m=int(f,0,10)m=quad(fx,0,10,1e-6)(1)符号函数积分。(2)数值积分。先建立一个函数文件：function fx=fx(x)fx=6+0.3*x;m=75再在MATLAB命令窗口，输入命令：m=75syms a b c z;f=pi*a*b*(c2-z2)/c2;V=int

41、(f,z,-c,c)例例36 求椭球的体积。V=4/3*pi*a*b*csyms t;x=3*t;y=3*t2;z=2*t3;f=diff(x,y,z,t)%求x,y,z对参数t的导数g=sqrt(f*f)%计算一型积分公式中的根式部分l=int(g,t,0,1)%计算曲线c的长度例例37 求空间曲线c从点(0，0，0)到点(3，3，2)的长度。求曲线c的长度是曲线一型，命令如下：f=3,6*t,6*t2 g=3*(1+4*t*conj(t)+4*t2*conj(t)2)(1/2)l=5（四）（四）级数的符号级数的符号求和求和求无穷级数的和需要符号表达式求和函数symsum，其调用格式为：sy

42、msum(s,v,n,m)其中s表示一个级数的通项，是一个符号表达式。v是求和变量，v省略时使用系统的默认变量。n和m是求和的开始项和末项。n=sym(n);s1=symsum(1/n2,n,1,inf)s2=symsum(-1)(n+1)/n,1,inf)%未指定求和变量,缺省为ns3=symsum(n*xn,n,1,inf)%此处的求和变量n不能省s4=symsum(n2,1,100)%计算有限级数的和 例例38 求级数之和。s1=1/6*pi2 s2=log(2)s3=3*t/(3*t-1)2 s4=338350（五）函数的泰勒展开（五）函数的泰勒展开 Matlab中提供了将函数展开为幂

43、级数的函数taylor，其调用格式为：该函数将函数f按变量v展开为泰勒级数，展开到第n项(即变量v的n-1次幂)为止，n的缺省值为6。v的缺省值与diff函数相同。参数a指定将函数f在自变量v=a处展开，a的缺省值是0。taylor(f,v,n,a)例例39 求函数在指定点的泰勒展开式。ans=1+2*x+2*x2-2*x4 ans=1/6*x2+x3+119/72*x4+239/72*x5x=sym(x);f1=(1+x+x2)/(1-x+x2);f2=sqrt(1-2*x+x3)-(1-3*x+x2)(1/3);taylor(f1,x,5)%展开到x的4次幂时应选择n=5taylor(f2

44、,6)x=sym(x);p=1+3*x+5*x2-2*x3;f=taylor(p,x,-1,4)例例40 将多项式表示成x+1的幂的多项式。f=-8-13*x+11*(1+x)2-2*(1+x)3x=sym(x);f=(1-x)(1/12);%定义函数，4000(1/12)=2f(96/212)g=taylor(f,4)%求f的泰勒展开式g,有4000(1/12)2g(96/212)b=96/212;a=1-b/12-11/288*b2-253/10368*b3%计算g(b)2*a%求4000(1/12)的结果4000(1/12)%用MATLAB的乘方运算直接计算例例41 应用泰勒公式近似计算

45、4000(1/12)。g=1-1/12*x-11/288*x2-253/10368*x3 a=ans=ans=（六）函数的傅立叶级数（六）函数的傅立叶级数function mfourier=mfourier(f,n)syms x a b c;mfourier=int(f,-pi,pi)/2;%计算a0for i=1:n a(i)=int(f*cos(i*x),-pi,pi);b(i)=int(f*sin(i*x),-pi,pi);mfourier=mfourier+a(i)*cos(i*x)+b(i)*sin(i*x);endreturn下述文件可以求任意函数的傅立叶级数。调用该函数时，需给出

46、被展开的符号函数f和展开项数n，不可缺省。x=sym(x);a=sym(a);f=x;mfourier(f,5)%求f(x)=x的傅立叶级数的前5项例例42 在-，区间展开函数为傅立叶级数。f=abs(x);mfourier(f,5)%求f(x)=|x|的傅立叶级数的前5项ans=2*pi*sin(x)-pi*sin(2*x)+2/3*pi*sin(3*x)-1/2*pi*sin(4*x)+2/5*pi*sin(5*x)ans=1/2*pi2-4*cos(x)-4/9*cos(3*x)-4/25*cos(5*x)f=cos(a*x);mfourier(f,6)%求f(x)=cos(ax)的傅立

47、叶级数的前6项ans=-2*sin(pi*a)/(a2-1)*sin(x)+4*sin(pi*a)/(a2-4)*sin(2*x)-6*sin(pi*a)/(a2-9)*sin(3*x)+8*sin(pi*a)/(a2-16)*sin(4*x)ans=sin(pi*a)/a-2*a*sin(pi*a)/(a2-1)*cos(x)+2*a*sin(pi*a)/(a2-4)*cos(2*x)-2*a*sin(pi*a)/(a2-9)*cos(3*x)+2*a*sin(pi*a)/(a2-16)*cos(4*x)-2*a*sin(pi*a)/(a2-25)*cos(5*x)+2*a*sin(pi*a

48、)/(a2-36)*cos(6*x)f=sin(a*x);mfourier(f,4)%求f(x)=sin(ax)的傅立叶级数的前4项四、四、符号方程的求解符号方程的求解（一）（一）符号符号代数方程（组）的求解代数方程（组）的求解 在Matlab中，求解用符号表达式表示的代数方程可由函数solve实现，其调用格式为：solve(s1,s2,sn,v1,v2,vn)：求解符号表达式s1,s2,sn组成的代数方程组，求解变量分别为v1,v2,vn。solve(s,v)：求解符号表达式s的代数方程,求解变量为vsolve(s)：求解符号表达式s的代数方程，求解变量为默认变量。u 当方程组不存在符号解时

49、，又无其他自由参数，则给出数值解。例例43 求线性方程组AX=b的解。此外，Matlab中提供了一个求解线性代数方程组求解线性代数方程组的函数linsolve，其调用格式为：linsolve(A,b)方程组符号求解符号求解的命令如下：syms a11 a12 a13 a21 a22 a23 a31 a32 a33 b1 b2 b3;A=a11,a12,a13;a21,a22,a23;a31,a32,a33;b=b1;b2;b3;X=linsolve(A,b)%调用linsolve函数求解XX=Ab%用左除运算求解方程组数值求解数值求解的命令如下：A=34,8,4;3,34,3;3,6,8;b=

50、4;6;2;X=linsolve(A,b)%调用linsolve函数求解Ab%用另一种方法求解X=ans=x=solve(1/(x+2)+4*x/(x2-4)=1+2/(x-2),x)f=sym(x-(x3-4*x-7)(1/3)=1);x=solve(f)x=solve(2*sin(3*x-pi/4)=1)x=solve(x+x*exp(x)-10,x)%仅标出方程的左端例例44 求非线性方程的解。x=1 x=3 x=5/36*pi x=一般格式ans=1/2/a*(-b+(b2-4*a*c)(1/2)1/2/a*(-b-(b2-4*a*c)(1/2)计算机格式例例45 求f=ax2+bx+

51、c 的解。f=a*x2+b*x+c;solve(f)例例46 求方程组的解。x y=solve(x+y-98,x(1/3)+y(1/3)-2,x,y)x y=solve(1/x3+1/y3=28,1/x+1/y=4,x,y)Warning:Explicit solution could not be found.In solve at 140 x=empty sym y=x=1 1/3 y=1/3 1 对第二个方程组Matlab给出了无解的结论，显然错误，请看完全与其同构的下述第一个方程组。输入命令如下：u,v=solve(u3+v3-98,u+v-2,u,v)x y=solve(x2+y2-

52、5,2*x2-3*x*y-2*y2)x=-1 1 2-2 y=2-2 1-1u=-3 5 v=5-3 例例47 解方程组 x+y+z=1 x-y+z=2 2x-y-z=1f=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)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:1x1 sym y:1x1 sym z:1x1 symx=2/3 y=-1/2z=5/6（二）（二）符号常微分方

53、程符号常微分方程的求解的求解dsolve(eqn1,eqn2,eqnN,condition1,conditionN,var1,varN)Matlab的符号运算工具箱中提供了功能强大的求解常微分方程（组）的函数dsolve。该函数的调用格式为：dsolve(eqn1,condition,var)该函数求解微分方程eqn1在初值条件condition下的特解。参数var描述方程中的自变量符号，省略时按缺省原则处理。若没有给出初值条件condition，则求方程的通解。函数求解微分方程组eqn1、eqnN在初值条件conditoion1、conditionN下的解，若不给出初值条件，则求方程组的通解

54、，var1、varN给出求解变量。在Matlab中，用大写字母D表示导数。例如，Dy表示y，D2y表示y，Dy(0)=5表示y(0)=5。D3y+D2y+Dy-x+5=0表示微分方程y+y+y-x+5=0。y=dsolve(Dy-(x2+y2)/x2/2,x)%方程的右端为0时可以不写y=dsolve(Dy*x2+2*x*y-exp(x),x)y=dsolve(Dy-x/y/sqrt(1-x2),x)例例48 求微分方程的通解。y=x*(-2+log(x)+C1)/(log(x)+C1)y=(exp(x)+C1)/x2 y=(-2*(1-x2)(1/2)+C1)(1/2)-(-2*(1-x2)

55、(1/2)+C1)(1/2)y=dsolve(Dy=2*x*y2,y(0)=1,x)y=dsolve(Dy-x2/(1+y2),y(2)=1,x)例例49 求微分方程的特解。y=-1/(-1+x2)y=-1/4*(4*x3-16+4*(20+x6-8*x3)(1/2)(1/3)+1/(4*x3-16+4*(20+x6-8*x3)(1/2)(1/3)-1/2*i*3(1/2)*(1/2*(4*x3-16+4*(20+x6-8*x3)(1/2)(1/3)+2/(4*x3-16+4*(20+x6-8*x3)(1/2)(1/3)1/2*(4*x3-16+4*(20+x6-8*x3)(1/2)(1/3)

56、-2/(4*x3-16+4*(20+x6-8*x3)(1/2)(1/3)-1/4*(4*x3-16+4*(20+x6-8*x3)(1/2)(1/3)+1/(4*x3-16+4*(20+x6-8*x3)(1/2)(1/3)+1/2*i*3(1/2)*(1/2*(4*x3-16+4*(20+x6-8*x3)(1/2)(1/3)+2/(4*x3-16+4*(20+x6-8*x3)(1/2)(1/3)y=dsolve(Dy+2*y/x-4*x,y(1)=2,x)%用符号方法得到方程的解析解。例例50 用微分方程的数值解法和符号解法解方程，并对结果进行比较。为了求方程的数值解，需要按要求建立一个函数文件

57、：function f=fxyy(x,y)f=(4*x2-2*y)/x;%只能是y=f(x,y)的形式，当不是这种形式时，要变形t,w=ode45(fxyy,1,2,2);%得到区间1，2中的数值解，以向量t、w存储y=(x4+1)/x2x=linspace(1,2,100);y=x.2+1./x.2;%为作图把符号解的结果离散化plot(x,y,b.,t,w,r-);为了对两种结果进行比较，在同一个坐标系中作出两种结果的图形。输入命令：x,y=dsolve(Dx=4*x-2*y,Dy=2*x-y,t)x,y=dsolve(D2x-y,D2y+x,t)例例51 求微分方程组的解。x=C1+C2

58、*exp(3*t)y=1/2*C2*exp(3*t)+2*C1 x=-C1*exp(1/2*2(1/2)*t)*sin(1/2*2(1/2)*t)-C2*exp(-1/2*2(1/2)*t)*sin(1/2*2(1/2)*t)+C3*exp(1/2*2(1/2)*t)*cos(1/2*2(1/2)*t)+C4*exp(-1/2*2(1/2)*t)*cos(1/2*2(1/2)*t)y=-C1*exp(1/2*2(1/2)*t)*cos(1/2*2(1/2)*t)+C2*exp(-1/2*2(1/2)*t)*cos(1/2*2(1/2)*t)-C3*exp(1/2*2(1/2)*t)*sin(1

59、/2*2(1/2)*t)+C4*exp(-1/2*2(1/2)*t)*sin(1/2*2(1/2)*t)x,y=dsolve(Dx=y,Dy=-x,x(0)=0,y(0)=1)dsolve(D2y=-a2*y,y(0)=1,Dy(pi/a)=0)x=sin(t)y=cos(t)ans=cos(a*t)y=dsolve(D2y+2*Dy+2*y=0,y(0)=1,Dy(0)=0)ans=exp(-t)*cos(t)+exp(-t)*sin(t)ezplot(y)%方程解y(t)的时间曲线图例例52 求下述方程的解例例53 求解两点边值问题：y=dsolve(x*D2y-3*Dy=x2,y(1)=

60、0,y(5)=0,x)ezplot(y,-1,6)hold onplot(1,5,0,0,.r,MarkerSize,20)text(1,1,y(1)=0)text(4,1,y(5)=0)title(x*D2y-3*Dy=x2,y(1)=0,y(5)=0)hold off （1）求解边值问题y=31/468*x4-1/3*x3+125/468（2）观察“解”的图形 两点边值问题的解曲线 五、五、符号数学的简易绘图函数符号数学的简易绘图函数 利用计算获得的符号符号表达式直接绘图绘图；据获得的符号表达式或符号数值结果，进而转换得到数值数据，再利用Matlab的数值绘图数值绘图指令绘制所需的图形符号

61、计算结果的可视化有两种方法：Matlab中有一组专门实现函数可视化的指令，它们的名称特点是：名称前两个字符都是“ez”。指令名含 义可 执 行 示 例ezcontour画等高线ezcontour(cos(x+sin(y)-sin(y),colormap(jet)ezcontourf画填色等高线colormap(flipud(cool),ezmesh(sin(x)*sin(y)hidden off,hold on,ezcontourf(sin(x)*sin(y),view(34,62),hold offezmesh画网线图ezmesh(exp(-s)*cos(t),exp(-s)*sin(t),

62、t,0,8,0,4*pi)ezmeshc画带等高线的网线图ezmeshc(y/(1+x2+y2),-5,5,-2*pi,2*pi)ezplot画二维曲线ezplot(1/y-log(y)+log(-1+y)+x-1)ezplot3画三维曲线ezplot3(sin(3*t)*cos(t),sin(3*t)*sin(t),t,animate)ezpolar画极坐标曲线ezpolar(sin(tan(t)ezsurf画曲面图ezsurf(x+8)*(y)2)/(x+8)2+(y)4+eps),circ)shading interp,colormap(flipud(hot),view(83,84)ez

63、surfc画带等高线的曲面图ezsurfc(sin(x)*sin(y)简捷作图指令汇总表（一）（一）二维绘图函数二维绘图函数ezplot(f)：绘制表达式f（x）的二维图形，轴坐标的近似范围为-2,2；ezplot(f,xmin,xmax)：绘制表达式f（x）的二维图形，轴坐标的近似范围为xmin,xmax；ezplot(f,xmin,xmax,ymin,ymax)：在指定x和y范围内，绘制y=f(x)描写的平面曲线。(1)对于函数对于函数f=f(x)，ezplot函数的调用格式为：函数的调用格式为：ezplot(f)：在默认区间-2x2和-2y2绘制f(x,y)=0的图形；ezplot(f,

64、xmin,xmax)：绘制f(x,y)=0的二维图形，轴坐标的近似范围为xmin,xmax；ezplot(f,xmin,xmax,ymin,ymax)：在指定x和y范围内，绘制f(x,y)=0的图形。(2)对于隐函数对于隐函数f=f(x,y)，ezplot函数的调用格式为：函数的调用格式为：ezplot(x,y):在默认区间0t2绘制x=x(t)和y=y(t)的图形ezplot(x,y,txmin,txmax)：在区间tmin t y=dsolve(Dy)2-x*Dy+y=0,x)clf,hold on,ezplot(y(1),-6,6,-4,8,1)cc=get(gca,Children);

65、set(cc,Color,r,LineWidth,5)for k=-2:0.5:2;ezplot(subs(y(2),C1,k),-6,6,-4,8,1);endhold off,title(fontname隶书fontsize16通解和奇解)y=1/4*x2 -C12+x*C1 通解和奇解曲线（二）（二）三维绘图函数三维绘图函数ezplot3(x,y,z)：绘制由表达式 x=x(t),y=y(t)和 z=z(t)定义的三维曲线，自变量 t 的变化范围为-2,2;ezplot3(x,y,z，tmin,tmax)：绘制由表达式 x=x(t),y=y(t)和 z=z(t)定义的三维曲线，自变量 t

66、 的变化范围为 tmin,tmax;ezplot3(,animate)：如果在函数中增加参数 animate，则绘制三维动态轨迹图。例例58 根据表达式 x=sin(t),y=cos(t)和 z=t,绘制三维曲线。syms t;ezplot3(sin(t),cos(t),t,0,6*pi,animate)syms t;ezplot3(sin(t),cos(t),t,0,6*pi)绘动态轨迹图（三）等高线绘图函数（三）等高线绘图函数ezcontour(f)：绘制由表达式f(x,y)定义的等高线，自变量 x和 y的变化范围为-2,2;ezcontour(f,domain)：绘制由表达式f(x,y)定义的等高线，自变量 x和 y的变化范围由domain确定。domain可以是41的矢量xmin,xmax,ymin,ymax；也可以是21的矢量min,max，此时，minxmax,minymax;ezcontour(,n)：按nn的网格密度绘制等高线图，n的缺省值为60。填充等高线的简易绘图函数为ezcontourf，其调用格式与ezcontour相同。例例59 根据表达式 绘制 f 的等高线