ch2 - 符号计算

第 2 章 符号计算

符号计算：

解算数学表达式、方程不是在离散化的数值点上进行，而是凭借一系列恒等式，数学定理，通过推理和演绎，获得解析结果。

符号计算建立在数值完全准确表达和推演严格解析的基础之上，所得结果完全准确。

特点：

一．相对于MATLAB的数值计算“引擎”和“函数库”而言，符号计算的“引

擎”和“函数库”是独立的。 二．在相当一些场合，符号计算解算问题的指令和过程，显得比数值计算更自然、

更简明。

三．大多数理工科的本科学生在学过高等数学和其他专业基础课以后，比较习惯

符号计算的解题理念和模式。

2.1 符号对象和符号表达式

MATLAB依靠基本符号对象（包括数字、参数、变量）、运算符及

一些预定义函数来构造和衍生符号表达式和符号方程。

2.1.1 一

符号对象的创建和衍生 生成符号对象的基本规则

? 任何基本符号对象都必须借助专门的符号函数指令sym或syms定义。 ? 任何包含符号对象的表达式或方程，将继承符号对象的属性。

1

二 符号数字

符号（类）数字的定义：

sym('Num') 创建一个符号数字Num

sc=sym('Num') 创建一个符号常数sc，该常数值准确等于Num 说明：Num代表一个具体的数字； Num必须处于（英文状态下的）单引号内，构成字符串（关于字符串参见附录A.1）；

Num须为整数、有理数（如321/1000）、预定义常数（如pi）；

Num应杜绝使用“十进制小数（如0.321）”、“科学法计述数（如321e-3）”，以便确保创建的符号数、符号表达式精准。此规定源于符号计算引擎的特性。

【例2.1-1】

（1）完全精准的符号数值或数值表达式的“推荐”创建法

clear

a=sin(3/10) % 创建一个双精度的数值类常数

sa=sym('sin(3/10)') % 利用有理数创建一个完全精准的符号数值表达式 vpa(sa-a,40) % 计算双精度与精准符号数间的40位精度误差 a =

0.2955 sa =

sin(3/10) ans =

0.00000000000000002892172294171061964864168785515107729858

（2）产生非精准的符号数值或数值表达式创建法（应杜绝使用！）

sa1=sym('sin(0.3)') sa2=sym('sin(3e-1)') d1=vpa(sa-sa1,40 ) sa1 =

0.29552020666133957510532074568503 sa2 =

0.29552020666133957510532074568503 d1 =

-6.349448172324523418097794487226753809579e-41

（3）

whos

Name Size Bytes Class Attributes

a 1x1 8 double ans 1x1 60 sym d1 1x1 60 sym sa 1x1 60 sym sa1 1x1 60 sym sa2 1x1 60 sym

2

三 基本符号变量

经典教科书里，表达式e-axsinbx中的a,b称为参数，x为变量。在MATLAB的符号计算中，a、b、x统称为基本符号变量，其中，x总被默认为自由（待解）符号变量，其他被作为符号参数处理。 定义基本符号变量的指令格式：

para=sym(' para') 定义单个复数域符号变量para

para=sym(' para', 'Flag') 定义单个Flag指定域符号变量para

syms para 定义单个复数域符号变量para的另一种方式

syms para Flag 定义单个Flag指定域符号变量para的另一种方式

syms para1 para2 paraN 定义多个复数域符号变量para1 para2 paraN syms para1 para2 paraN Flag 定义多个Flag指定域符号变量para1 para2 paraN

? 符号参数名不要用处于“字母表中小写字母x及其两侧的英文字母”开头。 ? Flag表示参数属性，可具体取以下词条：

positive 表示那些符号变量取正实数； real 表示那些符号变量限定为实数； unreal 表示那些符号变量为不限定的复数。 ? 默认值是复数域符号变量

? sym指令只能对单变量作用，syms不能用于对数值、常数相关的定义。

四 自由符号变量

解题通常是围绕自由符号变量进行。

确定自由符号变量的规则：

? 在专门指定变量名的符号运算中，解题一定围绕指定变量名进行。 ? 自动识别符号变量时，字母的优先次序为x,y,w,z,v等。

