第2章_符号计算

1、第第2章章 符号计算符号计算教学目标教学目标教学重点教学重点教学内容教学内容教学目标 一是讲述一是讲述MATLAB符号计算基本知识，包符号计算基本知识，包括符号对象的创建、符号数字、符号表达括符号对象的创建、符号数字、符号表达式的操作；式的操作； 二是介绍符号微积分的计算；二是介绍符号微积分的计算； 三是介绍符号矩阵分析和代数方程（组）三是介绍符号矩阵分析和代数方程（组）的符号解法；的符号解法； 四是介绍符号计算结果的可视化。四是介绍符号计算结果的可视化。教学重点 熟悉符号对象的创建、符号数字、符号熟悉符号对象的创建、符号数字、符号表达式的操作。表达式的操作。 熟悉符号微积分的基本计算函数指令

2、。熟悉符号微积分的基本计算函数指令。 熟悉代数方程（组）的符号解法。熟悉代数方程（组）的符号解法。 熟悉符号计算结果可视化的的基本指令。熟悉符号计算结果可视化的的基本指令。 了解符号计算帮助系统及其帮助指令。了解符号计算帮助系统及其帮助指令。 教学内容2.1符号对象和符号表达式符号对象和符号表达式2.2 符号数字及表达式的操作符号数字及表达式的操作 2.3 符号微积分符号微积分2.4 微分方程的符号解法微分方程的符号解法2.5 符号变换和符号卷积符号变换和符号卷积 2.6 符号矩阵分析和代数方程解符号矩阵分析和代数方程解 2.7 代数状态方程求符号传递函数代数状态方程求符号传递函数 2.8 符

3、号计算结果的可视化符号计算结果的可视化 2.9 符号计算资源深入利用符号计算资源深入利用符号运算与数值运算的区别：符号运算与数值运算的区别：符号运算中，解算数学表达式、方程时，不是在离散符号运算中，解算数学表达式、方程时，不是在离散化的数值点上进行，而是凭借一系列恒等式和数学定化的数值点上进行，而是凭借一系列恒等式和数学定理，通过推理和演绎，获得解析结果。这种计算建立理，通过推理和演绎，获得解析结果。这种计算建立在数值完全准确表达和推演严格解析的基础上，所得在数值完全准确表达和推演严格解析的基础上，所得结果是完全准确的。结果是完全准确的。符号运算符号运算-代数运算，公式推导代数运算，公式推导数

4、值运算数值运算-算术运算算术运算代值代值2.1 符号对象和符号表达式符号对象和符号表达式在在Matlab中中,数值和数值变量用于数值的存储和各种数值计算数值和数值变量用于数值的存储和各种数值计算.而而符号常量符号常量,符号变量符号变量,符号函数符号函数,符号操作等则是用来形成符号操作等则是用来形成符号符号表达式表达式,严格按照代数严格按照代数,微积分等课程中的规则微积分等课程中的规则,公式进行运算公式进行运算,并并尽可能给出尽可能给出解析表达式解析表达式.2.1.1 符号对象的创建和衍生符号对象的创建和衍生 数值计算数值计算-变量先赋值变量先赋值,再使用再使用. 符号计算符号计算-先定义基本的

5、符号对象先定义基本的符号对象(可以是可以是常量常量,变变 量量,表达式表达式),然后用这些基本符号对象去构成新的表达然后用这些基本符号对象去构成新的表达式式,再进行所需的符号运算再进行所需的符号运算2.1.1 符号对象的创建和衍生符号对象的创建和衍生1. 生成符号对象的基本规则生成符号对象的基本规则 任何基本符号对象（数字、参数、变量、表达式）任何基本符号对象（数字、参数、变量、表达式）都必须借助专门的符号函数指令都必须借助专门的符号函数指令sym或或syms定义。定义。 任何包含符号对象的表达式或方程，将继承符号任何包含符号对象的表达式或方程，将继承符号对象的属性。即任何包含符号对象的表达式

6、、方程对象的属性。即任何包含符号对象的表达式、方程也一定是符号对象。也一定是符号对象。2 符号数字的定义符号数字的定义格式：格式：sc=sym(num) % sc为值为为值为num的符号数字的符号数字注意：注意： i) 单引号必须在英文状态下输入，构成字符串单引号必须在英文状态下输入，构成字符串 ii) num为一个具体的数字为一个具体的数字如：如： sc=sym(2/3) sb=sym(pi+sqrt(5)sc=2/3sb =pi + 5(1/2)2 符号数字的定义符号数字的定义【例【例2.1-1】符号（类）数字与数值（类）数字之间的差异。】符号（类）数字与数值（类）数字之间的差异。a=pi

7、+sqrt(5)sa=sym(pi+sqrt(5)Ca=class(a)Csa=class(sa)vpa(sa-a) a = 5.3777sa = pi + 5(1/2)Ca = doubleCsa = symans =0.00000000000000001382237584108520004859354256418 本例表现符号数字总是被准确记录和运算，而数值数字并不本例表现符号数字总是被准确记录和运算，而数值数字并不总能保证被准确存储，运算时会引进截断误差。总能保证被准确存储，运算时会引进截断误差。3. 基本符号变量：基本符号变量：定义格式：定义格式：i) syms para para=s

8、ym(para) syms a; a=sym(a)ii) syms para flag para=sym(para, flag) syms a positive; a=sym(a, positive)flag为参数属性：为参数属性：positive-参数取正实数参数取正实数real-参数为实数参数为实数unreal-参数为限定的复数参数为限定的复数iii) syms a b c syms a b c flag无逗号无逗号符号参数符号参数（表达式中的参数）（表达式中的参数）, “待解符号变量待解符号变量”或或“自由符号变量自由符号变量” （表达式中的自变量（表达式中的自变量x，默认默认为为x）2

