MATLAB与数学建模2014-第6章符号计算.ppt

第6章MATLAB符号计算,主要内容：6.1符号计算基础6.2符号导数及其应用6.3符号积分6.4级数6.5代数方程的符号求解6.6常微分方程的符号求解,6.1符号计算基础,MATLAB中符号计算函数是数值计算函数的重载，符号计算工具箱采用的函数和数值计算的函数有一部分同名，为得到准确的在线帮助，应该用helpsym/函数名例如：helpsym/inv,6.1符号计算基础,6.1.1符号对象1.建立符号变量和符号常数(1)sym函数sym函数用来建立单个符号变量和符号表达式，例如，a=sym(a)建立符号变量a，此后，用户可以在表达式中使用变量a进行各种运算。rho=sym(1+sqrt(5)/2)f=rho2-rho-1rho=(1+sqrt(5)/2f=(1/2+1/2*5(1/2)2-3/2-1/2*5(1/2),例6.1考察符号变量和数值变量的差别。在MATLAB命令窗口，输入命令：a=sym(a);b=sym(b);c=sym(c);d=sym(d);%定义4个符号变量w=10;x=5;y=-8;z=11;%定义4个数值变量A=a,b;c,d%建立符号矩阵AB=w,x;y,z%建立数值矩阵Bdet(A)%计算符号矩阵A的行列式det(B)%计算数值矩阵B的行列式A=a,bc,dB=105-811,ans=a*d-b*cans=150,例6.2比较符号常数与数值在代数运算时的差别。pi1=sym(pi);k1=sym(8);k2=sym(2);k3=sym(3);%定义符号变量pi2=pi;r1=8;r2=2;r3=3;%定义数值变量sin(pi1/3)%计算符号表达式值sin(pi2/3)%计算数值表达式值sqrt(k1)%计算符号表达式值sqrt(r1)%计算数值表达式值sqrt(k3+sqrt(k2)%计算符号表达式值sqrt(r3+sqrt(r2)%计算数值表达式值,ans=1/2*3(1/2)ans=0.8660ans=2*2(1/2),ans=2.8284ans=(3+2(1/2)(1/2)ans=2.1010,6.1符号计算基础,(2)syms函数syms函数的一般调用格式为：symsvar1var2varn函数定义符号变量var1,var2,varn等。用这种格式定义符号变量时不要在变量名上加字符分界符()，变量间用空格而不要用逗号分隔。symsabcd,符号计算的结果是符号或符号表达式，如果其中的符号要用具体数值代替，可以用subs函数，例如将A中的符号a以数值5代替，可以用B=subs(A,a,5),clear;symsabcd;A1=a,b;c,d;A2=det(A1);A3=subs(A2,a,1);A3=subs(A3,b,2);A3=subs(A3,c,3);A3=subs(A3,d,4),例6.3计算行列式，当时计算A的数值结果。,A3=-2,6.1符号计算基础,2.建立符号表达式例6.4用两种方法建立符号表达式。在MATLAB窗口，输入命令：U=sym(3*x2+5*y+2*x*y+6)%定义符号表达式Usymsxy;%建立符号变量x、yV=3*x2+5*y+2*x*y+6%定义符号表达式V2*U-V+6%求符号表达式的值,U=3*x2+5*y+2*x*y+6V=3*x2+5*y+2*x*y+6ans=3*x2+5*y+2*x*y+12,例6.5计算3阶范得蒙矩阵行列式的值。设A是一个由符号变量a,b,c确定的范得蒙矩阵。symsabcU=a,b,cA=1,1,1;U;U.2%建立范得蒙符号矩阵det(A)%计算A的行列式值,U=a,b,cA=1,1,1a,b,ca2,b2,c2ans=b*c2-c*b2-a*c2+a*b2+a2*c-a2*b,6.1符号计算基础,例6.6建立x,y的一般二元函数。在MATLAB命令窗口，输入命令：symsxy;f=sym(f(x,y),f=f(x,y),6.1.2基本的符号运算1.符号表达式运算(1)符号表达式的四则运算例6.7符号表达式的四则运算示例。在MATLAB命令窗口，输入命令：symsxyz;f=2*x+x2*x-5*x+x3%符号表达式的结果为最简形式f=2*x/(5*x)%符号表达式的结果为最简形式f=(x+y)*(x-y)%符号表达式的结果不是x2-y2，而是%(x+y)*(x-y),f=-3*x+2*x3f=2/5f=(x+y)*(x-y),6.1符号计算基础,(2)因式分解与展开factor(S)对S分解因式，S是符号表达式或符号矩阵。expand(S)对S进行展开，S是符号表达式或符号矩阵。collect(S)对S合并同类项，S是符号表达式或符号矩阵。collect(S,v)对S按变量v合并同类项，S是符号表达式或符号矩阵。,例6.8对符号矩阵A的每个元素分解因式。命令如下：symsabxy;A=2*a2*b3*x2-4*a*b4*x3+10*a*b6*x4,3*x*y-5*x2;4,a3-b3;factor(A)%对A的每个元素分解因式,A=2*a2*b3*x2-4*a*b4*x3+10*a*b6*x4,3*x*y-5*x24,a3-b3,ans=2*a*b3*x2*(5*b3*x2-2*b*x+a),-x*(-3*y+5*x)4,(a-b)*(a2+a*b+b2),例6.9计算表达式S的值。symsxy;s=(-7*x2-8*y2)*(-x2+3*y2);expand(s)%对s展开collect(s,x)%对s按变量x合并同类项(无同类项)factor(ans)%对ans分解因式expand(cos(x+y),ans=7*x4-13*x2*y2-24*y4ans=7*x4-13*x2*y2-24*y4ans=(8*y2+7*x2)*(x2-3*y2)ans=cos(x)*cos(y)-sin(x)*sin(y),6.1符号计算基础,(3)表达式化简MATLAB提供的对符号表达式化简的函数有：simplify(S)应用函数规则对S进行化简。simple(S)调用MATLAB的其他函数对表达式进行综合化简，并显示化简过程。例：symsxy;s=(x2+y2)2+(x2-y2)2;simple(s)%MATLAB自动调用多种函数对s进行%化简，并显示每步结果,ans=2*x4+2*y4,6.1符号计算基础,2.符号矩阵运算transpose(S)返回S矩阵的转置矩阵。,6.1.3符号表达式中变量的确定MATLAB中的符号可以表示符号变量和符号常数。findsym可以帮助用户查找一个符号表达式中的符号变量。该函数的调用格式为：findsym(S,n)函数返回符号表达式S中的n个符号变量，若没有指定n，则返回S中的全部符号变量。在求函数的极限、导数和积分时，如果用户没有明确指定自变量，MATLAB将按缺省原则确定主变量并对其进行相应微积分运算。可用findsym(S,1)查找系统的缺省变量，事实上，MATLAB按离字符x最近原则确定缺省变量。,symsaxyztfindsym(sin(pi*t)ans=tfindsym(x+i*y-j*z)ans=x,y,z,findsym(a+y,1)returnsy.,6.2符号导数及其应用,6.2.1函数的极限limit函数的调用格式为：limit(f,x,a)limit函数的另一种功能是求单边极限，其调用格式为：limit(f,x,a,right)或limit(f,x,a,left),例6.10求极限。在MATLAB命令窗口，输入命令：symsamx;f=(x(1/m)-a(1/m)/(x-a);limit(f,x,a)%求极限(1)f=(sin(a+x)-sin(a-x)/x;limit(f)%求极限(2)limit(f,inf)%求f函数在x(包括+和-)%处的极限limit(f,x,inf,left)%求极限(3)f=(sqrt(x)-sqrt(a)-sqrt(x-a)/sqrt(x*x-a*a);limit(f,x,a,right)%求极限(4),6.2符号导数及其应用,6.2.2符号函数求导及其应用MATLAB中的求导的函数为：diff(f,x,n)diff函数求函数f对变量x的n阶导数。参数x的用法同求极限函数limit，可以缺省，缺省值与limit相同，n的缺省值是1。,例6.11求函数的导数。symsabtxyz;f=sqrt(1+exp(x);diff(f)%求(1)。未指定求导变量和阶数，按缺省规则处理f=x*cos(x);diff(f,x,2)%求(2)。求f对x的二阶导数diff(f,x,3)%求(2)。求f对x的三阶导数f1=a*cos(t);f2=b*sin(t);diff(f2)/diff(f1)%求(3)。按参数方程求导公式求y对x的导数(diff(f1)*diff(f2,2)-diff(f1,2)*diff(f2)/(diff(f1)3%求(3)。求y对x%的二阶导数f=x*exp(y)/y2;diff(f,x)%求(4)。z对x的偏导数diff(f,y)%求(4)。z对y的偏导数f=x2+y2+z2-a2;zx=-diff(f,x)/diff(f,z)%求(5)。按隐函数求导公式求z对x的偏导数zy=-diff(f,y)/diff(f,z)%求(5)。按隐函数求导公式求z对y的偏导数,6.2符号导数及其应用,ans=1/2/(1+exp(x)(1/2)*exp(x)ans=-2*sin(x)-x*cos(x)ans=-3*cos(x)+x*sin(x)ans=-b*cos(t)/a/sin(t)ans=-(a*sin(t)2*b+a*cos(t)2*b)/a3/sin(t)3ans=exp(y)/y2ans=x*exp(y)/y2-2*x*exp(y)/y3zx=-x/zzy=-y/z,6.2符号导数及其应用,例6.12在曲线y=x3+3x-2上哪一点的切线与直线y=4x-1平行。命令如下：x=sym(x);y=x3+3*x-2;%定义曲线函数f=diff(y);%对曲线求导数g=f-4;solve(g)%求方程f-4=0的根，即求曲线%何处的导数为4,ans=1/3*3(1/2)-1/3*3(1/2),6.3符号积分,6.3.1不定积分在MATLAB中，求不定积分的函数是int，其调用格式为：int(f,x)int函数求函数f对变量x的不定积分。参数x可以缺省，缺省原则与diff函数相同。,例6.13求不定积分。x=sym(x);f=(3-x2)3;int(f)%求不定积分(1)f=sqrt(x3+x4);int(f)%求不定积分(2)g=simple(ans)%调用simple函数对结果化简,ans=27*x-1/7*x7+9/5*x5-9*x3ans=-1/48*(x3+x4)(1/2)*(-16*(x2+x)(3/2)+12*(x2+x)(1/2)*x+6*(x2+x)(1/2)-3*log(1/2+x+(x2+x)(1/2)/x/(x+1)*x)(1/2)g=1/3*(x+1)*x)(1/2)*x2+1/12*(x+1)*x)(1/2)*x-1/8*(x+1)*x)(1/2)+1/16*log(1/2+x+(x+1)*x)(1/2),6.3.2符号函数的定积分定积分在实际工作中有广泛的应用。在MATLAB中，定积分的计算使用函数：int(f,x,a,b),例6.14求定积分。x=sym(x);t=sym(t);int(abs(1-x),1,2)%求定积分(1)f=1/(1+x2);int(f,-inf,inf)%求定积分(2)int(4*t*x,x,2,sin(t)%求定积分(3)f=x3/(x-1)100;I=int(f,2,3)%求定积分(4)double(I)%将上述符号结果转换为数值,ans=1/2ans=pians=2*t*(sin(t)2-4)I=97893129180187301565519001875382615/1192978373971185320372138406360121344ans=0.0821,6.3符号积分,例6.15求椭球的体积。命令如下：symsabcz;f=pi*a*b*(c2-z2)/c2;V=int(f,z,-c,c)V=4/3*pi*a*b*c,例6.16轴的长度为10米，若该轴的线性密度计算公式是f(x)=6+0.3x千克/米(其中x为距轴的端点距离)，求轴的质量。(1)符号函数积分。在MATLAB命令窗口，输入命令：symx;f=6+0.3*x;m=int(f,0,10)(2)数值积分。先建立一个函数文件fx.m：functionfx=fx(x)fx=6+0.3*x;再在MATLAB命令窗口，输入命令：m=quad(fx,0,10,1e-6),6.3符号积分,6.3符号积分,6.3.3积分变换1.傅立叶(Fourier)变换在MATLAB中，进行傅立叶变换的函数是：fourier(fx,x,t)求函数f(x)的傅立叶像函数F(t)。ifourier(Fw,t,x)求傅立叶像函数F(t)的原函数f(x)。,例：symsxt;y=abs(x)Ft=fourier(y,x,t)%求y的傅立叶变换fx=ifourier(Ft,t,x)%求Ft的傅立叶逆变换,y=abs(x)Ft=-2/t2fx=x*(2*heaviside(x)-1),6.3符号积分,2.拉普拉斯(Laplace)变换在MATLAB中，进行拉普拉斯变换的函数是：laplace(fx,x,t)求函数f(x)的拉普拉斯像函数F(t)。ilaplace(Fw,t,x)求拉普拉斯像函数F(t)的原函数f(x)。,6.3符号积分,例6.19计算y=x2的拉普拉斯变换及其逆变换.命令如下：x=sym(x);y=x2;Ft=laplace(y,x,t)%对函数y进行拉普拉斯变换fx=ilaplace(Ft,t,x)%对函数Ft进行拉普拉斯逆变换,Ft=2/t3fx=x2,6.3符号积分,3.Z变换对数列f(n)进行z变换的MATLAB函数是：ztrans(fn,n,z)求fn的Z变换像函数F(z)iztrans(Fz,z,n)求Fz的z变换原函数f(n)例6.20求数列fn=e-n的Z变换及其逆变换。命令如下：symsnzfn=exp(-n);Fz=ztrans(fn,n,z)%求fn的Z变换f=iztrans(Fz,z,n)%求Fz的逆Z变换,Fz=z/exp(-1)/(z/exp(-1)-1)f=exp(-1)n,6.4级数,6.4.1级数的符号求和函数symsum用于对符号表达式求和，调用格式为：symsum(a,n,n0,nn)对表达式a的符号变量n从n=n0到n=nn求和例6.22求级数之和。命令如下：n=sym(n);s1=symsum(1/n2,n,1,inf)%求s1s2=symsum(-1)(n+1)/n,1,inf)%求s2。未指定求和变量，缺省为ns3=symsum(n*xn,n,1,inf)%求s3。此处的求和变量n%不能省略。s4=symsum(n2,1,100)%求s4。计算有限级数的和,6.4级数,6.4.2函数的泰勒级数MATLAB中提供了将函数展开为幂级数的函数taylor，其调用格式为：taylor(f,x,a,k)将函数f对x在x=a处泰勒展开到k-1项。例6.23求函数在指定点的泰勒展开式。命令如下：x=sym(x);f1=(1+x+x2)/(1-x+x2);f2=sqrt(1-2*x+x3)-(1-3*x+x2)(1/3);taylor(f1,x,5)%展开到x的4次幂时应选择k=5,ans=1+2*x+2*x2-2*x4,6.4级数,例6.24将多项式表示成x+1的幂的多项式。命令如下：x=sym(x);p=1+3*x+5*x2-2*x3;f=taylor(p,x,-1,4),f=-8-13*x+11*(1+x)2-2*(1+x)3,6.5代数方程的符号求解,求解方程组的函数是solve，调用格式为：solve(eqn1,eqn2,eqnN,var1,var2,varN)例6.28解方程。命令如下：x=solve(1/(x+2)+4*x/(x2-4)=1+2/(x-2),x)%解方程(1)f=sym(x-(x3-4*x-7)(1/3)=1);x=solve(f)%解方程(2)x=solve(2*sin(3*x-pi/4)=1)%解方程(3)x=solve(x+x*exp(x)-10,x)%解方程(4)。仅标出方程的左端,例6.29求方程组的解。xy=solve(1/x3+1/y3-28,1/x+1/y-4,x,y)%解方程组(1)xy=solve(x+y-98,x(1/3)+y(1/3)-2,x,y)%解方程组(2)Warning:Explicitsolutioncouldnotbefound.InC:MATLABR11toolboxsymbolicsolve.matline136x=emptysymy=对方程组(2)MATLAB给出了无解的结论，显然错误，请看完全与其同构的方程组(3)。输入命令如下：u,v=solve(u3+v3-98,u+v-2,u,v)%解方程组(3),6.6常微分方程的符号求解,MATLAB的符号运算工具箱中提供了功能强大的求解常微分方程的函数dsolve。该函数的调用格式为：dsolve(eqn1,condition,var)该函数求解微分方程eqn1在初值条件condition下的特解。参数var描述方程中的自变量符号，省略时按缺省原则处理，若没有给出初值条件condition，则求方程的通解。dsolve在求微分方程组时的调用格式为：dsolve(eqn1,eqn2,eqnN,condition1,conditionN,var1,varN)函数求解微分方程组eqn1、eqnN在初值条件conditoion1、conditionN下的解，若不给出初值条件，则求方程组的通解，var1、varN给出