自动识别表达式中自由、独立的符号变量的指令： symvar(expression) 列出表达式中所有自由符号变量

symvar(expression, n) 列出表达式中距离x最近的n个自由符号变量

原来版本是

findsym(EXPR) 确认表达式EXPR中所有自由符号变量

findsym(EXPR, N) 确认表达式EXPR中距离x最近的N个自由符号变量

2【例2.1-2】用符号计算研究方程sin(3)uz?vz?3w-a5?0的解。 （1）产生符号表达式

syms u v w z a5 %定义符号参数和变量 f=sym('3'); %定义符号常数 Eq=sin(f)*u*z^2+v*z+f*w-a5;

3

（2）基本符号变量和自由符号变量的认定

symvar(Eq) % 按字母表顺序认定基本符号变量，注意没有f ans =

[ a5, u, v, w, z]

symvar(Eq,100) %按离x的距离认定所有自由符号变量 ans =

[ w, z, v, u, a5]

symvar(Eq,1) %仅认定一个自由符号变量 ans = w

（3）不指定变量情况

result_1=solve(Eq) result_1 =

a5/3 - (v*z)/3 - (u*sin(3)*z^2)/3

（4）指定变量情况

result_2=solve(Eq,z) result_2 =

-(v - (v^2 + 4*a5*u*sin(3) - 12*u*w*sin(3))^(1/2))/(2*u*sin(3)) -(v + (v^2 + 4*a5*u*sin(3) - 12*u*w*sin(3))^(1/2))/(2*u*sin(3))

【例2.1-3】对符号变量的辨认。

（1）

syms a b x X Y k=sym('3'); z=sym('c*sqrt(d)+y*sin(t)'); EXPR=a*z*X+(b*x^2+k)*Y;

% 定义基本符号变量 % 定义符号常数

% 创建“元”符号表达式 % 构成衍生符号表达式

（2）

symvar(EXPR) ans =

[ X, Y, a, b, c, d, t, x, y]

%列出EXPR中全部基本符号变量，k、z除外

（3）

symvar(EXPR,10) %列出EXPR中全部自由符号变量 ans =

[ x, y, t, d, c, b, a, X, Y]

（4）

symvar(EXPR,1) %列出EXPR中最优先的一个自由符号变量 ans = x

（5） %列出EXPR中最优先的3个自由符号变量 symvar(EXPR,3) ans =

[ x, y, t]

（6）

E3=sym('a*sqrt(theta)') % theta为sqrt函数唯一输入量引起问题 ??? Error using sym>convertExpression (line 2256)

Conversion to 'sym' returned the MuPAD error: Error: Expecting an arithmetical expression. [sqrt]

4

Error in sym>convertChar (line 2167) s = convertExpression(x);

Error in sym>convertCharWithOption (line 2150) s = convertChar(x);

Error in sym>tomupad (line 1881)

S = convertCharWithOption(x,a); Error in sym (line 108)

S.s = tomupad(x,'');

E4=sym('a*sqrt(theta123)') % 创建符号表达式时，回避MATLAB自用关键词 E4 =

a*theta123^(1/2)

【例2.1-4】symvar确定自由变量是对整个矩阵进行的。

syms a b t u v x y

A=[a+b*x,sin(t)+u;x*exp(-t),log(y)+v] symvar(A,1) A =

[ a + b*x, u + sin(t)] [ x*exp(-t), v + log(y)] ans = x

2.1.2

符号计算中的算符

? 与数值计算中的算符在形状、名称和使用方法上几乎完全相同。

2.1.3

符号计算中的函数指令

表2.1-1 MATLAB中可调用的符号计算函数指令 类 别 情况描述 与数值计算对应关系 基本函数 三角函数、双曲函数及反函数；除atan2外 名称和使用方法相同 指数、对数函数（如exp，expm） 名称和使用方法相同 复数函数（注意：没有幅角函数angle） 名称和使用方法相同 矩阵分解函数（如eig，svd） 名称和使用方法相同 方程求解函数solve 不同 微积分函数(如diff，int) 不完全相同 积分变换和反变换函数（如laplace，ilaplace） 只有离散Fourier变换 绘图函数（如ezplot，ezsurf） 数值绘图指令更丰富 经典特殊函数 如误差函数erf、贝塞尔函数besselj、第一类完全无对应函数 椭圆积分EllipticK等；通过mfunlist可以看到所有经典函数名 MuPAD库函数 MuPAD库函数在符号计算的扩展目录上；可通过无对应函数 mhelp index 看到各子函数库的名称；函数的数量很大；使用库函数，需要具备MuPAD语言知识

