MATLAB中的符号运算.ppt

《MATLAB中的符号运算.ppt》由会员分享，可在线阅读，更多相关《MATLAB中的符号运算.ppt（61页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、MATLABMATLAB中的符号运算中的符号运算2004.8.42004.8.4MATLAB所具有的符号数学工具箱与其它所有工具不同，它适用于广泛的用途，而不是针对一些特殊专业或专业分支。另外，MATLAB符号数学工具箱与其它的工具箱区别还因为它使用字符串来进行符号分析，而不是基于数组的数值分析。符号数学工具箱是操作和解决符号表达式的符号数学工具箱(函数)集合，有复合、简化、微分、积分以及求解代数方程和微分方程的工具。另外还有一些用于线性代数的工具，求解逆、行列式、正则型式的精确结果，找出符号矩阵的特征值而无由数值计算引入的误差。工具箱还支持可变精度运算，即支持符号计算并能以指定的精度返回结果

2、。符号数学工具箱中的工具是建立在功能强大的Maple之上。它最初是由加拿大的滑铁卢（Waterloo）大学开发的。当要求MATLAB进行符号运算时，它就请求Maple去计算并将结果返回到MATLAB命令窗口。因此，在MATLAB中的符号运算是MATLAB处理数字的自然扩展。符号表达式是代表数字、函数、算子和变量的MATLAB字符串，或字符串数组。不要求变量有预先确定的值，符号方程式是含有等号的符号表达式。符号算术是使用已知的规则和给定符号恒等式求解这些符号方程的实践，它与代数和微积分所学到的求解方法完全一样。符号矩阵是数组，其元素是符号表达式。符号表达式1/(2*xn)y=1/sqrt(2*x

3、)cos(x2)-sin(2*x)M=sym(a，b；c，d)符号表达式MATLAB表达式f=int(x3/sqrt(1-x)，a，b)MATLAB符号函数可用多种方法来操作这些表达式diff(cos(x)ans=-sin(x)%differentiate cos(x)with respect to xM=sym(a，b；c，d )M=a，b c，d%create a symbolic matrix Mdeterm(M)ans=a*d-b*c%find the determinant of the symbolic matrix上面的第一个例子的符号表达式是用单引号以隐含方式定义的。它告诉MA

4、TLAB cos(x)是一个字符串并说明diff(cosx )是一个符号表达式而不是数字表达式；而在第二个例子中，用函数sym显式地告诉MATLAB M=sym(a，b；c，d )是一符号表达式。在MATLAB可以自己确定变量类型的场合下，通常不要求显式函数sym。符号变量 当字符表达式中含有多于一个的变量时，只有一个变量是独立变量。如果不告诉MATLAB哪一个变量是独立变量，MATLAB将基于以下规则选择一个：在符号表达式中缺省的独立变量是唯一的，除去i和j的小写字母，不是单词的一部分。如果没有这种字母，就选择x作为独立变量。如字符不是唯一的，就选择在字母顺序中最接近x的字母。如果有相连的字

5、母，就选择在字母表中较后的那一个。缺省的独立变量，有时称作自由变量，在表达式 1/(5+cos(x)中是 x x ；在 3*y+z 中是 y y ；在 a+sin(t)是 t t 。在表式 sin(pi/4)-cos(3/5)中自由符号变量是 x x ，因为此式是一个符号常数无符号变量。可利用函数symvar询问MATLAB在符号表达式中哪一个变量它认为是独立变量。symvar(a*x+y*)%find the default symbolic variableans=xsymvar(a*t+s/(u+3)%u is the closest to x ans=usymvar(sin(omega

6、)%omega is not a singlee character。ans=xsymvar(3*i+4*j )%i and j are equel to sqrt(-1)ans=xsymvar(y+3*s ，t )ans=s%find the variable closest to t rather than x如果利用规则symvarsymvar不能找到一个缺省独立变量，它便假定无独立变量并返回x。这一结论对含有由多个字母组成的变量，如：alpha或s2的表达式，或不含变量的符号常数均成立。如果需要，绝大多数命令都使用用户选项以指定独立变量。diff(xn )%differentiate

7、with respect to the default variable x ans=xn*n/x diff(xn ，n )%differentiate xn with respect to n ans=xn*log(x)diff(sin(omega)%differentiate using the default variables(x)ans=0 diff(sin(omega)，omega )%specify the independent variable ans=cos(omega)符号表达式运算 一旦创建了一个符号表达式，或许想以某些方式改变它；也许希望提取表达式的一部分，合并两个表

8、达式或求得表达的数值。有许多符号工具可以帮助完成这些任务。所有符号函数(很少特殊例外的情况)作用到符号表达式和符号数组，并返回符号表达式或数组。其结果有时可能看起来象一个数字，但事实上它是一个内部用字符串表示的一个符号表达式。可以运用MATLAB函数isstrisstr来找出像似数字的表达式是否真是一个整数或是一个字符串。如果表达式是一个有理分式(两个多项式之比)，或是可以展开为有理分式(包括哪些分母为1的分式)，可利用numdennumden来提取分子或分母。例如，给定如下的表达式：提取分子和分母在必要时，numdennumden将表达式合并、有理化并返回所得的分子和分母。进行这项运算的MA

9、TLAB语句是：m=x2 m=x2%create a simple expression m=m=x2x2nn，d=d=numdennumden(m)(m)%extract the numerator and denominatorn=x2n=x2d=1d=1f=a*x2/(b-x)f=a*x2/(b-x)%create a rational expressionf=a*x2/(b-x)f=a*x2/(b-x)n n，d=d=numdennumden(f)(f)%extract the numerator and denominatorn=n=a*x2a*x2d=d=b-xb-xg=3/2*x

10、2+2/3*x-3/5%rationalize and extract the partsg=3/2*x2+2/3*x-3/5n，d=numden(g)n=45*x2+20*x-18d=30h=(x2+3)/(2*x-1)+3*x/(x-1)%the sum of rational polynomialsh=(x2+3)/(2*x-1)+3*x/(x-1)n，d=numden(h)%rationalize and extractn=x3+5*x2-3d=(2*x-1)*(x-1)很多标准的代数运算可以在符号表达式上执行，函数symadd、symsub、symlnul和symdiv为加、减、乘、

11、除两个表达式，sympow将一个表达式上升为另一个表达式的幂次。标准代数运算例如：给定两个函数f=2*x2+3*x-5%definethesymbolicexpressionf=2*x2+3*x-5g=x2-x+7g=x2-x+7symadd(f，g)%findanexpressionforf+gans=3*x2+2*x+2symsub(f，g)%findanexpressionforf-gans=x2+4*x-12symmul(f，g)%findanexpressionforf*gans=(2*x2+3*x-5)*(x2-x+7)symdiv(f，g)%findanexpressionfor

12、f/gans=(2*x2+3*x-5)/(x2-x+7)sympow(f，3*x)%findanexpressionforans=(2*x2+3*x-5)3*变换函数函数symsym可获取一个数字参量并将其转换为符号表达式。函数numnericnumneric的功能正好相反，它把一个符号常数（无变量符号表达式）变换为一个数值。phi=(1+sqrt(5)/2%the golden ratiophi=(1+sqrt(5)/2%convert to a numeric valuenumeric(phi)ans=1.6180符号函数sym2polysym2poly将符号多项式变换成它的MATLAB等

13、价系数向量。函数poly2syrnpoly2syrn功能正好相反，并让用户指定用于所得结果表达式中的变量。f=2*x2+x3-3*x+5%f is the symbolic polynomialsf=2*x2+x3-3*x+5n=sym2poly(f)%extract eht numeric coefficient vectorn=1 2 -3 5poly2sym(n)%recreate the polynomials in x(the default)ans=2*x2+x3-3*x+5poly2sym(n，s)%recreate the polynomials in sans=s3+2*s2

14、-3*s+5变量替换 假设有一个以x为变量的符号表达式，并希望将变量转换为y。MATLAB提供一个工具称作subssubs，以便在符号表达式中进行变量替换。其格式为subssubs（f f，newnew，old)old)，其中f是符号表达式，newnew和oldold是字符、字符串或其它符号表达式。新字符串将代替表达式f中各个旧字符串。f=a*x2+b*x+c%create a function f(x)f=a*x2+b*x+csubs(f，s，x)%substitute s for x in the expression f ans=a*s2+b*s+c subs(f，alpha，a)%su

15、bstitute alpha for a in f ans=alpha*x2+b*x+cg=3*x2+5*x-4%create another functiong=3*x2+5*x-4h=subs(g，2，x)%substitute 2 for x in g h=18isstr(h)%show that the result is a symbolic expression ans=1微分和积分 微分和积分是微积分学研究和应用的核心，并广泛地用在许多工程学科。MATLAB符号工具能帮助解决许多这类问题。微分导数（包括偏导数）函数 diffdiff格式 diff(S,v)、diff(S,sym(

16、v)%对表达式S中指定符号变量v计算S的1阶导数。diff(S)%对表达式S中的符号变量v计算S的1阶导数，其中v=findsym(S)。diff(S,n)%对表达式S中的符号变量v计算S的n阶导数，其中v=findsym(S)。diff(S,v,n)%对表达式S中指定的符号变量v计算S的n阶导数。syms x y tD1=diff(sin(x2)*y2,2)%计算D2=diff(D1,y)%计算D3=diff(t6,6)计算结果为：D1=-4*sin(x2)*x2*y2+2*cos(x2)*y2D2=-8*sin(x2)*x2*y+4*cos(x2)*yD3=720符号表达式的微分以四种形式

17、利用函数diffdiff：f=a*x3+x2-b*x-c%define a symbolic expressionf=a*x3+x2-b*x-cdiff(f)%differentiate with respect to the default variable xans=3*a*x2+2*x-bdiff(f，a)%differentiate with respect to aans=x3diff(f，2)%differentiate twice with respect to xans=6*a*x+2diff(f，a，2)%differentiate twice with respect to

18、 aans=0函数diffdiff也可对数组进行运算。如果F是符号向量或数组，diff(F)diff(F)对数组内的各个元素进行微分。F=sym(a*x，b*x2；c*x3，d*s)%create a symbolic arrayF=a*x，b*x2c*x3，d*sdiff(F)%differentiate the element with respect to x ans=a，2*b*x 3*c*x2，0积分 积分函数intint(f)(f)，其中f是一符号表达式，它力图求出另一符号表达式F使diff(F)=fdiff(F)=f。正如从研究微分学所了解的，积分比微分复杂得多。积分或逆求导不一

19、定是以封闭形式存在，或许存在但软件也许找不到，或者软件可明显地求解，但超过内存或时间限制。当MATLAB不能找到逆导数时，它将返回未经计算的命令。函数 intint 格式 R=int(S,v)%对符号表达式S中指定的符号变量v计算不定积分。注意的是，表达式R只是函数S的一个原函数，后面没有带任意常数C。R=int(S)%对符号表达式S中的符号变量v计算不定积分，其中v=findsym(S)。R=int(S,v,a,b)%对表达式s中指定的符号变量v计算从a到b的定积分 R=int(S,a,b)%对符号表达式s中的符号变量v计算从a到b的定积分，其中v=findsym(S)。int(log(x)

20、/exp(x2)%attempt to integrateans=log(x)/exp(x2)同微分一样，积分函数有多种形式。形式intint(f)(f)相对于缺省的独立变量求逆导数；形式(f f，s)s)相对于符号变量s积分；形式intint(f(f，a a，b)b)和intint(f(f，s s，a a，b)b)，a a，b b是数值，求解符号表达式从a a到b b的定积分；形式intint(f(f，m m ，n)n)和形式intint(f(f，s s，m m，n)n)，其中mm，n n是符号变量，求解符号表达式从mm到n n的定积分。f=sin(s+2*x)%crate a symbol

21、ic functionf=sin(s+2*x)int(f)%integrate with respect to xans=-1/2*cos(s+2*x)int(f，s)%integrate with respect to sans=-cos(s+2*x)int(f，pi/2，pi)%integrate with respect to x from pi/2 to pi ans=-cos(x)int(f，s，pi/2，pi)%integrate with respect to s from/2 to ans=cos(2*x)-sin(2*x)int(f，m，n)%integrate with r

22、espect to x from m to nans=-1/2*cos(s+2*n)+1/2*cos(s+2*m)正如函数diffdiff一样，积分函数intint对符号数组的每一个元素进行运算。F=sym(a*x，b*x2；c*x3，d*s )%create a symbolic arrayF=a*x，b*x2c*x3，d*sdiff(F)%ubtegrate the array elements with respect to x ans=1/2*a*x2，1/3*b*x3 1/4*c*x4，d*s*xsyms x z t alphaINT1=int(-2*x/(1+x3)2)INT2=i

23、nt(x/(1+z2),z)INT3=int(INT2,x)INT4=int(x*log(1+x),0,1)INT5=int(2*x,sin(t),1)INT6=int(exp(t),exp(alpha*t)当xa时：limit(F,a)%用命令findsym(F)确定F中的自变量，设为变量x，再计算F的极限值；当xa时：limit(F)%用命令findsym(F)确定F中的自变量，设为变量x，再计算F的极限值；当x0时：limit(F,x,a,right)或limit(F,x,a,left)%计算符号函数F的单侧极限：左极限xa-或右极限xa+。极限函数 limitlimit格式 limit

24、(F,x,a)：%计算符号表达式F=F(x)的极限值，syms x a t h n;L1=limit(cos(x)-1)/x)L2=limit(1/x2,x,0,right)L3=limit(1/x,x,0,left)L4=limit(log(x+h)-log(x)/h,h,0)v=(1+a/x)x,exp(-x);L5=limit(v,x,inf,left)L6=limit(1+2/n)(3*n),n,inf)计算结果为：L1=0、L2=inf、L3=-inf、L4=1/xL5=exp(a),0、L6=exp(6)符号表达式画图 在许多的场合，将表达式可视化是有利的。MATLAB提供了函数e

25、zplotezplot来完成该任务。y=16*x2+64*x+96%expression to ploty=16*x2+64*x+96ezplot(y)ezplotezplot绘制了定义域为-2x2的给定符号函数，并相应地调整了y轴比例，还加了网格栅和标志。在这个例子中，我们感兴趣的时间是从0到6。让我们再试一下，并指定时间范围.ezplot(y，0 6)%plot y for 0 x6 符号表达式简化 有时MATLAB返回的符号表达式难以理解，有许多工具可以使表达式变得更易读懂。第一个就是函数prettypretty，该命令以类似于数学课本上的形式来显示符号表达式。f=taylor(log(

26、x+1)/(x-5)%6 terms is the defaultf=-1/5*x+3/50*x2-41/750*x3+293/7500*x4-1207/37500*x5+O(x6)pretty(f)2 41 3 293 4 1207 5 6 -1/5 x+3/50 x -x +-x -x +O(x)750 7500 37500符号表达式可用许多等价形式来提供。在不同的场合，某种形式可能胜于另一种。MATLAB用许多命令来简化或改变符号表达式。f=sym(x2-1)*(x-2)*(x-3)%create a functionf=(x2-1)*(x-2)*(x-3)collect(f)%coll

27、ect all like termsans=x4-5*x3+5*x2+5*x-6hornor(ans)%change to Horner or nested representationans=-6+(5+(5+(-5+x)*x)*x)*xfactor(ans)%express as a product of polynomialsans=(x-1)*(x-2)*(x-3)*(x+1)expand(f)%distribute products over sumsans=x4-5*x3+5*x2+5*x-6 SimplifySimplify是功能强大、通用的工具。它利用各种类型代数恒等式，包括求

28、和、积分和分数幂、三角、指数和log函数、Bessel函数、超几何函数和函数，来简化表达式。下面例子说明函数的乘幂。simplify(log(2*x/y)ans=log(2)+log(x)-log(y)simplify(sin(x)2+3*x+cos(x)2-5 )ans=3*x-4 simplify(-a2+1)/(1-a)ans=a+1方程求解 用MATLAB所具有的符号工具可以求解符号方程。有一些工具已经在前面介绍过，更多将在本节予以检验。求解单个代数方程：MATLAB具有求解符号表达式的工具。若表达式不是一个方程式(不含等号)，则在求解之前函数solvesolve将表达式置成等于0。s

29、olve(a*x2+b*x+c )%solve for the roots of the quadratic equtionans=1/2/a*(-b+(b2-4*a*c)1/2)1/2/a*(-b-(b2-4*a*c)1/2)结果是符号向量，其元素是方程的2个解。如果想对非缺省x变量求解，solvesolve必须指定变量。solve(a*x2+b*x+c ，b )%solve for bans=-(a*x2+c)/x 带有等号的符号方程也可以求解。f=solve(cos(x)=sin(x)%solve for xf=1/4*pit=solve(tan(2*x)=sin(x)t=0acos(1

30、/2+1/2*3(1/2)acos(1/2=1/2*3(1/2)并得到数值解。numeric(f)ans=0.7854numeric(t)ans=0 0+0.8314i1.9455注意在求解周期函数方程时，有无穷多的解。在这种情况下，solvesolve对解的搜索范围限制在接近于零的有限范围，并返回非唯一的解的子集。代数方程组求解可以同时求解若干代数方程，语句solve(s1solve(s1，s2s2，.，snsn)对缺省变量求解n n个方程，语句solve(s1solve(s1，s2s2，.，snsn，v1 v1，v2v2，.，vnvn)对n n个个 v1v1，v2v2，.vnvn 的未知数

31、求解n n个方程。如果电影票价为3.00美元，爆米花为1.00美元，糖棒为50美分，她有足够的钱去买这三样东西？如何处理中小学典型的代数问题？黛安娜(Diane)想去看电影，她从小猪存钱罐倒出硬币并清点，她发现：10美分的硬币数加上5美分的硬币总数的一半等于25美分的硬币数。1美分的硬币数比5美分、10美分以及25美分的硬币总数多10。25美分和10美分的硬币总数等于1美分的硬币数加上1/4的5美分的硬币数25美分的硬币数和1美分的硬币数比5美分的硬币数加上8倍的10美分的硬币数多1。首先，根据以上给出的信息列出一组线性方程，假如p，n，d和q分别表示1美分，5美分，10美分，和25美分的硬币

32、数然后，建立MATLAB符号方程并对变量求解。eq1=d+(n+p)/2=q ；eq2=p=n+d+q-10 ；eq3=q+d=p+n/4 ；eq4=q+p=n+8*d-1 ；pennies，nickles，dimes，quarters=solve(equ1，equ2，equ3，equ4，p，n，d，q )pennies=16 nickles=8 dimes=3 quarters=15所以，黛安娜有16枚1美分的硬币，8枚5美分的硬币，3枚10美分的硬币，15枚25美分的硬币，这就意味着money=.01*16+.05*8+.10*3+.25*15money=4.6100她就有足够的钱去买电影

33、票，爆米花和糖棒并剩余11美分。单个微分方程 常微分方程有时很难求解，MATLAB提供了功能强大的工具，可以帮助求解微分方程。函数dsovledsovle计算常微分方程的符号解。因为我们要求解微分方程，就需要用一种方法将微分包含在表达式中。所以，dsovledsovle句法与大多数其它函数有一些不同，用字母D D来表示求微分，D2D2，D3D3等等表示重复求微分，并以此来设定方程。任何D D后所跟的字母为因变量。方程=0用符号表达式D2y=0来表示。独立变量可以指定或由symvarsymvar规则选定为缺省。函数 dsolvedsolve格式 r=dsolve(eq1,eq2,cond1,co

34、nd2,v)说明:对给定的常微分方程（组）eq1,eq2,中指定的符号自变量v，与给定的边界条件和初始条件cond1,cond2,.求符号解（即解析解）r；若没有指定变量v，则缺省变量为t；在微分方程（组）的表达式eq中，大写字母D表示对自变量(设为x)的微分算子：D=d/dx，D2=d2/dx2，。微分算子D后面的字母则表示为因变量，即待求解的未知函数。，初始和边界条件由字符串表示：y(a)=b，Dy(c)=d，D2y(e)=f，等等，分别表示 ，；若边界条件少于方程（组）的阶数，则返回的结果r中会出现任意常数C1，C2，；dsolve命令最多可以接受12个输入参量（包括方程组与定解条件个数

35、，当然我们可以做到输入的方程个数多于12个，只要将多个方程置于一字符串内即可）。若没有给定输出参量，则在命令窗口显示解列表。若该命令找不到解析解，则返回一警告信息，同时返回一空的sym对象。这时，用户可以用命令ode23或ode45求解方程组的数值解。例如，一阶方程dy/dx=1+y2的通解为：dsolve(Dy=1+y2 )%find the general solutionans=-tan(-x+C1)其中，C1C1是积分常数。求解初值y(0)=1的同一个方程就可产生：dsolve(Dy=1+y2，y(0)=1)%add an initial conditiony=tan(x+1/4*pi

36、)独立变量可用如下形式指定：dsolve(Dy=1+y2，y(0)=1，v)%find solution to dy/dvans=tan(v+1/4*pi)让我们举一个二阶微分方程的例子，该方程有两个初始条件：=cos(2x)-y (0)=0 y(0)=1y=dsolve(D2y=cos(2*x)-y，Dy(0)=0，y(0)=1)y=-2/3*cos(x)2+1/3+4/3*cos(x)y=simple(y)%y looks like it can be simplifiedy=-1/3*cos(2*x)+4/3*cos(x)通常，要求解的微分方程含有一阶以上的项，并以下述的形式表示：-2

37、-3y=0通解为：y=solve(D2y-2Dy-3*y=0)y=C1*exp(-x)+C2*exp(3*x)加上初始条件：y(0)=0和y(1)=1可得到：y=solve(D2y-2Dy-3*y=0 ，y(0)=0，y(1)=1 )y=1/(exp(-1)-exp(3)*exp(-x)-1/(exp(-1)-exp(3)*exp(3*x)y=simple(y)%this looks like a candidate for simplificationy=-(exp(-x)-exp(3*x)/(exp(3)-exp(-1)pretty(y)%pretty it up exp(-x)-exp(

38、3 x)-exp(3)-exp(-1)微分方程组：函数dsolvedsolve也可同时处理若干个微分方程式，下面有两个线性一阶方程。=3f+4g=-4f+3g通解为：f=exp(3*x)*sin(4*x)g=exp(3*x)*cos(4*x)f，g=dsolve(Df=3*f+4*g ，Dg=-4*f+3*g )f=C1*exp(3*x)*sin(4*x)+C2*exp(3*x)*cos(4*x)g=-C2*exp(3*x)*sin(4*x)+C1*exp(3*x)*cos(4*x)加上初始条件：f(0)=0和g(0)=1，我们可以得到：f，g=dsolve(Df=3*f+4*g ，Dg=-4*f+3*g ，f(0)=0，g(0)=1 )D1=dsolve(D2y Dy=exp(x)D2=dsolve(t*D2f=Df*log(Dy)/t)D3=dsolve(Dy)2+y2=1,s)D4=dsolve(Dy=a*y,y(0)=b)%带一个定解条件D5=dsolve(D2y=-a2*y,y(0)=1,Dy(pi/a)=0)%带两个定解条件x,y=dsolve(Dx=y,Dy=-x)%求解线性微分方程组u,v=dsolve(Du=u+v,Dv=u-v)