1、第三讲 MATLAB 的符号运算 matlab 不仅具有数值运算功能，还开发了在 matlab 环境下实现符号计算的工具包Symbolic Math Toolbox 符号运算的功能 符号表达式、符号矩阵的创建 符号线性代数 因式分解、展开和简化 符号代数方程求解 符号微积分 符号微分方程一、创建符号变量1. 什么是符号运算 与数值运算的区别 数值运算中必须先对变量赋值，然后才能参与运算。 符号运算无须事先对独立变量赋值，运算结果以标准的符号形式表达。 特点： 运算对象可以是没赋值的符号变量 可以获得任意精度的解 Symbolic Math Toolbox符号运算工具包通过调用 Maple 软件
2、实现符号计算的。 maple 软件主要功能是符号运算，它占据符号软件的主导地位。2. Sym 函数定义符号变量(1) S=sym(arg)Construct symbolic numbers, variables and objects.S = SYM(A) constructs an object S, of class sym, from A.If the input argument is a string, the result is a symbolic number or variable.If the input argument is a numeric scalar or m
3、atrix,the result is a symbolic representation of the given numeric valuesx = sym(x) creates the symbolic variable with name x and stores the result in x. x = sym(x,real) also assumes that x is real, so that conj(x) is equal to x. alpha = sym(alpha) and r = sym ( Rho , real) are other examples. Simil
4、arly, k =sym(k,positive) makes k a positive (real) variable. x = sym(x,unreal) makes x a purely formal variable with no additional properties (i.e., insures that x is NEITHER real NOR positive). See also: SYMS.Statements like pi = sym(pi) and delta = sym(1/10) create symbolic numbers which avoid the
5、 floating point approximations inherent in the values of pi and 1/10. The pi created in this way temporarily replacesthe built-in numeric function with the same name.S = sym(A,flag) converts a numeric scalar or matrix to symbolic form.The technique for converting floating point numbers is specified
6、by the optional second argument, which may be f, r, e or d. The default is r.f stands for floating point. All values are represented in the form 1.F*2(e) or -1.F*2(e) where F is a string of 13 hexadecimal digits and e is an integer. This captures the floating point values exactly, but may not be con
7、venient for subsequent manipulation.For example, sym(1/10,f) is 1.999999999999a*2(-4) because 1/10 cannot be represented exactly in floating point.r stands for rational. Floating point numbers obtained by evaluating expressions of the form p/q, p*pi/q, sqrt(p), 2q and 10q for modest sized integers p
8、 and q are converted to the corresponding symbolic form. This effectively compensates for the roundoff error involved in the original evaluation, but may not represent the floating point value precisely. If no simple rational approximation can be found, an expression of the form p*2q with large inte
9、gers p and q reproduces the floating point value exactly. For example, sym(4/3,r) is 4/3, but sym(1+sqrt(5),r) is 286977268806824*2(-51)e stands for estimate error. The r form is supplemented by a term involving the variable eps which estimates the difference between the theoretical rational express
10、ion and its actual floating point value. For example, sym(3*pi/4) is 3*pi/4-103*eps/249.d stands for decimal. The number of digits is taken from the current setting of DIGITS used by VPA. Fewer than 16 digits looses some accuracy, while more than 16 digits may not be warranted.For example, with digi
11、ts(10), sym(4/3,d) is 1.333333333, while with digits(20), sym(4/3,d) is 1.3333333333333332593,which does not end in a string of 3s, but is an accurate decimal representation of the floating point number nearest to 4/3.3.syms 函数定义符号变量Short-cut for constructing symbolic objects.SYMS arg1 arg2 .is shor
12、t-hand notation forarg1 = sym(arg1);arg2 = sym(arg2); .SYMS arg1 arg2 . realis short-hand notation forarg1 = sym(arg1,real);arg2 = sym(arg2,real); .SYMS arg1 arg2 . positiveis short-hand notation forarg1 = sym(arg1,positive);arg2 = sym(arg2,positive); .SYMS arg1 arg2 . unrealis short-hand notation f
13、orarg1 = sym(arg1,unreal);arg2 = sym(arg2,unreal); .Each input argument must begin with a letter and must contain onlyalphanumeric characters.By itself, SYMS lists the symbolic objects in the workspace.Examples:syms x beta realis equivalent to:x = sym(x,real);beta = sym(beta,real);syms k positiveis
14、equivalent to:k = sym(k,positive);To clear the symbolic objects x and beta of real or positive status, typesyms x beta unreal二、创建符号 1、创建符号表达式和符号方程 2、创建符号矩阵 3、数字矩阵和符号矩阵的转换 4、符号矩阵的引用和修改 5、建立符号数学函数 6、数据类型间的相互转换1、创建符号表达式和符号方程 （1） 采用 sym 函数 f=sym(a*x2+b*x+c) 符号表达式 f=sym(a*x2+b*x+c0) 符号方程 （2）直接法 f= a*x2+b
15、*x+c 符号表达式 f= a*x2+b*x+c0 符号方程2、创建符号矩阵 （1）由 sym 函数直接创建 A=sym(4+x, x2, x; x3,5*x-3,x*a)或者 syms x,aA=4+x, x2, x; x3,5*x-3,x*a（2） 由字符串直接创建a=4+x x2 x ;x3 5*x-3 x*a注意：字符串的长度相等，否则出错3、数字矩阵和符号矩阵的转换 用 sym 函数 例： A=2/5,4/0.78, sqrt(23)/3;0.33,0.333,log(4) A = 0.4000 5.1282 1.5986 0.3300 0.3330 1.3863 FA=sym(A)
16、 FA = 2/5, 200/39, sqrt(23/9) 33/100, 333/1000, 6243314768165359*2(-52)4、符号矩阵的引用和修改 在符号计算中，引用和修改只能对符号矩阵的元素一个一个地进行 例：syms a b c d A=a,a+c,d+b;c,d,a+c;a+c+d,c,c+d*a; syms eee ddd A(1，3）=eee A(3,2)=ddd5、建立符号数学函数 （1） 建立一般数学函数 sym x y z Fsin(x+y)/(x-y) G=sqrt(x2+y2+z2) 建立 m 文件 （2） 建立抽象的数学函数 sym var1 var2
17、 var3 F=sym(f(var1 var2 var3) 例： 建立一阶差分函数 hxffd)( Syms x h F=sym(f(x) Df=(subs(f, x,x+h)-f)/h6、数据类型之间的相互转换 数值型、字符型、符号型数据类型，符号型级别最高，数值型级别最低 （1）转换为数值变量 X=double(S) 将符号变量变为数值变量 （2）转换为符号变量 S=sym(f) 对变量 f 没有限制符号矩阵与数值矩阵的转换 将数值矩阵转化为符号矩阵函数调用格式：sym(A)A=1/3,2.5;1/0.7,2/5A =0.3333 2.50001.4286 0.4000sym(A)ans
18、= 1/3, 5/210/7, 2/5 将符号矩阵转化为数值矩阵函数调用格式： numeric(A)A = 1/3, 5/210/7, 2/5 numeric(A)ans =0.3333 2.50001.4286 0.4000三、符号运算1. 符号矩阵运算所有涉及符号运算的操作都有专用函数来进行，5.0 以上也可以用运算符，矩阵运算操作指令都比较直观、简单。例如：a=b+c; a=a*b ；A=2*a2+3*a-5 等。符号矩阵运算的函数：symadd(a,d) 符号矩阵的加symsub(a,b) 符号矩阵的减symmul(a,b) 符号矩阵的乘symdiv(a,b) 符号矩阵的除sympow
19、(a,b) 符号矩阵的幂运算symop(a,b) 符号矩阵的综合运算例 1：f= 2*x2+3*x-5; g= x2+x-7;h= symadd(f,g)h= 3*x2+4*x-12例 2：f=cos(x);g= sin(2*x);symop(f,/,g,+,f,*,g) ans =cos(x)/sin(2*x)+cos(x)*sin(2*x)例 1：f= 2*x2+3*x-5; g= x2+x-7; syms x f=2*x2+3*x-5; g= x2+x-7; h=f+gh = 3*x2+4*x-12例 2：f=cos(x);g= sin(2*x); syms x f=cos(x);g=s
20、in(2*x); f/g+f*gans =cos(x)/sin(x)+cos(x)*sin(x) 符号运算函数：symsize 求符号矩阵维数charploy 特征多项式determ 符号矩阵行列式的值eigensys 特征值和特征向量inverse 逆矩阵transpose 矩阵的转置jordan 约当标准型simple 符号矩阵简化2. 任意精度的数学运算 在 symbolic 中有三种不同的算术运算：1. 数值类型 matlab 的浮点算术运算2. 有理数类型 maple 的精确符号运算3. vpa 类型 maple 的任意精度算术运算 浮点算术运算1/2+1/3 (定义输出格式 for
21、mat long)ans =0.83333333333333 符号运算sym(1/2)+(1/3)ans =5/6 精确解 任意精度算术运算digits(n) 设置可变精度，缺省 16 位vpa(x,n) 显示可变精度计算digits(25)vpa(1/2+1/3)ans =.8333333333333333333333333vpa(5/6,40) ans =.8333333333333333333333333333333333333333a=sym(1/4,exp(1);log(3),3/7)a = 1/4,exp(1)log(3), 3/7vpa(a,10)ans =.2500000000
22、, 2.7182818281.098612289, .42857142863. 符号微积分与积分变换 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（被积表达式， 积分变量 ， 积分上限， 积分下限） 定积分积分上限与积分下限缺省时为不定积分mtaylor(f,n) 泰勒级数展开ztrans(f) Z 变换Invztrans(f
23、) 反 Z 变换Laplace(f) 拉氏变换Invlaplace(f) 反拉氏变换fourier(f) 付氏变换Invfourier(f) 反付氏变换例 1.计算二重不定积分 dxyeF=int(int(x*exp(-x*y),x),y)F=1/y*exp(-x*y)例 2.计算 f=x*exp(-x*10)的 Z 变换F=ztrans(f)F=z*exp(-10)/(z-exp(-10)2 syms x y F=int(int(x*exp(-x*y),x),y)F =1/y*exp(-x*y) syms x f=x*exp(-x*10); F=ztrans(f) F=ztrans(x*ex
24、p(-x*10);F =z*exp(-10)/(z-exp(-10)2例 3. 计算指数函数 eAt。用拉氏反变换法计算 eAt 的公式为：eAt = L-1(SI-A)-1系统矩阵 3210A结果： ttttt eee22 a=0 1;-2 -3; syms s b=(s*eye(2)-a)b = s, -1 2, s+3 B=inv(b) (s+3)/(s2+3*s+2), 1/(s2+3*s+2) -2/(s2+3*s+2), s/(s2+3*s+2) b11=ilaplace(B(1,1);b(1,1)=b11; b12=ilaplace(B(1,2);b(1,2)=b12; b21=
25、ilaplace(B(2,1);b(2,1)=b21; b22=ilaplace(B(2,2);b(2,2)=b22; bb = -exp(-2*t)+2*exp(-t), exp(-t)-exp(-2*t) -2*exp(-t)+2*exp(-2*t), 2*exp(-2*t)-exp(-t)4.符号代数方程求解matlab 符号运算能够解一般的线性方程、非线性方程及一般的代数方程、代数方程组。当方程组不存在符号解时，又无其他自由参数，则给出数值解。命令格式：solve(f) 求一个方程的解Solve(f1,f2, fn) 求 n 个方程的解例 1. f = ax2+bx+c 求解f=a*x
26、2+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)以上为计算机格式此式为一般格式acb24 solve(f , b ) 对指定变量 b 求解ans =-(a*x2+c)/x例 2. 符号方程 cos(x)=sin(x) ；tan(2*x)=sin(x)求解f1=solve(cos(x)=sin(x),f1 =1/4*pif2=solve(tan(2*x)=sin(x)f2 = matlab4.2 的解 0acos(1/2+1/2*3(1/2)acos(1/2 -1/2*3(1/2
27、)f3= matlab5.3 的解 0 pi atan(1/2*(-2*3(1/2)(1/2),1/2+1/2*3(1/2) atan(-1/2*(-2*3(1/2)(1/2),1/2+1/2*3(1/2) atan(1/2*2(1/2)*3(1/4)/(1/2-1/2*3(1/2)+pi -atan(1/2*2(1/2)*3(1/4)/(1/2-1/2*3(1/2)-pi numeric(f2)ans =0 0 + 0.8314i1.9455 numeric(f3) matlab4.2 与 5.3 的对比ans =0 3.1416 0 + 0.8314i0 - 0.8314i1.9455 -
28、1.9455 例 3. 解方程组 x+y+z=1 x-y+z=22x-y-z=1g1=x+y+z=1,g2=x-y+z=2,g3=2*x-y-z=1f=solve(g1,g2,g3)f=solve(x+y+z=1,x-y+z=2,2*x-y-z=1)f =f.z = 5/6,f. y = -1/2, f.x = 2/3f=solve(x+y+z=1,x-y+z=2,2*x-y-z=1)f = x: 1x1 sym f.xans =2/3y: 1x1 sym f.yans =-1/2z: 1x1 sym f.zans =5/6x,y,z=solve(x+y+z=1,x-y+z=2,2*x-y-z
29、=1)x = 2/3 y =-1/2z =5/65. 符号微分方程求解 用一个函数可以方便地得到微分方程的符号解符号微分方程求解指令：dsolve命令格式：dsolve(f,g) f 微分方程，可多至 12 个微分方程的求解；g 为初始条件 默认自变量为 x,可任意指定自变量t, u等 微分方程的各阶导数项以大写字母 D 表示 yydxt的 一 阶 导 数或t 222的 二 阶 导 数或 Dnynydxtn阶 导 数的或y1,y2=dsolve(x1,x2,xn) 返回微分方程的解一阶微分方程dsolve(Dx=y,Dy=x,x(0)=0,y(0)=1)ans =x(t) = sin(t),
30、y(t) = cos(t)二阶微分方程dsolve(D2y=-a2*y,y(0)=1,Dy(pi/a)=0)ans =cos(a*x)例 3. 求该方程的解0)(;1)0(;22 dxyydxyy=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)的时间曲线图-6 -4 -2 0 2 4050100150200250300xexp(-x)*cos(x)+exp(-x)*sin(x) İ s (t)四、maple 函数符号运算的扩展maple是专门进行数学运算的软件工具
31、，具有超强的符号运算能力，提供了几乎包括所有数学领域的专用函数matlab依赖于 maple 的内核与函数库，扩展了自己的符号运算功能。matlab 还设计了对 maple 库函数的调用功能使得已有的 maple 数学功能，可以扩充 matlab 中,作为自身符号运算能力的扩展。1. maple 内核访问函数可以访问 maple 内核的 matlab 函数:maple 访问 maple 内核函数mapleinit maple 函数初始化mpa maple 函数定义mhelp maple 函数帮助命令procread maple 函数程序安装 maple 的调用格式maple(表达式） 将表达式
32、送至 maple 内核，返回符号表达式结果。maple （ 函数 ，变量 1，变量 2）调用 maple 函数，传递给定变量。例 1. 展开 5 阶 bernoulli 多项式，计算 x=3 时 bernoulli 数。a=maple(bernoulli(5,x)a =-1/6*x+5/3*x3+x5-5/2*x4a=maple(bernoulli(5,3)a =85例 2. 化简三角函数式 sin2x+cos2xa=maple(simplify(sin(x)2+cos(x)2);)a =1例 4. 求 f(t)=e-3tsint 的拉式变换f=maple(laplace(exp(-3*t)*
33、sin(t),t,s);)f =1/(s+3)2+1)例 4. 寻找二次多项式的完全平方f (x) = x2+2x+2a=maple(completesquare(x2+2*x+2)a =completesquare(x2+2*x+2) maple(with(student);) 将工具包装入内存a=maple(completesquare(x2+2*x+2)a =(x+1)2+1 maple 软件中的所有函数，在初始化时并没有完全装入内存，可用 readlib 指令把库函数读入内存，或用 with 指令将应用工具包装入内存。 调用格式maple(readlib(函数名);)maple(wit
34、h(工具包名);)例 5.求 sin(x2+y2)在 x=0,y=0 处泰勒级数展开式，8 阶截断。maple(mtaylor(sin(x2+y2),x=0,y=0,8)ans =mtaylor(sin(x2+y2),x = 0, y = 0,8)maple(readlib(mtaylor);)maple(mtaylor(sin(x2+y2),x=0,y=0,8)ans =x2+y2-1/6*x6-1/2*y2*x4-1/2*y4*x2-1/6*y62. mpa maple 变量定义 任何一个 matlab 定义的函数 f，可使用 mpa 语句直接调用，还可把 f 定义成 maple变量 v。
35、 maple 的工作空间与 matlab 工作空间是相互独立的， 所以 f 与 v 是属于不同工作空间中的变量 mpa 的调用格式：mpa(v,f)mpa v ff 为 matlab 工作空间中已存在的变量例. 电磁力计算公式为 204)(),(xNISIF试 I=0.5，x=0.1 邻域展开泰勒级数，3 阶截断，令常数 ，420NSk2),(xkIF1.直接调用maple(readlib(mtaylor);)maple(mtaylor(k*I2/x2,I=0.5,x=0.1,3);)2.定义符号函数 f(matlab6.1 无 map 函数)f=k*I2/x2;maple(mtaylor(f
36、,I=0.5,x=0.1,3);)ans =mtaylor(f,I = .5, x = .1,3)mpa(u,f)maple(mtaylor(u,I=0.5,x=0.1,3);)ans =25.*k-.50e3*k*(x-.1)+.10e3*k*(I-.5)+7500.000000000000*k*(x-.1)2+.1e3*k*(I-.5)2-.20e4*k*(I-.5)*(x-.1)注意：matlab 符号运算时，可以识别 matlab 定义的符号变量，但在调用 maple 函数时，需将 matlab 变量定义为 maple 变量后，所调用的函数方可识别和执行3.mhelp maple 函数
37、帮助命令 mhelp 是协助检索 maple 库函数的专用命令调用格式：mhelp 相关词条（工具词条，函数词条）例如：mhelp intro maple 介绍mhelp maple maple 命令格式mhelp tutorial maple 入门mhelp index maple 检索mhelp index 用于工具包检索library maple 标准库函数packages 应用工具包libmisc 其它库函数statements maple 语句描述expressions maple 表达式datatypes maple 数据格式tables maple 表格和阵列procedures
38、 maple 程序misc maple 其它应用一般帮助文本主要包括以下部分 FUNCTION 函数功能说明 CALLING SEQUENCE 调用格式 PARAMETERS 调用参数说明 SYNOPSIS 语法说明 EXAMPLES 应用举例 SEE ALSO 相关词条4.maple 库函数maple 库函数共分四类maple 内部函数：驻留函数任何条件下都可调用mhelp indexinternal maple 的外部函数读库定义部分：调用时先执行读库命令，因此与内部函数一样可直接调用mhelp indexexternal maple 的外部函数读库装入部分maple 其余外部函数需要在使用前执行 maple(readlib(函数名);)命令将其装入内存mhelp indexlibmisc maple 的惰性函数不能直接调用，还需一些函数如 mod,evala,evalf 等才能调用mhelp indexintert 小 结本节介绍了 matlab 语言的符号运算功能，通过学习应该掌握： 掌握如何创建、修改符号矩阵 掌握符号运算功能 maple 函数调用 mhelp 检索