5

2.1.4 符号对象的识别

为了函数指令与数据对象的适配，MATLAB提供了用于识别数据对象属性的指令： class(var) 给出变量var的数据类别（如double，sym等） isa(var, 'Obj') 若变量var是Obj代表的类型，给出1，表示“真” whos 给出所有MATLAB内存变量的属性

【例2.1-5】数据对象及其识别指令的使用。 （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 2 3 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）

CMn=class(Mn) CMc=class(Mc) CMs=class(Ms) CMn = double CMc = char CMs = sym

（4）

isa(Mn,'double') isa(Mc,'char') isa(Ms,'sym') ans = 1 ans = 1 ans = 1

（5）

whos Mn Mc Ms

Name Size Bytes Class Attributes

% 产生4个数值变量

利用已赋值变量构成数值矩阵

% 字符串中的a,b,c,d与前面输入的数值变量无关 % Ms是一个符号矩阵，它与前面各变量无关

6

%

Mc 1x9 18 char Mn 2x2 32 double Ms 2x2 312 sym

2.1.5 符号运算的工作机理

一、符号运算的工作机理

? sym或syms启动MuPAD引擎，开启MuPAD单独的内存空间 ? 被定义变量保存至Matlab内存空间

? 对变量的限定性假设，保存在MuPAD空间中

二、sym和syms指令可作的限定性假设

x=sym('x') % 定义“复数”符号变量x，默认 x=sym('x','real') % 定义“实数”符号变量x x=sym('x','positive') % 定义“正数”符号变量x

syms x y z % 定义“复数”符号变量x，y，z, 默认 syms x y z real % 定义“实数”符号变量x，y，z syms x y z positive % 定义“正数”符号变量x，y，z

三、清除变量和撤销假设

clear x % 清除MATLAB内存中的x变量

syms x clear % 撤销MuPAD内存中对变量x的任何假设，而恢复为“复数”变量 syms（'x','clear'） % 功能同上

assumptions(x) % 显示符号变量x的限定性假设

assumptions % 显示MuPUD内存中已带限定性假设的全部符号变量

reset(symengine) % 重启MuPUD引擎，清空MuPUD内存中所有内容

【例2.1-7】

clear all reset(symengine) Da=1.2;Dw=1/3; syms sa sw sx sy sz syms A B positive syms C real

（2）

whos

（3）

assumptions

（4）

clear A who

assumptions

（5）

syms B clear who

7

assumptions

2.1.6 符号帮助体系

1、“直接调用符号计算指令”的求助。 不依赖mfun，evalin，feval指令调用。 doc symname

doc laplace

help symname

help laplace

2、借助mfun调用的“特殊函数指令”的求助

doc mfunlist help mfunlist

3、借助evalin，feval指令调用的“MuPAD指令”的求助,需要借助MuPUD帮助浏览器 doc(symengine)

【例2.1-8】 （1）

图2.1-1

8

（2）erfc

图2.1-2

? （3）doc(symengine)

图2.1-3

9

2.2

2.2.1 一

符号数字及表达式的操作

双精度数字与符号数字之间的转换 双精度数字向符号数字的转换

借助sym函数：

sym(Num,'r') 数值类数字Num的“有理分数”表达的符号数字 sym(Num) sym(Num,'r')的简略形式

sym(Num,'f') 数值类数字Num的“两个二进制数”近似表达符号数字 sym(Num,'e') 数值类数字Num的带估计误差的“广义有理表达”符号数字 sym(Num,'d') 数值类数字Num的“十进制浮点”近似表达符号数字

1、Num是“双精度数字”。

2、sym(Num)中Num和sym('Num')中'Num'的差别

Num：用作符号计算函数的输入量时，是理论数字的双精度近似。

'Num'：字符串数字，用作符号计算函数的输入量时，是数字的理论真值。 3、在符号运算中，“双精度数字”会自动地按sym（Num，'r'）格式转换为符号数字。