9、5sin(3)30uzxzwa4. 自由符号变量自由符号变量symvar(expression) 列出表达式中所有基本符号变量列出表达式中所有基本符号变量symvar(expression,n) 列出表达式中认定列出表达式中认定n个自由符号变量个自由符号变量expression是符号表达式或符号表达式矩阵，是符号表达式或符号表达式矩阵，x是是首选自由符号变量，认定优先次序为首选自由符号变量，认定优先次序为x, y, w, z, v等等u解题结果是解题结果是“用符号参数构成的表达式表述自由符号用符号参数构成的表达式表述自由符号变量变量”。解题时自由符号变量可。解题时自由符号变量可“人为指定人为指

10、定”，也可也可“默认地自动认定默认地自动认定”：与小写字母与小写字母 x 的的ASII码距离最码距离最小的变量。小的变量。syms u v w z a5f=sym(3);Eq=sin(f)*u*z2+v*z+f*w-a5;u【例【例2.1-2】 用符号计算研究方程用符号计算研究方程 的解的解。25sin(3)30uzvzwasymvar(Eq) %按字母表顺序列出基本符号变量按字母表顺序列出基本符号变量, 无无 fans = a5, u, v, w, z symvar(Eq,100) %按离按离x的距离列出所有自由符号变量的距离列出所有自由符号变量ans = w, z, v, u, a5res

