1、第五章 符号计算,符号计算的特点: 一、运算以推理解析的方式进行,因此不受计算误差积累问题困扰; 二、符号计算,或给出完全正确的封闭解,或给出任意精度的数值解(当封闭解不存在时); 三、符号计算指令的调用比较简单,与经典教科书公式相近; 四、计算所需时间较长,有时难以忍受。,定义基本符号对象的指令有两个:sym,syms 格式如下: f=sym(arg) 把数字、字符串或表达式arg转换为符号对象f。 f=sym(argn,flagn)把数值或数值表达式argn转换为flagn格式的符号对象。 argv=sym(argv,flagv)按flagv指定的要求把字符串argv定义为符号对象argv
2、。,5.1 符号对象和符号表达式 5.1.1 符号对象的生成和使用,说明: f=sym(argn,flagn)中的argn是数值或数值表达式时,flagn可取以下选项: d最接近的十进制浮点精确表示; e带(数值计算时)估计误差的有理表示; f十六进制浮点表示 r最接近有理表示,这是缺省设置。 两个正整数p,q构成的p/q,p*pi/q,sqrt(p),2q,10q,syms(argv1,argv2,argvk) 把字符argv1,argv2,argvk定义为基本符号对象。 syms argv1 argv2 argvk上述格式的简洁形式。,argv=sym(argv,flagv)中的argv是
3、字符时,flagv可取以下“限定性”选项: positive 限定argv为“正、实”符号变量。 real 限定argv为“实”符号变量。 unreal argv为非符号变量。,符号常数形成中的差异 a1=1/3,pi/7,sqrt(5),pi+sqrt(5) a2=sym(1/3,pi/7,sqrt(5),pi+sqrt(5) a3=sym(1/3,pi/7,sqrt(5),pi+sqrt(5),e) a4=sym(1/3,pi/7,sqrt(5),pi+sqrt(5) a24=a2-a4,几种输入下产生矩阵的异同 a1=sym(1/3,0.2+sqrt(2),pi) a2=sym(1/3,
4、0.2+sqrt(2),pi) a3=sym(1/3 0.2+sqrt(2) pi) a1_a2=a1-a2,把字符表达式转换为符号变量 y=sym(2*sin(x)*cos(x) y=simple(y),用符号计算验证三角等式 sin(fai1)cos(fai2)-cos(fai1)sin(fai2)=sin(fai1-fai2) syms fai1 fai2; y=simple(sin(fai1)*cos(fai2)-cos(fai1)*sin(fai2),例 求矩阵A=a11,a12;a21,a22的行列式值、逆和特征根syms a11 a12 a21 a22;A=a11,a12;a21
5、,a22 DA=det(A),IA=inv(A),EA=eig(A),syms A t tao w; yf=int(A*exp(-i*w*t),t,-tao/2,tao/2); Yf=simple(yf),5.1.2 符号计算中的算符和基本函数,由于matlab采用了重载技术,使得用来构成符号计算表达式的算符和基本函数,无论在形状、名称上、还是在使用方法上,都与数值计算中的算符和基本函数几乎完全相同,5.1.3 识别对象类别的指令数据对象及其识别指令的使用 (1) clear,a=1;b=2;c=3;d=4; Mn=a,b;c,d Mc=a,b;c,d Ms=sym(Mc) (2) SizeM
6、n=size(Mn),SizeMc=size(Mc),SizeMs=size(Ms),(3) CMn=class(Mn); CMc=class(Mc); CMs=class(Ms) (4)用isa判断每种矩阵的类别(若返回1,表示判断正确) isa(Mn,double); isa(Mc,char); isa(Ms,sym) (5) whos Mn Mc Ms,5.1.4 符号表达式中自由变量的确定,findsym(EXPR) 确认表达式EXPR中所有“自由”符号“变量”; findsym(EXPR,N) 从表达式EXPR中确认出靠x最近的N个独立自变量。 说明 EXPR可以是符号矩阵。此时,该
7、指令对自由变量的确认是对整个矩阵进行的,而不是对矩阵元素逐个进行的。,5.1.4 符号表达式中自由变量的确定,对独立自由符号变量的自动辨认。 (1)生成符号变量 syms a b x X Y; k=sym(3);z=sym(c*sqrt(delta)+y*sin(theta); EXPR=a*z*X+(b*x2+k)*Y; (2)找出EXPR中的全部自由符号变量 findsym(EXPR),(3)在EXPR中确定一个自由符号变量 findsym(EXPR,1) (4)在EXPR中确定2个和3个自由变量时的执行情况 findsym(EXPR,2),findsym(EXPR,3)findsym确定
8、自由变量是对整个矩阵进行的 syms a b t u v x y; A=a+b*x,sin(t)+u;x*exp(-t),log(y)+v findsym(A,3),5.2 符号表达式和符号函数的操作 5.2.1 符号表达式的操作,collect(EXPR,v)对EXPR表达式中指定的符号对象v的同幂项系数进行合并 expand(EXPR)对EXPR表达式进行多项式、三角函数、指数对数函数等展开 factor(EXPR)对EXPR表达式(或正整数)进行因式(或因子)分解 horner(EXPR)对多项式EXPR分解成嵌套形式 n,d=numden(EXPR)提取表达式EXPR的最小分母公因式d
9、和相应的分子多项式n,simplify(EXPR) 运用多种恒等式转换对表达式EXPR进行综合简化; simple(EXPR) 运用包括simplify在内的各种指令把EXPR转换成最简短形式 pretty(EXPR) 以习惯的“书写”方式显示EXPR表达式 说明 以上指令中的EXPR也可以是符号矩阵。在这种情况下,这些指令将对该矩阵的元素逐个进行操作。,按不同的方式合并同幂项 EXPR=sym(x2+x*exp(-t)+1)*(x+exp(-t); expr1=collect(EXPR) expr2=collect(EXPR,exp(-t),factor 指令的使用,(1)除x外不含其他自由
10、变量的情况 syms a x; f1=x4-5*x3+5*x2+5*x-6;factor(f1) (2)含其他自由变量的情况 f2=x2-a2;factor(f2) (3)对正整数的质数分解 factor(1025),对多项式进行嵌套型分解 clear;syms a x; f1=x4-5*x3+5*x2+5*x-6;horner(f1),(1)求矩阵各元素的分子、分母多项式 syms x; A=3/2,(x2+3)/(2*x-1)+3*x/(x-1);4/x2,3*x+4; n,d=numden(A) pretty(simplify(A) (2) pretty(simplify(n./d),(
11、1)运用simplify简化 syms x;f=(1/x3+6/x2+12/x+8)(1/3); sfy1=simplify(f),sfy2=simplify(sfy1) (2)运用simple简化 g1=simple(f),g2=simple(g1); 说明 simple给出的简化式比simplify给出的更短;多次使用simple指令,可找到最少字母的简化式,syms x;ff=cos(x)+sqrt(-sin(x)2); ssfy1=simplify(ff),ssfy2=simplify(ssfy1) gg1=simple(ff),gg2=simple(gg1),5.2.2 符号函数的求
12、反和复合,g=finerse(f,v)对指定自变量为v的函数f(v),求反函数 g=finverse(f)对缺省自变量求反函数g; fg=compose(f,g,v,w,t)对f(v)和v=g(w)求复合函数fg=f(g(w)|w=t fg=compose(f,g)对f(.)和v=g(.)求复合函数fg=f(g(.)。,说明 v,w分别是对函数f和g指定的自变量。t是所得复合函数得自变量符号 finverse(f,g)是finverse(f,v)的缺省形式。此时默认自变量由findsym确定。 compose(f,g)是compose(f,g,v,w,t)的缺省形式。此时,v,w都是由find
13、sym自动对函数f和g进行辨认而定的自变量。,例 求f = x2 的反函数 syms x; f=x2;g=finverse(f) fg=simple(compose(g,f)说明: 事实上,x2的反函数有两个,finverse只给出“主”反函数,(1)自变量由机器确定 syms x y u fai t; f=x/(1+u2);g=cos(y+fai); fg1=compose(f,g) (2)指定自变量 fg2=compose(f,g,u,fai,t),5.2.3 置换及其应用 1、 自动执行的子表达式置换指令,RS,ssub=subexpr(S,ssub) 运用符号变量ssub置换子表达式,
14、重写S为RS。 子表达式的置换表示 syms a b c d W; V,D=eig(a b;c d) RVD,W=subexpr(V;D,W),被置换的子表达式是机器自动寻找的。其置换原则与pretty指令相同 RVD的上两行子阵是置换后的特征向量阵;下两行是特征值阵,2、 通用置换指令,RES=subs(ES,old,new) 用new置换ES中的old后产生RES RES=subs(ES,new) 用new置换ES中的自由变量后产生RES。 RES=subs(ES) 用当前内存中已知值置换ES中所有可能的同名变量后产生RES。 subs 的置换规则 (1)产生符号函数 syms a x;f
15、=a*sin(x)+5;,(2)符号变量置换 f1=subs(f,sin(x),sym(y) (3)符号常量置换 f2=subs(f,a,x,2,sym(pi/3) (4)双精度数值置换 f3=subs(f,a,x,2,pi/3) (5)数值数组置换 f4=subs(subs(f,a,2),x,0:pi/6:pi) f5=subs(f,a,x,0:6,0:pi/6:pi),digits 显示当前采用的数值计算的精度 digits(n) 设置今后数值计算以n位相对精度进行 xs=vpa(x) 在digits指定精度下,给出x的数值型符号结果xs。 xs=vpa(x,n) 在n位相对精度下,给出x
16、的数值型符号结果xs,5.2.4 符号数值精度控制和任意精度计算任意精度的符号数值,指令使用实例 p0=sym(1+sqrt(5)/2); p1=sym(1+sqrt(5)/2) e01=vpa(abs(p0-p1) p2=vpa(p0) e02=vpa(abs(p0-p2),64) digits,5.2.5 符号对象与其它数据对象间的转换,符号、数值间的转换 phi=sym(1+sqrt(5)/2) double(phi) 各种多项式表示形式之间的转换 syms x;f=x3+2*x2-3*x+5; sy2p=sym2poly(f) p2st=poly2str(sy2p,x) p2sy=po
17、ly2sym(sy2p) pretty(f,x),s = symsum(f,v,a,b) 求通式f在指定变量v取遍a,b中所有整数时的和。 说明 f是矩阵时,求和对元素逐个进行,但自变量定义在整个矩阵上。 v缺省时,f中的自变量由findsym自动辨认;b可以取有限整数,也可以取无穷大。 a,b可同时缺省,此时默认求和的自变量区间为0,v-1,5.3 符号微积分,5.3.1 符号序列的求和,例 求,syms k t;f1=t k3;f2=1/(2*k-1)2,(-1)k/k; s1=simple(symsum(f1) s2=simple(symsum(f2,1,inf) 通式中 的自变量只取整
18、数值 指令中的f可以是符号矩阵,此时求和操作将对矩阵中的“元素通式”逐个进行。但通式矩阵的自变量及其取值区间对各“元素通式”是相同的。 由于findsym对f符号矩阵的自变量进行确认时,不是对该矩阵元素分别进行的,而是把矩阵f作为一个统一的对象处理的。所以当f矩阵中含有一个以上非常数符号变量时要特别注意,求和对哪个符号变量进行。,5.3.2 符号微分和jacobian 矩阵,dfdvn=diff(f,v,n) 求,fjac=jacobianfv(f,v) 求多元向量函数f的jacobian矩阵,f是矩阵时,求导对元素逐个进行,但自变量定义在整个矩阵上。 v缺省时,自变量会自动由findsym确
19、认;n缺省时,默认n=1。 注意:在数值计算中,指令diff是用来求差分的。,syms a t x;f=a,t3;t*cos(x), log(x); df=diff(f) dfdt2=diff(f,t,2) dfdxdt=diff(diff(f,x),t),syms x1 x2 x3; f=x1*exp(x2);x2;cos(x1)*sin(x2); v=x1 x2; fjac=jacobian(f,v),通用积分指令intf=int(f,v) 给出f对指定变量v的(不带积分常数的)不定积分 intf =int(f,v,a,b) 给出f对指定变量v的定积分 说明: 当f是矩阵时,积分将对元素逐
20、个进行 v缺省时,积分对findsys确认的变量进行 a,b分别是积分的下、上限,准许他们取任何值或符号表达式,5.3.3 符号积分,syms a b x; f=a*x,b*x2;1/x,sin(x); disp(The integral of f is); int(f) pretty(int(f),syms x y z F2=int(int(int(x2+y2+z2,z,sqrt(x*y),x2*y),y,sqrt(x),x2),x,1,2) VF2=vpa(F2),5.4 符号积分变换,Fourier 变换及其反变换,Fw=fourier(ft,t,w) 求“时域”函数ft的Fourier
21、变换Fw ft=ifourier(Fw,w,t)求“频域”函数Fw的Fourier反变换ft,(1)求Fourier 变换 syms t w;ut=sym(Heaviside(t); UT=fourier(ut) UTS=simple(UT)(2)求Fourier 反变换进行验算 Ut=ifourier(UT,w,t) Uts=ifourier(UTS,w,t),5.5 符号代数方程的求解,5.5.1 线性方程组的符号解 例 求d+n/2+p/2=q, n+d+q-p=10, q+d-n/4=p,q+p-n-8d=1线性方程组的解,A=sym(1 1/2 1/2 -1;1 1 -1 1;1 -
22、1/4 -1 1;-8 -1 1 1); b=sym(0;10;0;1); X1=Ab,5.5.2 一般代数方程组的解,S=solve(eq1,eq2,eqn,v1,v2,vn) 求方程组关于指定变量的解(标准格式) S=solve(exp1,exp2,expn,v1,v2,vn) 求方程组关于指定变量的解(可用格式) 说明; eq1,eq2,eqn或是字符串表达的方程,或是字符串表达式; v1,v2,vn)是字符串表达的求解变量名。 exp1,exp2,expn只能是符号表达式;v1,v2,vn 是求解的符号变量,如eq1,eq2,eqn是不含“等号”的表达式,则指令是对eq1=0,eq2=
23、0,.,eqn=0求解 S是一个构架数组。 在得不到“封闭型解析解”时,给出数值解。,例 求方程组uy2 + vz + w = 0 , y +z+w= 0 关于y, z 的解。 S=solve(u*y2+v*z+w=0,y+z+w=0,y,z) disp(S.y),disp(S.y),disp(S.z),disp(S.z),5.6 符号微分方程的求解,求微分方程符号解的一般指令,S=dsolve(a_1,a_2,a_n) 求解符号常微分方程最完整、通用的指令调用格式 说明 输入参数包括3部分内容:微分方程、初始条件、指定独立变量。其中微分方程是必不可少的输入内容。输入参数必须以字符形式编写。
24、在本调用格式中,输出参数S是“结构对象”。,图示微分方程y = xy ( y)2 的通解和奇解的关系 y=dsolve(y=x*Dy-(Dy)2,x) clf,hold on,ezplot(y(2),-6,6,-4,8,1) cc=get(gca,Children); set(cc,Color,r,LineWidth,5) for k=-2:0.5:2;ezplot(subs(y(1),C1,k),-6,6,-4,8,1); end hold off, title(fontname隶书fontsize16通解和奇解),mfun 对maple中若干经典特殊函数实施数值计算 mfunlist 能被
25、mfun计算的maple经典特殊函数列表 mhelp 查阅maple中的库函数及其调用方法 maple 进入maple工作空间计算,结果送回matlab工作间 procread 把按maple格式写的源程序读入maple工作空间,5.7 可视化数学分析界面,funtool,5.7.1 单变量函数分析的交互界面,5.7.2 泰勒级数逼近分析界面taylortool,练习题,syms x ; f=sin(x);g=x2-1; f1 = finverse(f) g1 = finverse(g) g2=factor(g) fg3=compose(f,g) f4=subs(f,x,pi/6) g4=subs(g,x,2)syms x y z F2=int(int(int(2*x2+3*y2+z2,z,sqrt(y),y),y,sqrt(x),x),x,2,4) VF2=vpa(F2),