二 符号数字向双精度数字转换

double(Num_sym) 把符号数字Num_sym转换为双精度数字

2.2.2 符号数字的任意精度表达形式

为了兼顾计算精度和速度，提供的“变精度”算法指令：

digits 显示当前环境下符号数字“十进制浮点”表示的有效数字位数（默认32位） digits(n) 设定符号数字“十进制浮点”表示的有效数字位数 xs=vpa(x) 据表达式x得到digits指定精度下的符号数字xs

xs=vpa(x,n) 据表达式x得到n位有效数字的符号数字xs，仅在运行的当时起作用。

【例2.2-1】digits, vpa指令的使用。 （1）

reset(symengine) % 重新启动符号计算引擎

sa=sym('1/3+sqrt(2)') % 定义准确符号数 sa =

2^(1/2) + 1/3

（2）

digits

Digits = 32

format long

% 观察有效位数

10

a=1/3+sqrt(2) % 定义双精度数

sa_Plus_a=vpa(sa+a,20) % 给出20位有效数字结果 sa_Minus_a=vpa(sa-a,20) % a =

1.747546895706428 sa_Plus_a =

3.4950937914128567869 sa_Minus_a =

-0.000000000000000022658064826339973669

（3）

sa32=vpa(sa) % 采用默认设置的32位有效数字 digits(48) % 设置48位有效数字 sa5=vpa(sa,5) % 仅影响sa5数位，对其后无影响。 sa48=vpa(sa) % 仍为48位有效数字 sa32 =

1.747546895706428382135022057543 sa5 = 1.7475 sa48 =

1.74754689570642838213502205754303141190300520871

2.2.3 符号表达式的基本操作

collect（合并同类项）

factor（进行因式分解） numden（提取公因式）等 最常用：

simple(EXPR) 把EXPR（符号表达式或矩阵）转换成最简形式

【例2.2-2】简化f?31612???8。 x3x2xsyms x

f=(1/x^3+6/x^2+12/x+8)^(1/3) g1=simple(f) f =

(12/x + 6/x^2 + 1/x^3 + 8)^(1/3) g1 =

1/x + 2

2.2.4 一

表达式中的置换操作 公因子法简化表达

RS=subexpr(S) 从S中自动提取公因子sigma，并把采用sigma重写的S赋给RS RS=subexpr(S,'w') 从S中自动提取公因子，记为w，并把采用w重写的S赋给RS [RS,w]=subexpr(S,'w') 该调用格式的效果与RS=subexpr(S, 'w')相同 S：符号表达式、符号表达式矩阵

11

【例2.2-3】对符号矩阵?（１）

?ab?进行特征向量分解。 ??cd?clear % 清空所有内存变量

A=sym('[a b;c d]') % 经字符串直接定义符号矩阵

[V,D]=eig(A) % 符号矩阵的特征值、特征向量分解 A =

[ a, b] [ c, d] V =

[(a/2+d/2-a^2-2*a*d+d^2+4*b*c)^(1/2)/2)/c-d/c, (a/2+d/2+(a^2-2*a*d+ d^2+4*b*c)^(1/2)/2)/c-d/c] [ 1, 1] D = [a/2 + d/2 - (a^2 - 2*a*d + d^2 + 4*b*c)^(1/2)/2, 0] [ 0, a/2 + d/2 + (a^2 - 2*a*d + d^2 + 4*b*c)^(1/2)/2]

（2）

subexpr([V;D]) % 自动提取公因式 who % 列出工作内存中的变量 sigma =

(a^2 - 2*a*d + d^2 + 4*b*c)^(1/2) ans =

[ (a/2 + d/2 - sigma/2)/c - d/c, (a/2 + d/2 + sigma/2)/c - d/c] [ 1, 1] [ a/2 + d/2 - sigma/2, 0] [ 0, a/2 + d/2 + sigma/2]

Your variables are:

A D V ans sigma

（3）

Dw=subexpr(D,'w') % 把自动提取的公因式记为w,Dw是用w重记D后的表达 w =

(a^2 - 2*a*d + d^2 + 4*b*c)^(1/2) Dw =

[ a/2 + d/2 - w/2, 0] [ 0, a/2 + d/2 + w/2]