11、ult_1=solve(Eq) result_1 =a5/3 - (v*z)/3 - (u*sin(3)*z2)/3 result_2=solve(Eq,z) %指定的自由符号变量指定的自由符号变量z解方程解方程result_2 = -(v - (v2 + 4*a5*u*sin(3) - 12*u*w*sin(3)(1/2)/(2*u*sin(3) -(v + (v2 + 4*a5*u*sin(3) - 12*u*w*sin(3)(1/2)/(2*u*sin(3)syms a b x X Yk=sym(3);z=sym(c*sqrt(d)+y*sin(t);EXPR=a*z*X+(b*x2+k

12、)*Y;u【例【例2.1-3】元符号表达式、衍生符号表达式定义，基】元符号表达式、衍生符号表达式定义，基本符号变量、自由符号变量的机器认定。本符号变量、自由符号变量的机器认定。E3=sym(a*sqrt(theta) ? Error using = sym.sym convertExpression at 2515 E4=sym(a*sqrt(theta1)E5=sym(a*sqrt(theta*t)symvar(EXPR)ans = X, Y, a, b, c, d, t, x, y 无无 k zsymvar(EXPR,10) ans = x, y, t, d, c, b, a, X, Y

13、4. 自由符号变量自由符号变量syms a b t u v x yA=a+b*x, sin(t)+u; x*exp(-t), log(y)+vsymvar(A,1) A = a + b*x, u + sin(t) x/exp(t), v + log(y)ans = x u【例【例2.1-4】symvar确定自由变量是对整个矩阵进行的。确定自由变量是对整个矩阵进行的。 2.1.2 符号计算中的算符符号计算中的算符u由于新版由于新版matlab采用了重载采用了重载(Overload)技术技术, 使得使得用来构成符号计算表达式的算符和基本函数用来构成符号计算表达式的算符和基本函数,无论在形无论在形式

14、式,名称名称,还是使用方法上还是使用方法上,都与数值计算中的算符和基都与数值计算中的算符和基本函数几乎完全相同本函数几乎完全相同,这给编程带来极大的方便这给编程带来极大的方便.(1) 基本运算符基本运算符算符算符 ”+”, ”-”, ”*”, ”, “/”,“” 分别构成矩阵的加分别构成矩阵的加, 减减, 乘乘,左除左除,右除右除,求幂运算求幂运算. 算符算符 ”.*”, “./”, “.”, “.” 分别实现元素对元素的数组乘分别实现元素对元素的数组乘,除除,求幂运算求幂运算.算符算符” ”, “ . ” 分别实现矩阵的共轭转置分别实现矩阵的共轭转置,非共轭转置非共轭转置2.1.2 符号计算

15、中的算符符号计算中的算符(2) 关系运算符关系运算符在符号对象的比较中在符号对象的比较中,没有大于没有大于,大于等于大于等于,小于小于,小于等小于等于的概念于的概念,而而只有是否等于的概念只有是否等于的概念。 ”=“ “=“分别分别用来对算符两边的对象进行相等和不等的比较用来对算符两边的对象进行相等和不等的比较,返回返回为逻辑量。事实为真时，比较结果为为逻辑量。事实为真时，比较结果为1，事实为假时，事实为假时，结果为结果为0.2.1.3 符号计算中的函数指令符号计算中的函数指令符号计算中的函数分成三个层次：符号计算中的函数分成三个层次：1.与数值类函数和指令对应的同名符号类函数和指令。与数值类

16、函数和指令对应的同名符号类函数和指令。2.约约50个经典特殊函数（误差函数个经典特殊函数（误差函数erf、贝塞尔函数、贝塞尔函数besselj、椭圆积分、椭圆积分ElliptiK等），借助等），借助mfun调用，调用，用用mfunlist可列出。可列出。3.数量较大的数量较大的MuPAD库函数，借助库函数，借助evalin和和feval调调用。用。(2) 指数指数,对数函数对数函数 sqrt, exp, expm在两者中用法相同符号计算中只在两者中用法相同符号计算中只有自然对数有自然对数log，而没有数值计算中的，而没有数值计算中的log2, log10(3) 复数函数复数函数 conj, i

17、mag, real, abs在两者中用法相同在两者中用法相同.但在符号但在符号计算中没有求相角的指令计算中没有求相角的指令angle.(1) 三角函数三角函数,双曲线函数以及他们的反函数双曲线函数以及他们的反函数除除atan2只能用于数值计算外只能用于数值计算外,另外的在两种运算中另外的在两种运算中使用方法相同使用方法相同.与数值类函数和指令对应的同名符号类函数和指令与数值类函数和指令对应的同名符号类函数和指令(4) 矩阵代数指令矩阵代数指令在符号计算中在符号计算中, matlab提供的常用矩阵代数指令有：提供的常用矩阵代数指令有：diag, tril, inv, det, rank, eig

18、, svd（ Singular value decomposition奇异值分解）等奇异值分解）等与数值类函数和指令对应的同名符号类函数和指令与数值类函数和指令对应的同名符号类函数和指令(5) 方程求解指令方程求解指令solve，与数值类不同。，与数值类不同。(6) 微积分如微积分如diff，int，与数值类不完全相同。，与数值类不完全相同。(7) 积分变换和反变换函数如积分变换和反变换函数如laplace， ilaplace，数，数值类只有值类只有Fourier变换。变换。(8) 绘图函数如绘图函数如ezplot，ezsurf，数值类绘图指令更丰富。，数值类绘图指令更丰富。u数值计算对象数值

19、计算对象,符号计算对象符号计算对象,字符串是字符串是MATALB中最常用中最常用的数据对象的数据对象.他们遵循各自不同的运算法则他们遵循各自不同的运算法则,但有时在外形但有时在外形上却十分相似上却十分相似.MATLAB提供了一些识别不同数据对象的提供了一些识别不同数据对象的指令指令,常用的有常用的有class, isa, whos 例例2.1.3-1 数据对象及其识别指令的使用数据对象及其识别指令的使用(1) 生成三种不同类型的矩阵，给出不同的显示形式生成三种不同类型的矩阵，给出不同的显示形式clear, a=1;b=2;c=3;d=4 % 产生四个数值变量产生四个数值变量Mn=a,b; c,

20、d % 利用已赋值变量构成数值矩阵利用已赋值变量构成数值矩阵Mc=a,b; c,d % 字符串中的字符串中的a,b,c,d与前面输入的数值变量无关与前面输入的数值变量无关Ms=sym(Mc) % 符号变量，与前面的各变量无关符号变量，与前面的各变量无关2.1.4 符号对象的识别符号对象的识别Mn = 1 2 3 4Mc =a,b;c,dMs = a, b c, dSizeMn=size(Mn), SizeMc=size(Mc), SizeMs=size(Ms)SizeMn = 2 2SizeMc = 1 9SizeMs = 2 2(3) 用用class获得每种矩阵的类别获得每种矩阵的类别CMn

21、=class(Mn), CMc=class(Mc), CMs=class(Ms)CMn = doubleCMc = charCMs = sym(4) 用用isa判断每种矩阵的类别判断每种矩阵的类别isa(Mn,double), isa(Mc,char), isa(Ms,sym)ans = 1ans = 1ans = 1Mn = 1 2 3 4 Mc = a,b;c,dMs = a, b c, d(2) 三个矩阵的大小不同三个矩阵的大小不同(5) 利用利用whos观察内存变量的类别和其他属性观察内存变量的类别和其他属性 Name Size Bytes Class Mc 1x9 18 char a

22、rray Mn 2x2 32 double array Ms 2x2 312 sym objecta=0:1:6; theta=sym(30*pi/180*a)alfa=sym(30*pi/180)*atheta =30*pi/180*aalfa = 0, 1/6*pi, 1/3*pi, 1/2*pi, 2/3*pi, 5/6*pi, pia与数组与数组a无关无关whos P26 Name Size Bytes Class ans 1x1 8 double t 1x201 1608 double y 1x201 1608 double 2.1.5 符号运算机理和变量假设符号运算机理和变量假设1

23、.符号运算的工作机理符号运算的工作机理u Matlab借助借助sym或或syms定义符号变量时，启动定义符号变量时，启动MuPAD引擎并启动一个专供引擎并启动一个专供MuPAD使用的内存工作使用的内存工作空间执行符号运算；空间执行符号运算；u matlab内存空间只保存该符号变量和计算结果。内存空间只保存该符号变量和计算结果。u 该定义变量的限定性假设（该定义变量的限定性假设（assumption）被保存在）被保存在MuPAD的内存空间。的内存空间。u 若变量不带限定性假设，则若变量不带限定性假设，则MuPAD默认为复数。默认为复数。2.1.5 符号运算机理和变量假设符号运算机理和变量假设2.

24、 对符号变量的限定性假设对符号变量的限定性假设i) syms x para=sym(x)ii) syms x flag para=sym(x, flag) syms a positive; a=sym(a, positive)iii) syms a b c syms a b c flagflag为参数属性：为参数属性：positive-参数取正实数参数取正实数real-参数为实数参数为实数unreal-参数为限定的复数参数为限定的复数2.1.5 符号运算机理和变量假设符号运算机理和变量假设3. 清除变量和撤销假设清除变量和撤销假设u符号变量和其假设存放在不同的内存空间，因此删除符号变量和其假设

25、存放在不同的内存空间，因此删除符号变量和撤销关于变量的假设是两件需要分别处理符号变量和撤销关于变量的假设是两件需要分别处理的事情。的事情。 clear all 可同时删除可同时删除clear x 删除删除matlab内存中的内存中的x变量变量syms x clear 撤销撤销MuPAD内存中对变量内存中对变量x的假设的假设evalin(symengine,getprop(x) ) 获取获取x的限定性假设的限定性假设 evalin(symengine,anames (Properties) 列出列出MuPAD内存中带限定性假设的符号变量内存中带限定性假设的符号变量reset(symengine)

26、 重启重启MuPAD引擎，清空引擎，清空MuPAD内存内存clear all 删除删除matlab及及 MuPAD内存中的所有变量内存中的所有变量【例【例2.1-6】syms 对符号变量限定性假设的影响对符号变量限定性假设的影响syms x clearf=x3+4.75*x+2.5;rf=solve(f,x)rf = -1/2 1/4 - (79(1/2)*i)/4 (79(1/2)*i)/4 + 1/4 evalin(symengine,getprop(x) ans =C_ syms x realrfr=solve(f,x) rfr =-1/2 evalin(symengine,getpro

27、p(x) ans =R_ clear xsyms xg=x2+x+5;rg=solve(g,x)Warning: Explicit solution could not be found. In solve at 98rg = empty sym evalin(symengine,anames(Properties) ans =xsyms x clearrg=solve(g,x) rg = - (19(1/2)*i)/2 - 1/2 (19(1/2)*i)/2 - 1/2 reset(symengine)clear all2.1.6 符号帮助体系符号帮助体系 uMatlab符号指令帮助系统与符

28、号指令帮助系统与2008年前不同，因为用年前不同，因为用MuPAD取代取代maple作为符号计算的内核。作为符号计算的内核。help SymNamehelpwin SymNamedoc SymNamedoc mfunlistdoc(symengine)2.1.6 符号帮助体系符号帮助体系 uMatlab符号指令帮助系统与符号指令帮助系统与2008年前不同，因为用年前不同，因为用MuPAD取代取代maple作为符号计算的内核。作为符号计算的内核。help SymNamehelpwin SymNamedocsearch laplacedoc mfunlistdoc(symengine)2.1.6

29、符号帮助体系符号帮助体系 uMatlab符号指令帮助系统与符号指令帮助系统与2008年前不同，因为用年前不同，因为用MuPAD取代取代maple作为符号计算的内核。作为符号计算的内核。help SymNamehelpwin SymNamedocsearch laplacedoc mfunlistdoc(symengine)u开启开启MuPAD帮帮助浏览器。助浏览器。2.2 符号数字及表达式的操作符号数字及表达式的操作2.2.1双精度数字与符号数字之间的转换双精度数字与符号数字之间的转换sym(num, r) % “有理分数有理分数”表达的符号数字表达的符号数字 : p/q, n(p/q)sym

30、(num) % sym(num, r) 的简写形式的简写形式1. 双精度数字向符号数字的转换双精度数字向符号数字的转换 ( 52)22eeN*2Nesym(num, f) %以数值以数值 表示表示“浮点数浮点数”，N, e为整为整数数 sym(1/10,f) is 3602879701896397/36028797018963968 sym(num, e) % “有理分数有理分数”表达表达+机器实际浮点表达误差机器实际浮点表达误差e sym(3*pi/4,e) is 3*pi/4 - 103*eps/249sym(num, d) %十进制小数近似表示，有效数字位数受十进制小数近似表示，有效数字

31、位数受 digits指令控制。默认为指令控制。默认为digits(32)情况情况.2.2.1双精度数字与符号数字之间的转换双精度数字与符号数字之间的转换sym(num) % sym(num, r) 的简写形式的简写形式sym(num) % num是字符串数字是字符串数字 1. 双精度数字向符号数字的转换双精度数字向符号数字的转换 sym(num) 中中num 用作符号计算函数的输入用作符号计算函数的输入时，体现其理论时，体现其理论真值真值，对对sym(num)中的中的num用作符号计算函数的输入用作符号计算函数的输入时，体现其字面数子的时，体现其字面数子的双精度近似双精度近似。2.2.1双精度

32、数字与符号数字之间的转换双精度数字与符号数字之间的转换double(num_sym) % 把符号数字转换为双精度数字把符号数字转换为双精度数字2. 符号数字向双精度数字转换符号数字向双精度数字转换 注意：注意：double(num) 把字符串数字转换为字符的把字符串数字转换为字符的ASCII码码double(0.1)ans = 48 46 49f=sym(10/3)Df=double(f)class(Df) Df = 3.3333ans = doubledigits 显示当前环境下符号数字显示当前环境下符号数字“十进制浮点十进制浮点” 表表示的有效数字位数；示的有效数字位数；digits(n)

33、 设定设定 “十进制浮点十进制浮点”表示的有效数字位数；表示的有效数字位数； xs=vpa(x) 据表达式据表达式x得到得到digits指定精度下的符号数字指定精度下的符号数字xs xs=vpa(x,n) 据表达式据表达式x得到得到n位有效数字的符号数字位有效数字的符号数字xs2.2.2 符号数字的任意精度表达形式符号数字的任意精度表达形式 u数值计算与符号计算的最重要区别：数值计算存在截数值计算与符号计算的最重要区别：数值计算存在截断误差，且在计算中不断传播，产生累计误差；符号断误差，且在计算中不断传播，产生累计误差；符号计算完全准确，没有累计误差，但降低计算速度，增计算完全准确，没有累计误

34、差，但降低计算速度，增加内存。为兼顾精度和速度，采用加内存。为兼顾精度和速度，采用“变精度变精度”算法。算法。reset(symengine) sa=sym(1/3+sqrt(2) sa =2(1/2) + 1/3 digits Digits = 32 format longa=1/3+sqrt(2)sa_Plus_a=vpa(sa+a,20)sa_Minus_a=vpa(sa-a,20)a = 1.747546895706428sa_Plus_a =3.4950937914128567869sa_Minus_a =-0.000000000000000022658064826339973669

35、 例例2.2-1 digits, vpa, symengine指令演示指令演示sa32=vpa(sa)digits(48)sa5=vpa(sa,5)sa48=vpa(sa)sa32 =1.747546895706428382135022057543sa5 =1.7475sa48 =1.74754689570642838213502205754303141190300520871 2.2.3 符号表达式的基本操作符号表达式的基本操作collect合并同类项numden获取最小公分母和相应分子expand展开指定项simplify简化表达式，消除冗余项factor因式分解simple搜索寻找最简表

36、达horner转换为嵌套形式 pretty习惯方式显示syms x;f=(1/x3+6/x2+12/x+8)(1/3);g1=simple(f)g1 =(2*x + 1)3/x3)(1/3)验证验证:f2=g13f2 =(2*x + 1)3/x3expand(f2)ans =12/x + 6/x2 + 1/x3 + 8例例2.2-2 简化简化32611238xxxf 2.2.4 表达式中的置换操作表达式中的置换操作1.公共子表达式法简化表达 RS=subexpr(S); RS=subexpr(S,w); RS,w=subexpr(S,w) 其中，w为MATLAB自动寻找的子表达式syms a

37、b c d W;V,D=eig(a,b;c,d);RVD,W=subexpr(V;D,W)【例【例2.2-3】 对符号矩阵对符号矩阵 进行特征向量分解。进行特征向量分解。abAcdRVD = (a/2 + d/2 - w/2)/c - d/c, (a/2 + d/2 + w/2)/c - d/c 1, 1 a/2 + d/2 - w/2, 0 0, a/2 + d/2 + w/2w =(a2 - 2*a*d + d2 + 4*b*c)(1/2) V =(a/2+d/2-(a2-2*a*d+d2+4*b*c)(1/2)/2)/c-d/c, (a/2+d/2+(a2-2*a*d+ d2+4*b*c

38、)(1/2)/2)/c-d/c 1, 1D =a/2 + d/2 - (a2 - 2*a*d + d2 + 4*b*c)(1/2)/2, 0 0, a/2 + d/2 + (a2 - 2*a*d + d2 + 4*b*c)(1/2)/2 2.2.4 表达式中的置换操作表达式中的置换操作2. 通用置换指令通用置换指令RESsubs(ES,old,new)RESsubs(ES,new)例例2.2-4 subs置换规则示例置换规则示例clearsyms a b x;f=a*sin(x)+b f1=subs(f,sin(x),log(y)f2=subs(f,a,3.11)f3=subs(f,a,b,x

39、,2,5,sym(pi/3)class(f3) f1=a* log(y) +b f2 =b + (311*sin(x)/100f3 =3(1/2) + 5ans =sym (5)例例2.2-4 subs置换规则示例置换规则示例format %设为默认输出格式显示设为默认输出格式显示format compactf4=subs(f,a,b,x,2,5,pi/3)class(f4) f4 = 6.7321ans =double f5=subs(f,x,0:pi/2:pi)class(f5)f5 = b, a + b, bans =sym t=0:pi/10:2*pi;f6=subs(f,a,b,x,

40、2,3,t)plot(t,f6)0123456711.522.533.544.55k=(0.5:0.1:1);f6=subs(subs(f,a,b,k,2),x,t);size(f6) plot(t,f6) ans = 6 210123456711.21.41.61.822.22.42.62.83 f6=2*sin(t)+3 f=a*sin(x)+bdfdvn=diff(f,v,n) 求fjac=jacobian(f,v) 求多元向量函数f(v)的jacobian矩阵r=taylor(f,n,v,a) 把f(v)在v=a处进行泰勒展开2.3 符号微积分符号微积分2.3.1 极限和导数的符号计算

41、极限和导数的符号计算limit(f,v,a) 求极限limit(f,v,a,right) 求右极限limit(f,v,a,left) 求左极限)(limvfav)(limvfav)(limvfavnndvvfd)( )10( )()!knkkfaxak【例2.3-1】试求kxxx11limsyms x kf=(1-1/x)(k*x); Lf=limit(f,x,inf) Lf =1/exp(k) Lf1=vpa(subs(Lf,k,sym(-1),48) Lf1=2.7182818284590452353602874713526624977572470937 syms a t x f=a, t

42、3; t*cos(x), log(x); df=diff(f) dfdt2=diff(f,t,2) dfdxdt=diff( diff(f,x) ,t) xxttaflncos3dxdf22dtfddtdxfd2【例2.3-2】 求df = 0, 0 -t*sin(x), 1/xdfdt2 = 0, 6*t 0, 0dfdxdt = 0, 0 -sin(x), 0 例例2.3-5:设设cos(x+siny)=siny, 求求dy/dx(隐函数求导隐函数求导).%将将dy/dx表达式用表达式用x,y表达表达syms xg=sym(cos(x+sin(y(x)=sin(y(x)dgdx=diff(

43、g,x) g =cos(x + sin(y(x) = sin(y(x)dgdx = -sin(x + sin(y(x)*(cos(y(x)*diff(y(x), x) + 1) = cos(y(x)*diff(y(x), x) dgdx1=subs(dgdx,diff(y(x),x),dydx)dgdx1 =-sin(x + sin(y(x)*(dydx*cos(y(x) + 1) = dydx*cos(y(x) dydx=solve(dgdx1, dydx) dydx =-sin(x + sin(y(x)/(cos(y(x) + cos(y(x)*sin(x + sin(y(x)sym xr

44、=taylor(x*exp(x),9,x,0)pretty(r) 例例2.3-6:求求f(x)=xex在在x=0处展开的处展开的8阶阶Maclaurin级数，即忽略级数，即忽略9阶及以上小量的泰勒级数展开阶及以上小量的泰勒级数展开 。 r =x8/5040 + x7/720 + x6/120 + x5/24 + x4/6 + x3/2 + x2 + x 8 7 6 5 4 3 x x x x x x 2 - + - + - + - + - + - + x + x 5040 720 120 24 6 2 MATLAB求解通式求和求解通式求和 问题的指令为：问题的指令为： s=symsum(f,v

45、,a,b) 求通式求通式f在指定变量在指定变量v取遍取遍 a,b中所有中所有整数整数时的和。时的和。 2.3.2 序列序列/级数的符号求和级数的符号求和( )bv af v说明：说明： （1）f是矩阵时，求和对逐个元素进行，但自变量定义在是矩阵时，求和对逐个元素进行，但自变量定义在整个矩阵上。整个矩阵上。（2）v省缺时，省缺时，f中的自变量由中的自变量由findsym（symvar）自动辨自动辨认；认；b可以取有限整数也可以取无穷大。可以取有限整数也可以取无穷大。（3）a，b同时省缺时，默认的自变量求和区间为同时省缺时，默认的自变量求和区间为0,v-1syms n kf1=1/(k*(k+1)

46、;s1=symsum(f1,k,1,n)s1 =1 - 1/(n + 1) 11(1)nkk k12)1(,)12(1kkkk 【例2.3-8】求 ， 。f2=1/(2*k-1)2, (-1)k/k s2=symsum(f2,1,inf) s3 = pi2/8, -log(2) 2.3.3 符号积分【例2.3-9】求 。lnxxdx intf=int(f,v) 给出给出f 对指定变量对指定变量v的不定积分的不定积分 intf=int(f,v,a,b) 给出给出f对指定变量对指定变量v的定积分的定积分syms xf= x*log(x) s=int(f,x)s=simple(s)f = x*log

47、(x)s = (x2*(log(x) - 1/2)/2s = x2*(log(x)/2 - 1/4)2.3.3 符号积分 intf=int(f,v) 给出给出f 对指定变量对指定变量v的不定积分的不定积分 intf=int(f,v,a,b) 给出给出f对指定变量对指定变量v的定积分的定积分syms a b x;f2=a*x, b*x2; 1/x, sin(x)disp(The integral of f is)pretty(int(f2) 【例2.3-10】求 。21sina xb xd xxxThe integral of f is +- -+ | 2 3 | | a x b x | | -

48、 , - | | 2 3 | | | | log(x), -cos(x) | +- -+ f2 = a*x, b*x2 1/x, sin(x) 例例2.3-11求积分求积分syms x y z;f=x2+y2+z2;F2=int(int(int(f, z, sqrt(x*y), x2*y), y, sqrt(x), x2), x, 1, 2)第一重积分第一重积分第二重积分第二重积分2222221()xx yxxyxyz dzdydx Warning: Explicit integral could not be found. F2 =(14912*2(1/4)/4641 - (6072064*

49、2(1/2)/348075 + (64*2(3/4)/225 + 1610027357/6563700VF2=vpa(F2) %默认默认32位有效数字位有效数字 VF2 =224.92153573331143159790710032805教材中教材中48位有效数字位有效数字 22cos,sin,( )( )xryrdlxydsyms a r theta phi r=a*theta; x=r*cos(theta); y=r*sin(theta); dLdth=sqrt(diff(x,theta)2+diff(y,theta)2);L=simple(int(dLdth,theta,0,phi)例例

50、2.3-12 2.3-12 求阿基米德螺线求阿基米德螺线r=ar=a* *(a0)(a0)在在=0=0到到 间的曲线长间的曲线长度函数，并求出度函数，并求出a=1, a=1, 时的曲线长度。时的曲线长度。2220( )( )Lxyd解解:L_2pi=subs(L,a,phi,sym(1,2*pi) L_2pi_vpa=vpa(L_2pi) L =(phi*(a2*phi2 + a2)(1/2) + log(phi + (phi2 + 1)(1/2)*(a2)(1/2)/2 L_2pi_vpa =21.256294148209098800702512272566例例2.3-12 2.3-12 画

51、阿基米德螺线画阿基米德螺线r = a*及其及其曲线长度图形。曲线长度图形。cos,sinxryrL1=subs(L,a,sym(1);ezplot(L1*cos(phi),L1*sin(phi),0,2*pi)grid on; hold onx1=subs(x,a,sym(1);y1=subs(y,a,sym(1);h1=ezplot(x1,y1,0,2*pi);set(h1,Color,r,LineWidth,5) title( ); legend(螺线长度螺线长度-幅角曲线幅角曲线,阿基米德螺线阿基米德螺线)hold offL =(phi*(a2*phi2 + a2)(1/2) + log

52、(phi + (phi2 + 1)(1/2)*(a2)(1/2)/2 -505101520-16-14-12-10-8-6-4-2024xyx = cos() (log( + (2 + 1)1/2)/2 + ( (2 + 1)1/2)/2), y = sin() (log( + (2 + 1)1/2)/2 + ( (2 + 1)1/2)/2)见见P1982.4 微分方程的符号解法微分方程的符号解法S=dslove(eq1,eq2,eqn, cond1,cond2,condn,v)2.4.1 符号解法和数值解法有互补作用符号解法和数值解法有互补作用2.4.2 求微分方程符号解的一般指令求微分方程

53、符号解的一般指令从数值计算角度看，微分方程边值问题的求解比初值问题复杂困从数值计算角度看，微分方程边值问题的求解比初值问题复杂困难；而对符号计算来说，边值问题和初值问题的求解微分方程的难；而对符号计算来说，边值问题和初值问题的求解微分方程的指令形式相同且简单，不过计算时间较长且未必有封闭形式解。指令形式相同且简单，不过计算时间较长且未必有封闭形式解。S=dslove(eq1,eq2,eqn,cond1,cond2,condn,v)默认独立变量默认独立变量t输入量为输入量为字符串字符串必选必选可选可选条件：条件：y(a)=b, Dy(c)=d若应变量为若应变量为y ，用，用“Dny”表示表示“y

54、的的n阶导数阶导数”， Dy为一阶导数。解在为一阶导数。解在S.y中。中。说明：说明：（1 1）输入量包括三部分：微分方程、初始条件、指）输入量包括三部分：微分方程、初始条件、指定独立变量（不指定时，定独立变量（不指定时，默认为默认为t t）。输入量）。输入量必须必须以字符串的形式给出以字符串的形式给出。（2 2）微分方程的记述规定：当）微分方程的记述规定：当y y是应变量时，是应变量时， 用用“DnyDny”表示表示“y y的的n n阶导数阶导数”。（3 3）关于初始条件或边界条件的规定：应写成）关于初始条件或边界条件的规定：应写成y(a)=b, Dy(cy(a)=b, Dy(c)=d)=d

55、等。等。（4 4）如果）如果y y是应变量，则它的解在是应变量，则它的解在S.yS.y中。中。【例2.4-1】求 的解。dxdtydydtx ,2.4.3 微分方程符号解示例微分方程符号解示例clearS=dsolve(Dx=y,Dy=-x) disp(微分方程的解, blanks(2),x,blanks(22),y)disp(S.x,S.y) S = x: 1x1 sym y: 1x1 sym微分方程的解 x y C2*cos(t) + C1*sin(t), C1*cos(t) - C2*sin(t)写微分方程遵循原则：写微分方程遵循原则：导数在前函数在后，导数阶数降阶导数在前函数在后，导数

56、阶数降阶x, y= dsolve(Dx=y,Dy=-x) %按字母表输出顺序不变按字母表输出顺序不变例例2.4-2.2.4-2.图示微分方程图示微分方程y=xyy=xy-(y-(y) )2 2通解和奇解的关系通解和奇解的关系clear ally=dsolve(Dy)2-x*Dy+y=0,x)y = x2/4 C3*x - C32 clf, hold onhy1=ezplot(y(1),-6,6,-4,8,1);set(hy1,Color,r,LineWidth,5)for k=-2:0.5:2 y2=subs(y(2),C3,k); ezplot(y2,-6,6,-4,8,1)endhold

57、off;box onlegend(奇解奇解,通解通解,Location,Best)ylabel(y)title(fontsize14微分方程微分方程, (y )2 xy + y = 0 ,的解的解) 见见P199例例2.4-3.求解两点边值问题求解两点边值问题: xy-3y=x2,y(1)=0,y(5)=0y=dsolve(x*D2y-3*Dy=x2,y(1)=0,y(5)=0,x) y =(31*x4)/468 - x3/3 + 125/468 ezplot(y,-1,6)hold onplot(1,5,0,0,. .r,MarkerSize,20)text(1,1,y(1)=0)text(

58、4,1,y(5)=0)title(x*D2y-3*Dy=x2, y(1)=0,y(5)=0)hold off2.5 符号变换和符号卷积符号变换和符号卷积2.5.1 Fourier变换及其反变换变换及其反变换 Fw=fourier(ft,t,w) 求求“时域时域”函数函数ft的的Fourier变变换换 ft=ifourier(Fw,w,t) 求求“频域频域”函数函数Fw的的Fourier 反变换反变换2.5.2 Laplace变换及其反变换变换及其反变换 Fs=laplace(ft,t,s) 求求“时域时域”函数函数ft的的Laplace变变换换 ft=ilaplace(Fs,s,t) 求求“频

59、域频域”函数函数Fs的的Laplace 反变换反变换2.5.3 Z变换及其反变换变换及其反变换 FZ=ztrans(ft,n,z) 求求“时域时域”序列序列ft的的Z变换变换 fn=itrans(FZ,z,n) 求求“频域频域”序列序列FZ的的Z反变换反变换2.5.4 符号卷积符号卷积 利用拉氏变换求取。利用拉氏变换求取。【例 2.5-1】求 的Fourier变换。0001)(tttf（1）求Fourier变换syms t wut=heaviside(t);UT=fourier(ut) UT = pi*dirac(w)-i/w （2）求Fourier反变换Ut=ifourier(UT,w,t)

60、 Ut =heaviside(t)plot(t(kk),ut(kk),.r,MarkerSize,30)ut(kk)=NaN; hold onplot(t,ut,-r,LineWidth,3)plot(t(kk),t(kk),ut(kk-1), ut(kk+1),or,MarkerSize,10)hold offgrid onaxis(-2,2,-0.2,1.2)xlabel(fontsize14t),ylabel(fontsize14ut)title(fontsize14Heaviside(t) （3）单位阶跃曲线t=-2:0.01:2;ut=heaviside(t);kk=find(t=0