matlab课件课堂讲授_第1页
matlab课件课堂讲授_第2页
matlab课件课堂讲授_第3页
matlab课件课堂讲授_第4页
matlab课件课堂讲授_第5页
已阅读5页,还剩23页未读 继续免费阅读

下载本文档

版权说明:本文档由用户提供并上传,收益归属内容提供方,若内容存在侵权,请进行举报或认领

文档简介

1、PAGE PAGE 27第二章 TOC o 1-3 h z u HYPERLINK l _Toc414526443 第2章 符号计算 PAGEREF _Toc414526443 h 1 HYPERLINK l _Toc414526444 2.1符号对象和符号表达式 PAGEREF _Toc414526444 h 3 HYPERLINK l _Toc414526445 2.1.1符号对象的创建和衍生 PAGEREF _Toc414526445 h 3 HYPERLINK l _Toc414526446 2.1.2符号计算中的算符 PAGEREF _Toc414526446 h 6 HYPERLI

2、NK l _Toc414526447 2.1.3符号计算中的函数指令 PAGEREF _Toc414526447 h 7 HYPERLINK l _Toc414526448 2.1.4符号对象的识别 PAGEREF _Toc414526448 h 8 HYPERLINK l _Toc414526449 2.2符号数字及表达式的操作 PAGEREF _Toc414526449 h 8 HYPERLINK l _Toc414526450 2.2.1数值数字与符号数字之间的转换 PAGEREF _Toc414526450 h 9 HYPERLINK l _Toc414526451 2.2.2符号数字

3、的任意精度计算 PAGEREF _Toc414526451 h 10 HYPERLINK l _Toc414526452 2.2.3符号表达式的基本操作 PAGEREF _Toc414526452 h 11 HYPERLINK l _Toc414526453 2.2.4表达式中的置换操作 PAGEREF _Toc414526453 h 11 HYPERLINK l _Toc414526454 2.3符号微积分 PAGEREF _Toc414526454 h 14 HYPERLINK l _Toc414526455 2.3.1极限和导数的符号计算 PAGEREF _Toc414526455 h

4、14 HYPERLINK l _Toc414526456 2.3.2序列/级数的符号求和 PAGEREF _Toc414526456 h 17 HYPERLINK l _Toc414526457 2.3.3符号积分 PAGEREF _Toc414526457 h 18 HYPERLINK l _Toc414526458 2.4微分方程的符号解法 PAGEREF _Toc414526458 h 20 HYPERLINK l _Toc414526459 2.4.1符号解法和数值解法的互补作用 PAGEREF _Toc414526459 h 20 HYPERLINK l _Toc414526460

5、2.4.2求微分方程符号解的一般指令 PAGEREF _Toc414526460 h 20 HYPERLINK l _Toc414526461 2.4.3微分方程符号解示例 PAGEREF _Toc414526461 h 21 HYPERLINK l _Toc414526462 2.5符号变换和符号卷积 PAGEREF _Toc414526462 h 23 HYPERLINK l _Toc414526463 2.5.1Fourier变换及其反变换 PAGEREF _Toc414526463 h 23 HYPERLINK l _Toc414526464 2.5.2Laplace变换及其反变换 P

6、AGEREF _Toc414526464 h 25 HYPERLINK l _Toc414526465 2.5.3Z变换及其反变换 PAGEREF _Toc414526465 h 26 HYPERLINK l _Toc414526466 2.5.4符号卷积 PAGEREF _Toc414526466 h 27 HYPERLINK l _Toc414526467 2.6符号矩阵分析和代数方程解 PAGEREF _Toc414526467 h 29 HYPERLINK l _Toc414526468 2.6.1符号矩阵分析 PAGEREF _Toc414526468 h 29 HYPERLINK

7、l _Toc414526469 2.6.2线性方程组的符号解 PAGEREF _Toc414526469 h 31 HYPERLINK l _Toc414526470 2.6.3一般代数方程组的解 PAGEREF _Toc414526470 h 33 HYPERLINK l _Toc414526471 2.7符号计算结果的可视化 PAGEREF _Toc414526471 h 35 HYPERLINK l _Toc414526472 2.7.1直接可视化符号表达式 PAGEREF _Toc414526472 h 35 HYPERLINK l _Toc414526473 2.7.2符号计算结果的

8、数值化绘图 PAGEREF _Toc414526473 h 40 HYPERLINK l _Toc414526474 习题2 PAGEREF _Toc414526474 h 41 符号对象和符号表达式从经典教科书可知,数学表达式和方程的基本组成是:数字、参数、变量;运算符号(加减乘除等);数学函数(如三角函数、指数函数等)。MATLAB作为面向对象的科学计算语言,也是依靠基本符号对象(包括数字、参数、变量)、运算符及一些预定义函数来构造和衍生符号表达式、符号方程的。符号对象的创建和衍生生成符号对象的基本规则基本符号对象(数字、参数、变量、表达式、函数)都必须借助专门的符号函数指令sym或sym