（4）

[RVD,w]=subexpr([V;D],'w') % <7> % 给出合成矩阵[V;D]的公因式表达方式 RVD =

[ (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/2] w =

(a^2 - 2*a*d + d^2 + 4*b*c)^(1/2)

二 通用置换指令

RES=subs(ES,old,new) 用new置换ES中的old后产生RES

RES=subs(ES,new) 用new置换ES中的自由变量后产生RES

12

【例2.2-4】用简单算例演示subs的置换规则。 （1）产生符号函数

syms a x;f=a*sin(x)+5

f =

a*sin(x)+5

（2）符号表达式置换

f1=subs(f,'sin(x)',sym('y')) class(f1)

f1 = % 符号表达式 a*y+5 ans = sym

%<2>

（3）符号常数置换

f2=subs(f,{a,x},{2,sym('pi/3')}) % <3> a被双精度数字置换，x被符号数字置换 class(f2)

f2 = % 3^(1/2)+5 ans = sym

（4）双精度数值置换

f3=subs(f,{a,x},{2,pi/3}) class(f3)

f3 = % 双精度数值变量 6.7321 ans = double

%<4>

（5）数值数组置换之一

f4=subs(subs(f,a,2),x,0:pi/6:pi) %<5> class(f4)

f4 = % 双精度数值变量

5.0000 6.0000 6.7321 7.0000 6.7321 6.0000 5.0000 ans = double

（6）数值数组置换之二

f5=subs(f,{a,x},{0:6,0:pi/6:pi}) %<6> class(f5)

f5 = % 双精度数值变量

5.0000 5.5000 6.7321 8.0000 8.4641 7.5000 5.0000 ans = double

13

2.3

2.3.1

符号微积分

极限和导数的符号计算

大学本科高等数学中的大多数微积分问题，都能用符号计算解决，手工笔算演绎的烦劳都可以由计算机完成。

limit(f,v,a) 求极限 limf(v)

v?a

?1?【例2.3-1】试求lim?1??。

x???x?kxsyms x k

Lim_f=limit((1-1/x)^(k*x),x,inf) Lim_f = exp(-k)

dnfdiff(f,v,n) 求 (n缺省时，默认n=1) ndv?at3?dfd2fd2f【例2.3-2】求f???求dx， dt2，dtdx。

?tcosxlnx?syms a t x;f=[a,t^3;t*cos(x), log(x)]; df=diff(f) dfdt2=diff(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]

14

2.3.2 序列/级数的符号求和

symsum(f,v,a,b) 求f在变量v取遍[a, b]中所有整数时的和。a,b缺省时默认求和区间[0, v-1]。

?1(?1)k?,【例2.3-8】求?[t,k]，???。 2(2k?1)kk?1?t?0?t?13?syms k t;f1=[t k^3];f2=[1/(2*k-1)^2,(-1)^k/k];

s1=simple(symsum(f1)) % f1的自变量被确认为t s2=simple(symsum(f2,1,inf)) % f1的自变量被确认为k s1 =

[ 1/2*t*(t-1), k^3*t] s2 =

[ 1/8*pi^2, -log(2)]

2.3.3 符号积分

intf=int(f,v) 求f对变量v的不定积分 intf=int(f,v,a,b) 求f对变量v的定积分

【例2.3-9】求

11?x?xxdx。

clear syms x

f=sqrt((1+x)/x)/x s=int(f,x) f =

((x + 1)/x)^(1/2)/x s =

- atan((1/x + 1)^(1/2)*i)*2*i - 2*(1/x + 1)^(1/2)

【例2.3-11】求积分

??? 1 2 x2 x2y x xy(x2?y2?z2)dzdydx。

syms x positive syms y z

F2=int(int(int(x^2+y^2+z^2,z,sqrt(x*y),x^2*y),y,sqrt(x),x^2),x,1,2) VF2=vpa(F2) %积分结果用32位数字表示 F2 =

(14912*2^(1/4))/4641 - (6072064*2^(1/2))/348075 + (64*2^(3/4))/225 + 1610027357/6563700 VF2 =

224.92153573331143159790710032805

15

本文来源：https://www.bwwdw.com/article/2jhx.html