1、第四章 符号计算,4.1 符号对象和符号表达式,符号常数、符号变量、符号函数、符号操作等是用来形成符号表达式，严格按照代数、微积分等课程中的规则、公式进行运算，并尽可能给出解析表达式。,一、符号对象的生成和使用,syms argv1 argv2 argvk syms(argv1,argv2,argvk) f=sym(arg,flagn)： flagn可取 d, r f=sym(arg) ：缺省为r argv=sym(argv,flagv): flagv可取positive, real, unreal,【例】符号常数形成中的差异 a1=1/3,pi/7,sqrt(5),pi+sqrt(5) a2
2、=sym(1/3,pi/7,sqrt(5),pi+sqrt(5) a3=sym(1/3,pi/7,sqrt(5),pi+sqrt(5) a23=a2-a3 a1 =0.3333 0.4488 2.2361 5.3777 a2 = 1/3, pi/7,sqrt(5), 6054707603575008*2(-50) a3 = 1/3, pi/7, sqrt(5), pi+sqrt(5) a23 = 0,0,0,189209612611719/35184372088832-pi-5(1/2),注： a3是绝对准确的， a2是近似的,最接近的有理表示,【例】把字符表达式转换为符号变量 y=sym(2
3、*sin(x)*cos(x) y=simple(y) y = 2*sin(x)*cos(x) y = sin(2*x),【例】用符号计算验证三角等式syms fai1 fai2 y=simple(sin(fai1)*cos(fai2)-cos(fai1)*sin(fai2) y =sin(fai1-fai2),【例】求矩阵 的行列式值、逆和特征根syms a11 a12 a21 a22 A=a11,a12;a21,a22 DA=det(A),IA=inv(A),EA=eig(A) A = a11, a12 a21, a22 DA = a11*a22-a12*a21 IA = a22/(a11*
4、a22-a12*a21), -a12/(a11*a22-a12*a21) -a21/(a11*a22-a12*a21), a11/(a11*a22-a12*a21) EA = 1/2*a11+1/2*a22+1/2*(a112-2*a11*a22+a222+4*a12*a21)(1/2) 1/2*a11+1/2*a22-1/2*(a112-2*a11*a22+a222+4*a12*a21)(1/2),【例】验证积分syms A t tao w yf=int(A*exp(-i*w*t),t,-tao/2,tao/2);Yf=simple(yf) Yf = 2*A*sin(1/2*tao*w)/w
5、,三、识别对象类别的指令常用指令：class, isa, whos,【例】数据对象及其识别指令的使用。 (1)生成三种不同类型的矩阵，给出不同的显示形式 clear,a=1;b=2;c=3;d=4; Mn=a,b;c,d Mc=a,b;c,d Ms=sym(Mc) Mn =1 23 4 Mc =a,b;c,d Ms = a, b c, d,（2）三种矩阵的大小不同 SizeMn=size(Mn),SizeMc=size(Mc), SizeMs=size(Ms) SizeMn =2 2 SizeMc =1 9 SizeMs =2 2,（3）用class获得每种矩阵的类别 CMn=class(Mn
6、),CMc=class(Mc),CMs=class(Ms) CMn =double CMc =char CMs =sym,（4）用 isa判断每种矩阵的类别（若返回1，表示判断正确） isa(Mn,double),isa(Mc,char),isa(Ms,sym) ans =1 ans =1 ans =1,（5）利用 whos观察内存变量的类别和其他属性 whos Mn Mc Ms Name Size Bytes ClassMc 1x9 18 char arrayMn 2x2 32 double arrayMs 2x2 408 sym object Grand total is 21 eleme
7、nts using 458 bytes,【例】对独立自由符号变量的自动辨认。 （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) ans = X, Y, a, b, c, delta, theta, x, y,（3）在EXPR中确定一个自由符号变量 findsym(EXPR,1) ans = x,（4）在EXPR中确定2个和3个自由变量时的执行情况 findsym(EXPR,2),findsym
8、(EXPR,3) ans = x,y ans = x,y,theta,【例】findsym确定自由变量是对整个矩阵进行的。 syms a b t u v x y; A=a+b*x,sin(t)+u;x*exp(-t),log(y)+v findsym(A,1) A = a+b*x, sin(t)+u x*exp(-t), log(y)+v ans =x,4.2 符号对象的操作和转换,simple(EXPR) 把EXPR转换成最简短形式 simplify(EXPR) 对EXPR进行综合化简 pretty(EXPR) 以习惯的书写方式显示EXPR,一、 符号表达式的操作,【例】简化 syms x
9、f=(1/x3+6/x2+12/x+8)(1/3); g1=simple(f) g2=simple(g1),g1 = (2*x+1)/x g2 = 2+1/x,二、符号函数的求反和复合,【例】求 的反函数 syms x;f=x2;g=finverse(f) Warning: finverse(x2) is not unique. g = x(1/2),fg=simple(compose(g,f) %验算g(f(x)是否等于x fg = x,二、符号函数的求反和复合,【例】求 的反函数 syms x;f=x2;g=finverse(f) Warning: finverse(x2) is not
10、unique. g = x(1/2),fg=simple(compose(g,f) %验算g(f(x)是否等于x fg = x,【例】求 的复合函数,（1）syms x y u fai t;f=x/(1+u2);g=cos(y+fai);fg1=compose(f,g) fg1 = cos(y+fai)/(1+u2),（2）fg2=compose(f,g,u,fai,t) fg2 = x/(cos(y+t)2+1),fg=compose(f,g) : f(*), v=g(*), fg=f(g(*) fg=compose(f,g,v,w,t) : f(v), v=g(w), fg=f(g(w)|
11、w=t,三、置换操作,1、 子表达式置换操作RS,ssub=subexpr(S,ssub) 运用符号变量ssub置换子表达式，重写S为RS.,【例】把复杂表达式中所含的多个相同子表达式用一个符号代替，使表达简洁。 clear all syms a b c d W V,D=eig(a b;c d)； RVD,W=subexpr(V;D,W) RVD = -(1/2*d-1/2*a-1/2*W)/c, -(1/2*d-1/2*a+1/2*W)/c 1, 1 1/2*d+1/2*a+1/2*W, 0 0, 1/2*d+1/2*a-1/2*W W =(d2-2*a*d+a2+4*b*c)(1/2),2
12、、 通用置换指令 RES=subs(ES,old,new) ：用 new 置换ES中的old后产生 RES。 RES=subs(ES,new)：用new 置换 ES 中的自由变量后产生 RES。,【例】用简单算例演示subs的置换规则。 syms a x;f=a*sin(x)+5; f1=subs(f,sin(x),sym(y) f1 = a*y+5 f2=subs(f,a,x,2,sym(pi/3) f2 = 3(1/2)+5,f3=subs(f,a,x,2,pi/3) f3 =6.7321 f4=subs(subs(f,a,2),x,0:pi/6:pi) f4 =5.0000 6.0000
13、 6.7321 7.0000 6.7321 6.0000 5.0000 f5=subs(f,a,x,0:6,0:pi/6:pi) f5 =5.0000 5.5000 6.7321 8.00008.4641 7.5000 5.0000,4.3 符号微积分,一、 符号序列的求和 s=symsum(f,v,a,b) : a,b缺省时，自变量区间为0,v-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) s1 = 1/2*t*(t-1), k3*t s2 =
14、1/8*pi2, -log(2),二、 符号微分和 矩阵dfdvn=diff(f,v,n) ：求fjac=jacobian(f,v)：求多元向量函数的jacobian矩阵。,【例】求 的 矩阵。 syms x1 x2;f=x1*exp(x2);x2;cos(x1)*sin(x2); v=x1 x2;fjac=jacobian(f,v) fjac = exp(x2), x1*exp(x2) 0, 1 -sin(x1)*sin(x2), cos(x1)*cos(x2),【例】求 、 和syms a t x;f=a,t3;t*cos(x), log(x); df=diff(f) dfdt2=diff
15、(f,t,2) dfdxdt=diff(diff(f,x),t) df = 0, 0 -t*sin(x), 1/x dfdt2 = 0, 6*t 0, 0 dfdxdt = 0, 0 -sin(x), 0,三、 符号积分 intf=int(f,v) ： 给出f对指定变量v的不定积分。 intf=int(f,v,a,b) ： 给出f对指定变量v的定积分。,【例】求syms a b x;f=a*x,b*x2;1/x,sin(x); disp(The integral of f is); pretty(int(f) The integral of f is 2 31/2 a x 1/3 b x lo
16、g(x) -cos(x) ,【例】求（1）求一般积分结果 F1=int(1/log(t),t,0,x) F1 = -Ei(1,-log(x),（2）利用mfun指令求 x=0.5, 0.6, 0.7, 0.8, 0.9 时的定积分 x=0.5:0.1:0.9 F115=-mfun(Ei,1,-log(x) x =0.5000 0.6000 0.7000 0.8000 0.9000 F115 =-0.3787 -0.5469 -0.7809 -1.1340 -1.7758,【例】求syms x y z F2=int(int(int(x2+y2+z2,z,sqrt(x*y),x2*y), y,sq
17、rt(x),x2),x,1,2) VF2=vpa(F2) %积分结果用32位数字表示 F2 = 1610027357/6563700-6072064/348075*2(1/2)+14912/4641*2(1/4)+ 64/225*2(3/4) VF2 = 224.92153573331143159790710032805,4.4 符号代数方程的求解,一、线性方程组的符号解,【例 】求线性方程组 的解。 A=sym(1 1/2 1/2 -1;1 1 -1 1;1 -1/4 -1 1;-8 -1 1 1); b=sym(0;10;0;1);X1=Ab X1 = 1 8 8 9,二、一般代数方程组的
18、解 S=solve(eq1,eq2,eqn,v1,v2,vn) S=solve(exp1,exp2,expn,v1,v2,vn),【例】求方程组 ， 关于 的解。 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) S = y: 2x1 symz: 2x1 sym S.y -1/2/u*(-2*u*w-v+(4*u*w*v+v2-4*u*w)(1/2)-w -1/2/u*(-2*u*w-v-(4*u*w*v+v2-4*u*w)(1/2)-w S.z 1/2/u*(-2*u*w-v+(4*u*w*v+
19、v2-4*u*w)(1/2) 1/2/u*(-2*u*w-v-(4*u*w*v+v2-4*u*w)(1/2),【例】求 的解。 clear all,syms x;s=solve(x+2)x=2,x) s = .69829942170241042826920133106081,【例】求 ， ， 构成的“欠定”方程组解。 syms d n p q; eq1=d+n/2+p/2-q;eq2=n+d+q-p-10; eq3=q+d-n/4-p; S=solve(eq1,eq2,eq3,d,n,p,q); S.d,S.n,S.p,S.q,，,ans = d ans = 8,ans = 4*d+4 ans = 3*d+6,4.5 符号微分方程的求解,求微分方程符号解的一般指令S=dsolve(a_1,a_2,a_n),【例】求 的解。S=dsolve(Dx=y,Dy=-x); disp(blanks(12),x,blanks(21),y), disp(S.x,S.y) x y cos(t)*C1+sin(t)*C2, -sin(t)*C1+cos(t)*C2,