9、s、symfun定义。任何包含符号对象的表达式或方程,将继承符号对象的属性。符号数字完全准确的数字是所谓的符号(类)数字。它的定义格式如下sym(Num)创建一个符号数字Numsc=sym(Num)创建一个符号常数sc,该常数值准确等于Num说明具体数字Num必须处于单引号内。【例2.1-1】符号(类)数字与数值(类)数字之间的差异。本例演示:它们在创建方式、显示形式不同、类别判断、具体数值上的不同。a=pi+sqrt(5)%创建数值类常数sa=sym(pi+sqrt(5)%创建符号类常数Ca=class(a)%判断a的数据类别Csa=class(sa)%判别sa的数据类别vpa(sa-a)%

10、在32位精度意义上计算两类数字之间的差 a = 5.3777sa =pi + 5(1/2)Ca =doubleCsa =symans =0.000000000000000013822375841085200048593542564188 ? 思考:Sa1=sym(pi+sqrt(5)Sa2=sym(pi+sqrt(5) Sa1 =189209612611719/35184372088832Sa2 =pi + 5(1/2) vpa(Sa1-Sa2) ans =-0.000000000000000013822375841085200048593542564188 符号变量和自由符号变量几种定义基本

11、符号变量的指令格式:syms Para定义符号参数ParaPara=sym(Para)(含义同上)syms Para Flag定义具有Flag指定属性的符号参数ParaPara=sym(Para, Flag)(含义同上)syms Para1 Para2 paraN定义Para1, Para2, ParaN等多个符号参数syms Para1 Para2 paraN Flag定义Para1, Para2, ParaN等多个具有Flag指定属性的符号参数说明参数名Para, Para1, Para2, ParaN分别代表符号的命名规则没有什么特殊。建议:符号参数名不用“x两侧的字母”开头。Flag表

12、示参数属性。它可具体取以下词条:positive限定取正实数;real 限定为实数;unreal不限定的复数。思考:下列表达式谁是变量?式中的a,b称谓参数,它们事先并不指定具体数值。它们是符号参数,是符号运算的基本组成元素。确定自由符号变量的规则:在指定变量名的符号运算中,解题围绕指定变量名进行。在不指定变量名的符号运算中,将按照与小写x的ASCII码距离,自动识别自由符号变量。此后的解题将围绕那被自动识别的变量进行。可实现独立自变量的自动认定的指令如下symvar(EXPR)确认EXPR中所有自由符号变量symvar(EXPR,N)从EXPR中确认出靠x最近的N个自由符号变量说明EXPR可

13、以是符号矩阵。确认对整个矩阵进行。x是首选符号变量,其后的次序规则是:与x的ASCII码值之差的绝对值小的字母优先;差绝对值相同时,ASCII码值大的字母优先。例如,符号变量字母的优先次序为x, y, w, z, v。【例2.1-2】用符号计算研究方程的解。本例演示:自动识别变量的能力;指定变量不同,解也不同;熟悉符号法解题基本指令;专用指令findsym。(1)不指定变量情况syms u v w z%定义符号参数和变量Eq=u*z2+v*z+w;%构成符号表达式result_1=solve(Eq)%采用自动识别变量解方程u*z2+v*z+w=0result_1 =- u*z2 - v*z f

14、indsym(Eq,1)%自动识别出的符号变量 ans =w symvar(Eq,1) ans =w (2)把z指定为变量的情况result_2=solve(Eq,z)%对于指定变量z解方程 result_2 = -(v + (v2 - 4*u*w)(1/2)/(2*u) -(v - (v2 - 4*u*w)(1/2)/(2*u) 【例2.1-4】symvar确定自由变量是对整个矩阵进行的。本例演示:符号矩阵的一种创建方法;矩阵在MATLAB中是作为一个整体看待的。syms a b t u v x yA=a+b*x,sin(t)+u;x*exp(-t),log(y)+vsymvar(A,1)%

15、确定矩阵A中的自由符号变量 A = a + b*x, u + sin(t) x*exp(-t), v + log(y)ans =x 符号计算中的算符由于MATLAB采用了重载技术,使得用来构成符号计算表达式的算符,无论在形状、名称上,还是在使用方法上,都与数值计算中的算符完全相同。下面就符号计算中的基本算符作简扼归纳。(1)基本运算符“+”、“-”、“*”、“”、“/”、“”“.*”、“./”、“.”、“.”表示“数组对应元素间”的乘、除、求幂。“”、“.”分别实现矩阵的共轭转置、非共轭转置。(2)关系运算符在符号对象的比较中,没有“大于”、“大于等于”、“小于”、“小于等于”的概念,而只有是

16、否“等于”的“= =”,“ =”概念。符号计算中的函数指令MATLAB提供面向对象的软件环境。对于不同的数据对象(如数值类和符号类),它借助重载(Overload)技术,把具有相同函数计算功能的文件采用同函数名加以保存。表2.1-1 MATLAB中可调用的符号计算函数指令类别情 况 描 述与数值计算对应关系基本函数三角、双曲、及反函数;除atan2外同名、同用法指数、对数函数(如exp, expm)同名、同用法复数函数(没有幅角函数angle)同名、同用法矩阵分解函数(如eig, svd)同名、同用法方程求解函数solve不同微积分函数(如diff, int)不完全相同积分变换、反变换函数只有

17、离散Fourier绘图函数(如ezplot, ezsurf)数值绘图指令更丰富特殊函数如误差函数erf、贝塞尔函数besselj、第一类完全椭圆积分EllipticK等;通过mfunlist可以看到所有经典函数名无对应函数库函数MuPAD库函数,借助evalin和feval命令可调用比mfunlist所列更广的MuPAD库函数;需要具备MuPAD语言知识。无对应函数说明虽然数值计算与符号计算中有许多同名函数,但是读者在使用函数时,还是要十分注意函数对数据类型的要求,否则容易出错。举例来说,就数字而言,有双精度数字和符号类数字之分。这两种数字,显示形式有时非常相似。但假如把符号类数字输入到某个只

18、对数值数据适用的函数(如plot)中,就一定产生错误。符号对象的识别数据对象属性的指令:class(var)给出变量var的数据类别(如double, sym等)isa(var,Obj)若变量var是Obj代表的类别,给出1,表示“真”whos给出所有MATLAB内存变量的属性【例2.1-5】数据对象及其识别指令的使用。本例演示:不同数据对象在创建、显示形式、大小性质、内存占用上的差异;演示class, isa, whos的具体使用。(1)生成三种不同类型的矩阵,给出不同的显示形式cleara=1;b=2;c=3;d=4;%产生四个数值变量Mn=a,b;c,d%利用已赋值变量构成数值矩阵Mc=

19、a,b;c,d%字符串中的a,b,c,d与前面输入的数值变量无关Ms=sym(Mc)%Ms是一个符号变量,它与前面各变量无关。 (2)三种矩阵的大小不同SizeMn=size(Mn)SizeMc=size(Mc)SizeMs=size(Ms)(3)用class获得每种矩阵的类别CMn=class(Mn)CMc=class(Mc)CMs=class(Ms) (4)用isa判断每种矩阵的类别(若返回1,表示判断正确)isa(Mn,double)isa(Mc,char)isa(Ms,sym) (5)利用whos观察内存变量的类别和其它属性whos Mn Mc Ms%观察三个变量的类别和属性 符号数字

20、及表达式的操作数值数字与符号数字之间的转换数值数字向符号数字的转换sym(Num, r )数值类数字Num 的广义有理表达sym(Num, d )数值类数字Num 的“十进制浮点”近似表达sym(Num, e )数值类数字Num 的带eps误差的理性近似表达sym(Num, f )数值类数字Num 的“十六进制浮点”近似表达说明第一个输入量只能是“数值类数字”。sym(Num, r ) 的简略形式是sym(Num)。该格式下的输出是广义有理表达的符号数字。所谓广义有理表达是指:通过p/q , p*pi/q , sqrt(p) , p*2q , p*10q 各项的累加表达式表示(这里p, q都是

21、正整数)。注意:sym(Num)与sym(Num) 的差别。Num 是数字字符串。sym(Num)体现数字的理论真值。Num的字面表示数字,实际上以双精度近似方式保存。sym(Num)体现理论数字的双精度近似。【例2.2-1T】数值类数组的四种符号化形态和精度比较.a=0.1333333331/3pi/7sqrt(5)7(-3)1/3+1/74(1/3)cos(0.2)pi+sqrt(5); %生成数值类数值数组asr=sym(a,r)%广义有理化形态ase=sym(a,e)%顾及最接近浮点表示的广义有理化形态asd=sym(a,d)%十进制定点化形态asf=sym(a,f)%十六进制浮点化形

22、态 符号数字向双精度数字转换double(Num_sym)把符号数字Num_sym转换为双精度数字符号数字的任意精度计算digits当前环境下符号数字“十进制浮点”表示的有效位数digits(n)设定符号数字“十进制浮点”表示的有效位数xs=vpa(x)据表达式x得到digits指定精度下的符号数字xsxs=vpa(x,n)据表达式x得到n为有效数字的符号数字xs说明变精度函数vpa(x)的运算精度受它之前运行的digits(n)控制。MATLAB对digits指令的缺省精度设置是32位。vpa(x,n) 只在运行的当时起作用。x可以是符号对象,也可以是数值对象,但指令运作后所得结果xs一定是

23、符号数字。【例2.2-1】digits, vpa指令的使用。本例演示:digits观察符号变精度算法的默认精度;sym(Num)与sym(Num)的区别;vpa给出在指定精度上的运算结果。digits%观察当前“十进制浮点”表示符号数字的有效位数p0=sym(1+sqrt(5)/2) %近似、准确?字面数字的完全准确表达pr=sym(1+sqrt(5)/2) %近似、准确?16位“广义有理表示”形式pd=sym(1+sqrt(5)/2,d) %准确?16位精度的32位十进制符号表达e32r=vpa(abs(p0-pr) %32位变精度计算误差e16=vpa(abs(p0-pd),16) %16

24、位变精度计算误差,结果?为0 e32d=vpa(abs(p0-pd) %32位变精度计算P0与pd间的误差符号表达式的基本操作collect(合并同类项)expand(对指定项展开)factor(进行因式或因子分解)horner(转换成嵌套形式)numden(提取公因式)simplify(恒等式简化)pretty(习惯方式显示)等其中最常用的是simple(EXPR)把EXPR转换成最简短形式simplify(恒等式简化)EXPR可以是符号表达式或矩阵。在这种情况下,这些指令将对该矩阵的元素逐个进行操作。【例2.2-2】简化。重点演示:simple的用法。syms xf=(1/x3+6/x2+

25、12/x+8)(1/3);g1=simple(f) g1 =(2*x + 1)3/x3)(1/3) g2=simple(g1) g2 =(2*x + 1)3/x3)(1/3) 说明为找到最少字母的简化式,可能要多次使用simple指令。表达式中的置换操作子表达式置换操作RS,ssub=subexpr(S,ssub)运用符号变量ssub置换子表达式,并重写S为RS【例2.2-3】对符号矩阵进行特征向量分解。本例演示:子表达式置换的作用;subexpr的使用格式。clear all%清空内存和调用过的函数syms a b c d WV,D=eig(a b;c d)RVD,W=subexpr(V;D

26、,W)%注意第一输入量是矩阵 说明被置换子表达式由机器自动寻找。只有比较长的子表达式才被置换。指令 的第一输出量的理解应参照第一输入量的形式。通用置换指令RES=subs(ES,old,new)用new置换ES中的old后产生RESRES=subs(ES,new)用new置换ES中的自由变量后产生RES说明subs指令的输出结果的属性取决于new的属性。当ES中的符号参数和变量,全部被双精度数字置换时,所得结果RES就是双精度数据。对于ES表达式中的标量符号参数或变量old,被 new置换时,new可以是数组。【例2.2-4】用简单算例演示subs的置换规则。(1)产生符号函数syms a b

27、 xf=a*sin(x)+b f =b + a*sin(x) (2)符号表达式置换f1=subs(f,sin(x),sym(y)%class(f1) f1 =b + a*yans =sym (3)符号常数置换(a,b被双精度数字2,5置换,x被符号数字置换)f2=subs(f,a,b,x,2,5,sym(pi/3)%注意:采用花括号的输入量是“元胞数组”,允许各元素的类型不同。class(f2) f2 =3(1/2) + 5ans =sym (4)双精度数值置换f3=subs(f,a,b,x,2,5,pi/3)%注意:方括号中元素数据类型相同class(f3) f3 =3(1/2) + 5an

28、s =sym (5)数值数组置换之一(取a=2,x=0:pi/6:pi)f4=subs(subs(f,a,b,2,5),x,0:pi/6:2*pi)class(f4)plot(0:pi/6:2*pi,f4) f4 = 5, 6, 3(1/2) + 5, 7, 3(1/2) + 5, 6, 5, 4, 5 - 3(1/2), 3, 5 - 3(1/2), 4, 5ans =sym (6)数值数组置换之二(取a=0:6,x=0:pi/6:2*pi)f5=subs(f,a,b,x,0:12,0:12,0:pi/6:2*pi)%注意:再次提醒花括号;自动执行数组乘class(f5)plot(f5) f

29、5 = 0, 3/2, 3(1/2) + 2, 6, 2*3(1/2) + 4, 15/2, 6, 7/2, 8 - 4*3(1/2), 0, 10 - 5*3(1/2), 11/2, 12ans =sym 表达式不能和符号变量同时包含在old内被置换。例如subs(f, a, sin(x),2,sym(y)的运行结果将不是所希望的2*y+5。符号微积分求极限、求导、求和、求积分极限和导数的符号计算limit(f,v,a)limit(f,v,a)求极限limit(f,v,a,right)求右极限limit(f,v,a,left)求左极限dfdvn=diff(f,v,n)求r=taylor(f,

30、n,v,a)把在处展开为幂级数说明f是矩阵时,求极限和求导操作对元素逐个进行,但自变量定义在整个矩阵上。注意:v省缺时,自变量自动由findsym确认;n省缺时,默认n=1。【例2.3-1】试求。syms x kLim_f=limit(1-1/x)(k*x),x,inf) 【例2.3-2】求求, ,。本例展示:求导运算是对矩阵元素逐个进行的。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) 【例2.3-4】,求,。本例演示:(A)在理论层面上:对问题本身的分析;导数的极限

31、定义;区间端点处的导数。(B)在符号计算层面上:subs的变量置换用法;limit的左极限、右极限用法;diff, legend, char的用法。(C)数值绘图指令如何用符号计算结果绘制曲线。首先对问题进行理性审视。观察可知:除处存疑外,是处处光滑可导的。(1)据导数定义,对于,借助极限计算“右导数”clearsyms x df_p=sin(x);%当时,改写如左dfdx_p=limit(subs(f_p,x,x+d)-f_p)/d,d,0,right)%当时,求右导数 dfdx_p_0=subs(dfdx_p,x,0)%求处的导数 dfdx_p =cos(x)dfdx_p_0 =1 (2)

32、据导数定义,对于,借助极限计算“左导数”f_n=sin(-x);%当时,改写如左dfdx_n=limit(subs(f_n,x,x+d)-f_n)/d,d,0,left)dfdx_n_0=subs(dfdx_n,x,0) %求处的导数 dfdx_n =-cos(x)dfdx_n_0 =-1 (3)对直接求导数f=sin(abs(x);dfdx=diff(f,x)%直接求导 dfdx =cos(abs(x)*sign(x) (4)图形观察xn=-3/2*pi:pi/30:0;xp=0:pi/30:3/2*pi;xnp=xn,xp(2:end);hold onplot(xnp,subs(f,x,x

33、np),k,LineWidth,2.5)plot(xn,subs(dfdx_n,x,xn),-r,LineWidth,2.5)plot(xp,subs(dfdx_p,x,xp),:r,LineWidth,2.5)legend(char(f),char(dfdx_n),char(dfdx_p),Location,NorthEast)grid onxlabel(x)hold off 图 2.3-1 函数及其导函数说明具体的导函数如下 或写为指令直接求导只能得到时的,而无法判断处的导数。指令所给结果中的abs(1,x)表示的一阶导数。对于“非0实数”,它等于sign(x);对于0,不定义。隐函数求导

34、【例2.3-5】设,求。本例演示:(A)如何实现隐函数求导。(B)maple指令的应用。(C)特殊指令isolate的用法。(1)对方程(即隐函数)求导clearsyms xg=sym(cos(x+sin(y(x)=sin(y(x) %y必须写为y(x)dgdx=diff(g,x) %对方程求导 (2)用符合规则的新变量名dydx替代dgdx中的diff(y(x),x)dgdx1=subs(dgdx,diff(y(x),x),dydx) (3)对变量dgdx1代表的符号方程关于dydx的求解,使通过x,y表达出来dydx=solve(dgdx1,dydx) 注意解释“此题结果表达方式与一般不同

35、”【例2.3-6】求在处展开的8阶Maclaurin级数。本例演示:一元函数的Taylor级数展开。留作课后练习。r=taylor(f,n,v,a)r=taylor(f,n,v,a)把在处展开为幂级数syms xr=taylor(x*exp(x),9,x,0)%忽略9阶小量的展开 序列/级数的符号求和s=symsum(f,v,a,b) s=symsum(f,v,a,b) 求通式f在指定变量v取遍a,b内所有整数时的和。符号积分intf=int(f,v)intf=int(f,v)给出f对指定变量v的(不带积分常数的)不定积分intf=int(f,v,a,b)给出f对指定变量v的定积分。说明与di

36、ff指令一样,当f是矩阵时,积分将对元素逐个进行。v省缺时,积分对symvar确认的变量进行。a, b分别是积分的下、上限,允许它们取任何值或符号表达式。与数值积分相比,符号积分指令简单,适应性强,但可能占用机器时间很长。积分限非数值时,符号积分有可能给出相当冗长而生疏的“闭”符号表达式,也有可能给不出“闭”解。但此时,假若用户把积分限用具体数值替代,那末符号积分将能给出具有“任意精度”的定积分值。【例2.3-8】求。本例演示:积分指令对符号函数矩阵的作用。syms a b xf=a*x,b*x2;1/x,sin(x)S=int(f) disp(The integral of f is)pre

37、tty(S) 对pretty结果给与解读【例2.3-9】求积分。本例演示:内积分上下限都是函数的情况。syms x y zF2=int(int(int(x2+y2+z2,z,sqrt(x*y),x2*y),y,sqrt(x),x2),x,1,2) VF2=vpa(F2,16)%积分结果用32位数字表示 【例2.3-10】求阿基米德(Archimedes)螺线在到间的曲线长度函数,并求出时的曲线长度。本例演示:(A)数学问题本身的分析。(B)符号变量的属性限定。(C)subs, diff, int, vpa的综合使用。(D)绘图指令ezplot的参变量绘图法,及改变曲线色彩的方法。据数学分析知,

38、为弧长元素,而曲线长度。(1)求曲线长度函数syms a r theta phi %限定符号参数和变量取正值x=r*cos(theta);x=subs(x,r,a*theta);y=r*sin(theta);y=subs(y,r,a*theta);dLdth=sqrt(diff(x,theta)2+diff(y,theta)2);L=simple(int(dLdth,theta,0,phi) (2)时的曲线长度L_2pi=subs(L,a,phi,sym(1,2*pi)%获得完全准确值L_2pi_vpa=vpa(L_2pi)%计算32精度近似值 (3)螺线和螺线长度函数的绘制(图2.3-2)L

39、1=subs(L,a,sym(1);ezplot(L1*cos(phi),L1*sin(phi),0,2*pi)grid onhold onx1=subs(x,a,sym(1);%使参数a数字化y1=subs(y,a,sym(1);h1=ezplot(x1,y1,0,2*pi);%为改变ezplot绘制曲线色彩必须借助“图柄”set(h1,color,r,LineWidth,3)%通过图柄,改变曲线图形对象属性hold off 说明函数中的是角度为时螺线的长度,在极坐标上用“幅值”表示。微分方程的符号解法符号解法和数值解法的互补作用对于符号计算来说,不论是初值问题,还是边值问题,其求解微分方程

40、的指令形式相同,且相当简单。符号计算可能花费较多的机时,可能得不到简单的解析解,可能得不到封闭形式的解,甚至也可能无法求解。没有万能的微分方程一般解法,那么请读者记住:求解微分方程的符号法和数值法有很好的互补作用。求微分方程符号解的一般指令求解符号微分方程最常用的指令格式为如下两种:S=dsolve(S=dsolve(eq1, eq2, , eqn, cond1, cond2, , condn, v)S=dsolve(eq1, eq2, , eqn, cond1, cond2, , condn, v)说明输入量包括三部分:微分方程、初始条件、指定独立变量。其中微分方程是必不可少的输入内容。其余

41、视需要而定,可有可无。输入量必须以字符串形式编写。若不对独立变量加以专门的定义,则默认小写英文字母 t 为独立变量。(不像其他场合,把x默认为变量)微分方程的记述规定:当 y 是“应变量”时,用“Dny”表示“y的n 阶导数”。在t为默认独立变量时,Dy表示;Dny表示。关于初始条件或边界条件的规定:应写成 y(a)=b , Dy(c)=d 等。a, b, c, d 可以是变量使用符以外的其它字符。当初始条件少于微分方程数时,在所得解中将出现任意常数符C1,C2,。解中任意常数符的数目等于所缺少的初始条件数。在本调用格式中,输出量S是“构架对象”。如果y是应变量,那末关于它的解在S.y中。微分

42、方程符号解示例【例2.4-1】求的解。本例演示:t为默认独立变量时的最简单调用;输出量的格式。S=dsolve(Dx=y,Dy=-x)%从它的输出,只能看到2个“域”S.x和S.ydisp(blanks(12),x,blanks(21),y)disp(S.x,S.y) %注意解释 S.x, S.y是“构架”数据类型 说明这是dsolve一种较简略的调用格式。symvar规则确定输出“应变量”为x和y,默认独立变量是t 。C1,C2是任意常数。在写微分方程时,最好遵循“导数在前函数在后,导数阶数降阶”的次序。否则有可能运行出错。【例2.4-2】求微分方程的解。解微分方程y=dsolve(y=x*

43、Dy-(Dy)2,x) y=dsolve(Dy)2-x*Dy+y=0,x) (2)画“解”曲线clf,hold on%hy1=ezplot(y(1),-6,6,-4,8,1);% set(hy1,Color,r,LineWidth,5)Sv=symvar(y(2);for k=-2:0.5:2% y2=subs(y(2),Sv(1),k);%ezplot(y2,-6,6,-4,8,1)endhold off%box onlegend(奇解,通解,Location,Best)ylabel(y)title(fontsize14微分方程, (y )2 xy + y = 0 ,的解) 【例2.4-3】

44、求解两点边值问题: 。本例演示:dsolve解边值问题;可视化微分方程的解;ezplot和plot的混合使用。(1)求解边值问题y=dsolve(x*D2y-3*Dy=x2,y(1)=0,y(5)=0,x) (2)观察“解”的图形ezplot(y,-1,6)hold on%ezplot(1,5,0,0)plot(1,5,0,0,.r,MarkerSize,20)text(1,1,y(1)=0)text(4,1,y(5)=0)%title(x*D2y-3*Dy=x2, y(1)=0,y(5)=0)hold off 符号变换和符号卷积Fourier变换、Laplace变换、Z变换和卷积在信号处理和

45、系统动态特性研究中起着重要作用。本节将讨论这些变换和卷积符号算法实现。Fourier变换及其反变换时域中的与它在频域中的Fourier变换之间存在如下关系(2.5.1-1)(2.5.1-2)由计算机完成这种变换的途径有两条:一是直接调用指令fourier和ifourier进行;另是根据式(2.5.1-1)和(2.5.1-2)定义,利用积分指令int实现。下面只介绍指令fourier和ifourier的使用及注意事项。 Fw=fourier(ft,t,w)Fw=fourier(ft,t,w)求“时域”函数ft的Fourier变换Fwft=ifourier(Fw,w,t)求“频域”函数Fw的Fou

46、rier反变换ftft是以t为自变量的“时域”函数;Fw是以圆频率w为自变量的“频域”函数。【例2.5-2】根据Fourier变换定义,用积分指令求方波脉冲的Fourier变换。本例演示:变换、反变换;heaviside的调用;绘图指令的配用。求Fourier变换HELP:heaviside(x) has the value 0 for x 0, and NaN for x = 0. heaviside is not a function in the strict sense.syms A t wsyms tao %positive%若不限定,反变换失败yt=heaviside(t+tao/

47、2)-heaviside(t-tao/2);Yw=fourier(A*yt,t,w) (2)用反变换验算Yt=ifourier(Yw,w,t) Yt_e=simple(Yt) (3)曲线绘制(设tao=3,A=1)yt3=subs(yt,tao,3)Yw3=subs(Yw,A,tao,1,3)subplot(2,1,1)Ht=ezplot(yt3,-3,3);set(Ht,Color,r,LineWidth,3)%为醒目subplot(2,1,2),ezplot(Yw3) 【例2.5-3】求的Fourier变换,在此是参数,是时间变量。本例演示:fourier的缺省调用格式的使用要十分谨慎;在

48、被变换函数中包含多个符号变量的情况下,对被变换的自变量给予指明,可保证计算结果的正确。syms t x w;ft=exp(-(t-x)*heaviside(t-x)F1=simple(fourier(ft,t,w) %给出以w为频率变量的正确结果 F2=simple(fourier(ft,t)%误把x当作时间变量,又误把t当作频率变量 F3=simple(fourier(ft)%误把x当作时间变量 Laplace变换及其反变换Laplace变换和反变换的定义为(2.5.2-1)(2.5.2-2)与Fourier(2.5.2-1)(2.5.2-2)Fs=laplace(ft,t,s)求“时域”函

49、数ft的Laplace变换Fsft=ilaplace(Fs,s,t)求“频域”函数Fs的Laplace反变换ft说明ft是以t为自变量的“时域”函数;Fs是以复频率s为自变量的“频域”函数。在此给出的是laplace , ilaplace 指令的完整调用格式。【例2.5-4】求的Laplace变换。syms t s;syms a b %positiveDt=dirac(t-a);%单位脉冲函数Ut=heaviside(t-b);Mt=Dt,Ut;exp(-a*t)*sin(b*t),t2*exp(-t);MS=laplace(Mt,t,s) 说明若不对符号参数 a, b进行限定,不能确保hea

50、viside(t-b)被正确变换。dirac(t-a) 函数的数学含义是。Z变换及其反变换一个离散因果序列的Z变换及其反变换的定义为(2.5.3-1)(2.5.3-1)(2.5.3-2)FZ=ztrans(fn,n,z)FZ=ztrans(fn,n,z)求“时域”序列fn的Z变换FZfn=iztrans(FZ,z,n)求“频域”序列FZ的Z反变换fn这两个指令的调用格式及注意事项,与前两节指令情况相同。【例2.5-6】求序列 的Z变换,并用反变换验算。本例演示:离散单位函数的构成;ztrans, iztrans的用法。(1)单位函数及性质验证syms nDelta=sym(charfcn0(n

51、);%定义单位函数D0=subs(Delta,n,0);%计算D15=subs(Delta,n,15);%计算disp(D0,D15);disp(D0,D15) (2)求序列的Z变换syms zfn=2*Delta+6*(1-(1/2)n)%借助自定义的Delta构造序列FZ=simple(ztrans(fn,n,z);disp(FZ = )pretty(FZ)%以易读的形式显示FZ_n=iztrans(FZ,z,n) 说明指令中,charfcnA(x) 是Maple中的表达式(或集)的指征函数(Charactristic function),它满足 。符号卷积由于卷积在信号处理和系统动态特性

52、中的占有特殊的地位,本书在此专辟一节以算例方式讨论符号卷积的机器实现。【例2.5-7】已知系统冲激响应,求输入下的输出响应。本例演示:卷积的时域积分法;simple的反复简化。由系统分析可知,输出响应等于卷积。据此,可推出该题的理论计算结果是。下面是计算机实现的指令。syms t taosyms T ut=exp(-t);%定义系统输入ht=exp(-t/T)/T;%定义系统冲激响应uh_tao=subs(ut,t,tao)*subs(ht,t,t-tao);%运用变量替换指令形成被积函数yt=simple(simple(int(uh_tao,tao,0,t)%实施卷积 【例2.5-8】采用L

53、aplace变换和反变换求上例的输出响应。本例演示:通过变换法求卷积,即系统冲激响应。(p84)【例2.5-9】求函数和的卷积。本例演示:限定性符号变量;heaviside的使用;collect的作用。maple restartsyms tao t;t=sym(t,positive);%有无限定只影响结果表达形式ut=heaviside(t)-heaviside(t-1);ht=t*exp(-t);yt=int(subs(ut,t,tao)*subs(ht,t,t-tao),tao,0,t)yts=simple(yt)%进行同类项归并 yt1=collect(yt,heaviside(t-1)

54、%进行同类项归并 符号矩阵分析和代数方程解符号矩阵分析det(A)det(A)行列式diag(A)取对角元构成向量,或据向量构成对角阵V, D=eig(A)特征值分解,使expm(A)矩阵指数inv(A)矩阵逆 V, J=jordan(A)Jordan分解,使poly(A)矩阵的特征多项式rank(A)矩阵秩U, S, V=svd(A)矩阵的奇异值分解,使【例2.6-1】求矩阵的行列式、逆和特征根。%syms a11 a12 a21 a22%A=a11,a12;a21,a22A=sym(A%d%d,2,2)DA=det(A)IA=inv(A)EA=eig(A) (2)借助公因子表达的参数矩阵特

55、征值分析EA,D=eig(A)DA=subexpr (EA,D),w)例2.6-2】著名的Givens旋转(变换)对矩阵的旋转作用。本例演示:Givens旋转的几何意义;符号矩阵乘法;符号矩阵的数值化;绘图指令的配合。(1)产生符号矩阵,和旋转变换syms tA=sym(sqrt(3)/2,1/2;1/2,sqrt(3)/2)G=cos(t),-sin(t);sin(t),cos(t);%一般的Givens矩阵GA=G*A (2)图示Givens旋转(取转角t=110度)An=subs(GA,t,110/180*pi);%转角为110度的Givens数值矩阵Op=0;0;%原点坐标向量,为绘图

56、用。Ad=double(A);%符号矩阵数值化v1=Op,Ad(:,1);v2=Op,Ad(:,2);%形成A矩阵列向量v11=Op,A(:,1);v21=Op,A(:,2);u1=Op,An(:,1);u2=Op,An(:,2);%旋转后GA矩阵列向量plot(v1(:,1),v1(:,2),-k,v2(:,1),v2(:,2),b)axis(-1,1,-1,1),axis squarehold onhu=plot(u1(:,1),u1(:,2),-k,u2(:,1),u2(:,2),b);set(hu,LineWidth,4)title(Givens Rotation)legend( v1

57、, v2, u1, u2,Location,South)hold off grid on 思考:变量v1,v2,u1,u2,v11的数据类型?线性方程组的符号解【例 2.6-3】求线性方程组的解。本例演示,符号线性方程组的矩阵除求解法;solve 求代数方程解。(1)解法一:矩阵除该方程组的矩阵形式是 。该式简记为。求符号解的指令如下clear allA=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);ticX1=Abt1=toctic%左除,与inv(A)*b相当,但更有效X2=inv(A)*bt2=toc (2)解法二:代数方程通用求解指令solve法ticS=solve(d+n/2+p/2-q,d+n-p+q-10,d-n/4-p+q,-8*d-n+p+q

温馨提示

  • 1. 本站所有资源如无特殊说明,都需要本地电脑安装OFFICE2007和PDF阅读器。图纸软件为CAD,CAXA,PROE,UG,SolidWorks等.压缩文件请下载最新的WinRAR软件解压。
  • 2. 本站的文档不包含任何第三方提供的附件图纸等,如果需要附件,请联系上传者。文件的所有权益归上传用户所有。
  • 3. 本站RAR压缩包中若带图纸,网页内容里面会有图纸预览,若没有图纸预览就没有图纸。
  • 4. 未经权益所有人同意不得将文件中的内容挪作商业或盈利用途。
  • 5. 人人文库网仅提供信息存储空间,仅对用户上传内容的表现方式做保护处理,对用户上传分享的文档内容本身不做任何修改或编辑,并不能对任何下载内容负责。
  • 6. 下载文件中如有侵权或不适当内容,请与我们联系,我们立即纠正。
  • 7. 本站不保证下载资源的准确性、安全性和完整性, 同时也不承担用户因使用这些下载资源对自己和他人造成任何形式的伤害或损失。

评论

0/150

提交评论