MATLAB符号计算.docx
《MATLAB符号计算.docx》由会员分享，可在线阅读，更多相关《MATLAB符号计算.docx（32页珍藏版）》请在冰豆网上搜索。
MATLAB符号计算
第5章MATLAB符号计算
Matlab本身并没有符号计算功能，1993年通过购买Maple的使用权后，开始具备符号运算的功能．符号运算的类型很多，几乎涉及数学的所有分支．
符号运算使用一种特殊的数据类型，称为符号对象（SymbolicObject），用字符串形式表达，但又不同于字符串（CharArray）。
符号中的变量，函数和表达式都是符号对象。
Forpersonaluseonlyinstudyandresearch;notforcommercialuse
5.1符号对象的定义
5.1.1建立符号变量
MATLAB提供了两个建立符号对象的函数：
sym和syms，两个函数的用法不同。
Forpersonaluseonlyinstudyandresearch;notforcommercialuse
1.sym函数
sym函数用来建立单个符号对象，一般调用格式为：
符号对象名=sym（'符号字符串'）
该函数可以建立一个符号对象，符号字符串可以是常量、变量、函数或表达式。
符号对象名=sym（num,’d’）
该函数将数值表达式num转化为符号表达式，数值用十进制表示；
试比较分析下列命令:
a=pi^2,b=’pi^2’,c=sym（pi^2）,d=sym（pi^2,’d’）,e=sym（b）
>>a=pi^2%结果为一个数
a=
9.8696
>>b='pi^2'%结果为一个字符串
b=
pi^2
>>c=sym（pi^2）%数值转化为符号对象，有理表示！
c=
*2^（-49）
>>c=sym（pi^2,'d'）%数值转化为符号对象，十进制表示！
c=
9.
>>e=sym（b）%或e=sym（'pi^2'）字符串转化为符号对象
e=
pi^2
2.syms函数
函数sym一次只能定义一个符号变量，使用不方便。
MATLAB提供了另一个函数syms，一次可以定义多个符号变量（对象）。
syms函数的一般调用格式为：
syms符号变量名1符号变量名2…符号变量名n
用这种格式定义符号变量时不要在变量名上加字符串分界符（‘），变量间用空格而不要用逗号分隔。
>>symsxyz
>>h=x^3+2*y^2+e
h=
x^3+2*y^2+pi^2
5.1.2建立符号表达式
含有符号对象的表达式称为符号表达式。
建立符号表达式有以下3种方法：
（1）利用单引号来生成符号表达式。
（2）用sym函数建立符号表达式。
（3）使用已经定义的符号变量组成符号表达式。
>>A=[cb;e-hh-x^3]
A=
[*2^（-49）,pi^2]
[-x^3-2*y^2,2*y^2+pi^2]
5.1.3计算精度和数据类型转换
利用函数sym可以将数值表达式变换成它的符号表达式,函数numeric或eval可以将符号表达式变换成数值表达式。
double（s）将符号对象转化为双精度数值
char（s）将符号对象转化为字符串
digits（n）将数值计算精度设为n位；
x=vpa（s）求s的数值结果
x=vpa（s,n）采用n位计算精度求s的数值结果
>>2^10000
ans=
Inf
>>a=sym
（2）;b=a^10000
b=
216268…709376%很长的整数，准确的，而不是近似的
>>vpa（b）
ans=
.e3011
>>formatlong;pi^2,formatshort%用长格式显示pi的值，再恢复为短格式
ans=
9.108936
试比较以下C=PI^2的显示结果。
>>c=sym（'pi^2'）;
>>vpa（c,16）%16位
ans=
9.1089357
>>vpa（c）%32位，默认
ans=
9.9998761
>>vpa（c,100）%100位，高精度
ans=
9.02
>>vpa（c,2）%2位，低精度
ans=
9.6
>>double（c）%双精度，数值
ans=
9.8696
>>symsx;f=sym（x^3-1）
f=
x^3-1
>>fun=inline（f）
fun=
Inlinefunction:
fun（x）=x.^3-1
5.2符号矩阵和符号函数
5.2.1符号矩阵
MATLAB大部分矩阵和数组运算符及指令都可以用于符号如：
inv，./，.\，eig等
>>clear;A=sym（'[a,b;c,d]'）;
>>B=inv（A）
B=
[d/（a*d-b*c）,-b/（a*d-b*c）]
[-c/（a*d-b*c）,a/（a*d-b*c）]
>>A.\B,A\B
ans=
[d/（a*d-b*c）/a,-1/（a*d-b*c）]
[-1/（a*d-b*c）,a/（a*d-b*c）/d]
ans=
[（d^2+b*c）/（d^2*a^2-2*a*d*b*c+b^2*c^2）,-b*（a+d）/（d^2*a^2-2*a*d*b*c+b^2*c^2）]
[-c*（a+d）/（d^2*a^2-2*a*d*b*c+b^2*c^2）,（b*c+a^2）/（d^2*a^2-2*a*d*b*c+b^2*c^2）]
>>A*B
ans=
[a*d/（a*d-b*c）-b*c/（a*d-b*c）,0]
[0,a*d/（a*d-b*c）-b*c/（a*d-b*c）]
>>eig（A）
ans=
[1/2*d+1/2*a+1/2*（d^2-2*a*d+a^2+4*b*c）^（1/2）]
[1/2*d+1/2*a-1/2*（d^2-2*a*d+a^2+4*b*c）^（1/2）]
2.符号函数计算
大部分MATLAB的数学函数和逻辑关系运算也可以用于符号对象，另外还有：
factor（expr）对expr作因式分解
expand（expr）将expr展开
collect（expr,v）将expr按变量v合并同类项
simplify（s）：
应用函数规则对s进行化简。
simple（s）：
调用MATLAB的其他函数对表达式进行综合化简，并显示化简过程。
g=finverse（f,v）求函数f（v）的反函数g（v）
fg=compose（f,g）求函数f（v）和g（v）的复合函数f（g（v））
[n,d]=numden（expr）分式通分，n返回分子，d返回分母
如果符号表达式是一个有理分式或可以展开为有理分式，可利用numden函数来提取符号表达式中的分子或分母。
例1.已知f（x,y）=（x-y）3,g（x,y）=（x+y）3，s=x^6+1;考虑相关的运算，如求h=f*g并展开。
>>clear;symsxyz;f=（x-y）^3;g=（x+y）^3;
>>h=f*g%两函数相乘
h=
（x-y）^3*（x+y）^3
>>hs=expand（h）%函数展开
hs=
x^6-3*x^4*y^2+3*x^2*y^4-y^6
>>s=x^6+1;sf=factor（s）%因式分解
sf=
（x^2+1）*（x^4-x^2+1）
>>sz=subs（sf,x,z^2+z+1）%用z^2+z+1替换sf中的x
sz=
（（z^2+z+1）^2+1）*（（z^2+z+1）^4-（z^2+z+1）^2+1）
>>[n,d]=numden（x/y+y/x）%分式通分
n=
x^2+y^2
d=
x*y
下面是复合函数与反函数的例子：
>>clear;symsxyzt;
>>f=1/（1+x^2）;g=sin（y）;
>>A1=compose（f,g）%求复合函数
A1=
1/（sin（y）^2+1）
>>A2=compose（f,g,t）%求复合函数，变量用t表示
A2=
1/（sin（t）^2+1）
>>v=finverse（f）
Warning:
finverse（1/（x^2+1））isnotunique.
>InD:
\toolbox\symbolic\@sym\finverse.matline43
v=
1/x*（-x*（x-1））^（1/2）%结果不唯一。
5.3符号表达式中变量的确定
MATLAB中的符号可以表示符号变量和符号常量。
findsym可以帮助用户查找一个符号表达式中的的符号变量。
该函数的调用格式为：
findsym（s,n）
函数返回符号表达式s中的n个符号变量，若没有指定n，则返回s中的全部符号变量。
符号变量确定原则
（1）除了i和j之外，字母位置最接近x的字母；若距离相等，则取ASCII码大的；
（2）若没有除了i与j以外的字母，则视x为默认的符号变量；（微分方程视t为默认的符号变量）
>>symsabxy
>>findsym（3*a*b+y^2+1,1）
ans=
y
>>findsym（3*a*b+y^2+1）
ans=
a,b,y
5.4符号微积分
5.4.1符号极限
limit函数的调用格式为：
（1）limit（f,x,a）：
求符号函数f（x）的极限值。
即计算当变量x趋近于常数a时，f（x）函数的极限值。
（2）limit（f,a）：
求符号函数f（x）的极限值。
由于没有指定符号函数f（x）的自变量，则使用该格式时，符号函数f（x）的变量为函数findsym（f）确定的默认自变量，即变量x趋近于a。
（3）limit（f）：
求符号函数f（x）的极限值。
符号函数f（x）的变量为函数findsym（f）确定的默认变量；没有指定变量的目标值时，系统默认变量趋近于0，即a=0的情况。
（4）limit（f,x,a,'right'）：
求符号函数f的极限值。
'right'表示变量x从右边趋近于a。
（5）limit（f,x,a,‘left’）：
求符号函数f的极限值。
‘left’表示变量x从左边趋近于a。
例5-1求下列极限（改用课本例子）。
（1）
（2）
（3）
（4）
>>symshnx
>>L=limit（'（log（x+h）-log（x））/h',h,0）%单引号可省略掉
L=1/x
>>M=limit（'（1-x/n）^n',n,inf）
M=exp（-x）
>>limit（（1-cos（x））/（x*sin（x））,x,0）
ans=
1/2
>>symsx;
>>f=（sqrt（x）-sqrt
（2）+sqrt（x-2））/sqrt（x*x-4）;
>>limit（f,x,2,'right'）
ans=
1/2
例5-2（A）（B）（C）二元，其中C不存在！
1.（A）求二元函数的极限：
>>fun1='（x^3*y+x*y^3+x*y）/（x+y）^3';
>>limit（limit（fun1,x,2）,y,1）
ans=
4/9
（B）1.求极限:
>>symsmn
>>limit（symsum（m/n^2,m,1,n）,n,inf）
ans=
1/2
（C）
>>symsxy
>>limit（limit（x*y/（x^2+y^2）,x,0）,y,0）
ans=
0
注：
此二元函数的极限不存在，但计算结果为零，是不正确的。
要注意二元函数的极限与二次极限的区别。
故使用MATLAB计算极限，要在极限存在的前提下进行。
5.4.2符号导数
diff函数用于对符号表达式求导数。
该函数的一般调用格式为：
diff（s）：
没有指定变量和导数阶数，则系统按findsym函数指示的默认变量对符号表达式s求一阶导数。
diff（s,'v'）：
以v为自变量，对符号表达式s求一阶导数。
diff（s,n）：
按findsym函数指示的默认变量对符号表达式s求n阶导数，n为正整数。
diff（s,'v',n）：
以v为自变量，对符号表达式s求n阶导数。
例5-3
（1）设z=e2xcos（3y），求
（化简过程中只写最后一个最简结果）。
2.设z=e2xcos（3y），求
并求
（2）y=ex（sinx+cosx）
>>z=exp（2*x）*cos（3*y）;
>>f1=diff（diff（z,x）,y）
f1=
-6*exp（2*x）*sin（3*y）
%或f1=diff（diff（z,y）,x）结果同上!
>>ff2=subs（f1,x,1）;
>>ff2=subs（ff2,y,pi）
ff2=-1.6288e-014
（2）
>>clear;symsxy
>>y=exp（x）*（sin（x）+cos（x））;
>>dy=diff（y,x）
dy=
exp（x）*（sin（x）+cos（x））+exp（x）*（cos（x）-sin（x））
>>simple（dy）
ans=
2*exp（x）*cos（x）
>>simplify（dy）
ans=
2*exp（x）*cos（x）
例5-3求下列函数的导数。
（同前面合并！
）
（2）
，
>>clear;symsxy
>>s=diff（x^2*exp（-y）,x,2）
s=
2*exp（-y）
>>t=diff（x^2*exp（-y）,x）;
>>t=diff（t,y）;
>>t=subs（t,x,1）
>>t=subs（t,y,2）
t=
-0.2707
例5-4.
（1）求由方程y5+2y-x-3x7=0所确定的隐函数的导数
。
（2）.求函数z=ln（x2+y2）的全微分,并化简（化简过程中只写最后一个最简结果）。
（3）求摆线
，在t=π所对应的点的切线斜率
（4）求y=ln（sinx）的微分
>>symsxy;
>>fxy=y^5+2*y-x-3*x^7;
>>fx=diff（fxy,'x'）;fy=diff（fxy,'y'）;
>>dv=-fx/fy
dv=（1+21*x^6）/（5*y^4+2）
>>symsxydxdy
>>z=log（x^2+y^2）;
>>dz=diff（z,x）*dx+diff（z,y）*dy;
>>dz=simple（dz）
dz=2*（x*dx+y*dy）/（x^2+y^2）
>>symsat;
>>x=a*（t-sin（t））;y=a*（1-cos（t））;
>>f1=diff（y,t）;f2=diff（x,t）;
>>f=f1/f2
f=
sin（t）/（1-cos（t））
>>ff=inline（f）;
>>k=ff（pi）
k=
6.1232e-017
>>symsxydxdy
>>y=log（sin（x））;
>>f1=diff（y,x）;
>>dy=f1*dx
dy=
cos（x）/sin（x）*dx
5.4.3符号积分
符号积分由函数int来实现。
该函数的一般调用格式为：
int（s）：
没有指定积分变量和积分阶数时，系统按findsym函数指示的默认变量对被积函数或符号表达式s求不定积分。
int（s,v）：
以v为自变量，对被积函数或符号表达式s求不定积分。
int（s,v,a,b）：
求定积分运算。
a,b分别表示定积分的下限和上限。
该函数求被积函数在区间[a,b]上的定积分。
a和b可以是两个具体的数，也可以是一个符号表达式，还可以是无穷（inf）。
当函数f关于变量x在闭区间[a,b]上可积时，函数返回一个定积分结果。
当a,b中有一个是inf时，函数返回一个广义积分。
当a,b中有一个符号表达式时，函数返回一个符号函数。
当系统求不出解析解，会自动求原点附近的一个近似解。
例5-5求下列积分。
>>symsx
>>f=（x^2+1）/（x^2-2*x+2）^2;
>>g=cos（x）/（sin（x）+cos（x））;
>>h=exp（-x^2）;
>>I=int（f）
>>J=int（g,0,pi/2）
>>K=int（h,0,inf）
结果为：
I=1/4*（2*x-6）/（x^2-2*x+2）+3/2*atan（x-1）
J=1/4*pi
K=1/2*pi^（1/2）
4.求不定积分：
并化简（化简过程中只写最后一个最简结果）。
（4）
>>f2=int（x*exp（a*x））
f2=
1/a^2*（a*x*exp（a*x）-exp（a*x））
>>simple（f2）
ans=
exp（a*x）*（a*x-1）/a^2
例5-5计算积分（无解析解）
（3）
（4）
（5）重积分
>>symsx;t3=int（3*sin（x^2）/x,1,4）
t3=
3/2*sinint（16）-3/2*sinint
（1）%特殊函数sinint（a）表示积分int（sin（a*x）/x,0,1）
>>t3=vpa（t3,5）%用vpa求数值解
t3=
1.0279
>>t4=int（exp（-x^sin（x））,0,1）%求不出解析解
Warning:
Explicitintegralcouldnotbefound.
>InD:
\toolbox\symbolic\@sym\int.matline58
t4=
int（exp（-x^sin（x））,x=0..1）
>>t4=vpa（t4,5）%也可用vpa求数值解
t4=
.45491
>>symsxy
>>iy=int（2*sqrt（1-x^2）,y,-sqrt（1-x^2）,sqrt（1-x^2））;
>>t5=int（iy,x,-1,1）
t5=
16/3
5.5级数
5.5.1级数符号求和
求无穷级数的和需要符号表达式求和函数symsum，其调用格式为：
symsum（s,v,n,m）
其中s表示一个级数的通项，是一个符号表达式。
v是求和变量，v省略时使用系统的默认变量。
n和m是求和的开始项和末项。
例5-6求级数
的和S,以及前十项的部分和S1．
>>symsnx
>>S=symsum（1/n^2,1,inf）
S=
1/6*pi^2
>>S1=symsum（1/n^2,1,10）
S1=
1968329/1270080
重要说明：
当求函数项级数
的和S2时，可用命令：
>>symsnx
>>S2=symsum（x/n^2,n,1,inf）
S2=
1/6*x*pi^2
>>S2=symsum（x/n^2,n,1,10）
S2=
1968329/1270080*x
两点说明：
（1）注意观察S2与S1的细微区别！
（2）当通项公式的Matlab表达式较长时，表达式要加上单引号．后面的练习中会遇到此问题．
例5-7求下列级数之和。
3.求幂级数
的和函数（注意最后结果要化简!
）
（3）symsnx;%（1分）
s1=symsum（1/（2*n+1）/（2*x+1）^（2*n+1）,n,0,inf）;%（4分）
s1=simple（s1）%（3分）
s1=1/2*log（（x+1）/x）
5.5.2函数的泰勒级数
MATLAB提供了taylor函数将函数展开为幂级数，其调用格式为：
taylor（f,v,n,a）
该函数将函数f按变量v展开为泰勒级数，展开到第n项（即变量v的n-1次幂）为止，n的缺省值为6。
v的缺省值与diff函数相同。
参数a指定将函数f在自变量v=a处展开，a的缺省值是0。
例5-8求函数在指定点的泰勒级数展开式。
4.将函数
在x=0点展开成x的幂级数（取前9项，实际结果为4项）
（4）f4=log（（1+x）/（1-x））;%（2分）
taylor（f4,x,0,9）%（6分）
ans=
2*x+2/3*x^3+2/5*x^5+2/7*x^7
4.将函数
在x=-1点展开成x的幂级数（采用默认阶数）
（4）taylor（1/x^2,-1）%（8分）
ans=3+2*x+3*（x+1）^2+4*（x+1）^3+5*（x+1）^4+6*（x+1）^5
5.6符号代数方程求解
在MATLAB中，求解用符号表达式表示的代数方程可由函数solve实现，其调用格式为：
solve（s）：
求解符号表达式s的代数方程，求解变量为默认变量。
solve（s,v）：
求解符号表达式s的代数方程，求解变量为v。
solve（s1,s2,…,sn,v1,v2,…,vn）：
求解符号表达式s1,s2,…,sn组成的代数方程组，求解变量分别v1,v2,…,vn。
注：
线性方程可用linsolve函数
当系统求不出解析解，会自动求原点附近的一个近似解。
例5-9解方程
（1）二次方程ax2+bx+c=0;
（2）无解析解x2-3x+ex=2;
（3）方程组
，
这里y,z是未知量。
>>symsabcx;solve（a*x^2+b*x+c,x）
ans=
[1/2/a*（-b+（b^2-4*a*c）^（1/2））]
[1/2/a*（-b-（b^2-4*a*c）^（1/2））]
>>solve（x^2-3*x+exp（x）-2,x）%无解析解，得一近似解
ans=
1.8160
>>symsyzuvw;s=solve（u*y^2+v*z+w,y+z+w,y,z）
s=
y:
[2x1sym]
z:
[2x1sym]
>>s.y,s.z%察看结构的具体内容，结果略，也可用下列方法更直接
>>[y,z]=solve（u*y^2+v*z+w,y+z+w,y,z）
y=
[-1/2/u*（-2*u*w-v+（4*u*w*v+v^2-4*u*w）^（1/2））-w]
[-1/2/u*（-2*u*w-v-（4*u*w*v+v^2-4*u*w）^（1/2））-w]
z=
[1/2/u*（-2*u*w-v+（4*u*w*v+v^2-4*u*w）^（1/2））]
[1/2/u*（-2*u*w-v-（4*u*w*v+v^2-4*u*w）^（1/2））]
例5-9解下列方程。
***4.问k为何值时，下面的方程组有非零解？
请写出相应的MATLAB命令
x1-3x3=0
x1+2x2+kx3=0
2x1+kx2-x3=0
（4）
>>symsk%（1分）
>>A=[10-3;12k;2k-1];%（1分）
>>D=det（A）%（2分）
D=10-k^2-3*k
>>solve（D）
ans=
[-5]
[2]
当k=-5或k=2时,det（A）=0,从而有非零解!
。
5.7符号常微分方程求解
在MATLAB中，用大写字母D表示导数。
例如，Dy表示y'，D2y表示y''，Dy
升级会员 