精华区文章阅读

发信人: bage (最近比较烦), 信区: Matlab
标  题: 精通Matlab(十八)
发信站: 哈工大紫丁香 (Sun Feb  4 14:04:40 2001), 转信

发信人: Security (淼水), 信区: MathTools
发信站: BBS 水木清华站 (Tue Jun  1 22:43:29 1999)

第22章　符号数学工具

MATLAB所具有的符号数学工具箱与其它所有工具不同，它适用于广泛的用途，而不是针对一些特殊专业或专业分支。另外，MATLAB符号数学工具箱与其它的工具箱区别还因为它使用字符串来进行符号分析，而不是基于数组的数值分析。为此，本章包含了该工具箱的教学辅导材料。

22.1  引言

符号数学工具箱是操作和解决符号表达式的符号数学工具箱(函数)集合，有复合、简化、微分、积分以及求解代数方程和微分方程的工具。另外还有一些用于线性代数的工具，求解逆、行列式、正则型式的精确结果，找出符号矩阵的特征值而无由数值计算引入的误差。工具箱还支持可变精
度运算，即支持符号计算并能以指定的精度返回结果。
符号数学工具箱中的工具是建立在功能强大的称作Maple软件的基础上。它最初是由加拿大的滑铁卢（Waterloo）大学开发的。当要求MATLAB进行符号运算时，它就请求Maple去计算并将结果返回到MATLAB命令窗口。因此，在MATLAB中的符号运算是MATLAB处理数字的自然扩展。

22.2  符号表达式

符号表达式是代表数字、函数、算子和变量的MATLAB字符串，或字符串数组。不要求变量有预先确定的值，符号方程式是含有等号的符号表达式。符号算术是使用已知的规则和给定符号恒等式求解这些符号方程的实践，它与代数和微积分所学到的求解方法完全一样。符号矩阵是数组，其元
素是符号表达式。
MATLAB在内部把符号表达式表示成字符串，以与数字变量或运算相区别；否则，这些符号表达式几乎完全象基本的MATLAB命令。表22.1列有几则符号表达式例子以及MATLAB等效表达式。

表22.1
符号表达式      MATLAB表达式
         ' 1/(2*x^n) '
        y= ' 1/sqrt(2*x) '
         ' cos(x^2)-sin(2*x) '
        M=sym( ' [a，b；c，d] ' )
        f=int( ' x^3/sqrt(1-x) ' ， ' a ' ， ' b ' )

MATLAB符号函数可让用户用多种方法来操作这些表达式，比如，

>> diff( ' cos(x) ' )  %  differentiate cos(x) with respect to x
ans=
        -sin(x)

>> M=sym( ' [a，b；c，d] ' )  %  create a symbolic matrix M
M=
        [a，b]
        [c，d]

>> determ(M)  %  find the determinant of the symbolic matrix M
ans=
        a*d-b*c

请注意，上面的第一个例子的符号表达式是用单引号以隐含方式定义的。它告诉MATLAB  ' cos(x) ' 是一个字符串并说明diff( ' cosx ' )是一个符号表达式而不是数字表达式；然而在第二个例子中，用函数sym显式地告诉MATLAB M=sym( ' [a，b；c，d] '
)是一符号表达式。在MATLAB可以自己确定变量类型的场合下，通常不要求显式函数sym。
正如在第八章所阐述，MATLAB中函数function argument形式是与function( ' argument ' )等价的。其中，function是一个函数，argument是一字符串。例如，MATLAB可以构造diff cos(x)和diff( ' cos(x) ' )两者都意味diff (sym ' cos(x) '
)。但第一种形式显然更便于输入。然而，很多时候sym是必要的。在上述的第二个例子中，

>> M=[a，b；c，d]  %  M is a numeric matrix using value of a through d
???Uundefine function or variable a.

>> M= ' [a，b；c，d] '   %  M is a character string， but not a symbolic matrix
M=
        [a，b；c，d]

>> M=sym( ' [a，b；c，d] ' )  %  M is a symbolic matrix
M=
        [a，b]
        [c，d]

M以三种方式定义：数字型(如果a、b、c、d已预先确定)、字符串型或符号矩阵型。许多符号函数非常巧妙能够自动将字符转变为符号表达式。但在某些情况下，尤其是建立符号数组时，必须用函数sym，特别地将字符串变为符号表达式。隐含形式，例如diff
cos(x)，对于那些不需要参考先前结果的简单任务，最有用。但是最简单形式(无引号)要求一个参量，它是一个单字符的字符串、不包含插入的空格。

>> diff x^2+3*x+5  %  the argument is equivalent to  ' x^2+3*x+5 '
ans=
        2*x+3

>> diff x^2 + 3*x + 5  %  spaces break the argument into separate strings
???Error using==>diff
Too manyinput arguments

无变量的符号表达式称作符号常量。符号常量常常与整数很难区别，例如

>> f=symop( ' (3*4-2)/5+1 ' )  %  reduce a symbolic constant to its simplest form
f =
        3

>> isstr(f) %  is f a string? (1=yes， 0=no)
ans=
        1

在这个例子中，f代表符号常数 ' 3 ' ；而不是数字3。正如第六章所阐述的，MATLAB是以字符ASCII码形式来存储字符串的。所以，如果对字符串进行数字运算，则在运算中，采用各字符串的ASCII码值。因为数字51是字符 ' 3 ' 的ASCII表示，所以f加1在数值上不能得到期望的结果

>> f+1
ans=
        52

符号变量

当字符表达式中含有多于一个的变量时，只有一个变量是独立变量。如果不告诉MATLAB哪一个变量是独立变量，MATLAB将基于以下规则选择一个：

在符号表达式中缺省的独立变量是唯一的，除去i和j的小写字母，不是单词的一部分。如果没有这种字母，就选择x作为独立变量。如字符不是唯一的，就选择在字母顺序中最接近x的字母。如果有相连的字母，就选择在字母表中较后的那一个。

缺省的独立变量，有时称作自由变量，在表达式 ' 1/(5+cos(x)) ' 中是 ' x ' ；在 ' 3*y+z ' 中是 ' y ' ；在 ' a+sin(t) ' 是 ' t ' 。在表式 ' sin(pi/4)-cos(3/5) ' 中自由符号变量是 ' x '
，因为此式是一个符号常数无符号变量。可利用函数symvar询问MATLAB在符号表达式中哪一个变量它认为是独立变量。

>> symvar( ' a*x+y*)  %  fi
f=
        a*x^2+b*x+c

>> subs(f，' s '，' x ')  %  substitute  ' s '  for  ' x '  in the expression f
ans=
        a*s^2+b*s+c

>> subs(f，' alpha '，' a ')  %  substitute  ' alpha '  for  ' a '  in  f
ans=
        alpha*x^2+b*x+c

>> g=' 3*x^2+5*x-4 '  %  create another function
g=
        3*x^2+5*x-4

>> h=subs(g，' 2 '，' x ')  %  substitute  ' 2 ' for  ' x '  in  g
h=
        18

>> isstr(h)  %  show that the result is a symbolic expression
ans=
        1

最后一个例子表明subs如何进行替换，并力图简化表达式。因为替换结果是一个符号常数，MATLLAB可以将其简化为一个符号值。注意，因为subs是一个符号函数，所以它返回一个符号表达式。尽管看似数字，实质上是一个符号常数。为了得到数字，我们需要使用函数numeric或eval来转换
字符串。

>> numeric(h)  %  convert a symbolic expression to a number
ans=
        18

>> isstr(ans)  %  show that the result is a numeric value
ans=
        0

22.4  微分和积分

微分和积分是微积分学研究和应用的核心，并广泛地用在许多工程学科。MATLAB符号工具能帮助解决许多这类问题。

微分

符号表达式的微分以四种形式利用函数diff：

>> f= ' a*x^3+x^2-b*x-c '   %  define a symbolic expression
f=
        a*x^3+x^2-b*x-c

>> diff(f)  %  differentiate with respect to the default variable x
ans=
        3*a*x^2+2*x-b

>> diff(f，'a ')  %   differentiate with respect to a
ans=
        x^3

>> diff(f，2)  %   differentiate twice with respect to x
ans=
        6*a*x+2

>> diff(f，' a '，2)  %   differentiate twice with respect to a
ans=
        0

函数diff也可对数组进行运算。如果F是符号向量或数组，diff(F)对数组内的各个元素进行微分。

>> F=sym(' [a*x，  b*x^2；  c*x^3，  d*s] ')  %  create a symbolic array
F=
        [   a*x，  b*x^2]
        [c*x^3，     d*s]

>> diff(F)    %   differentiate the element with respect to x
ans=
        [           a，2*b*x]
        [3*c*x^2，       0]

注意函数diff也用在MATLAB，计算数值向量或矩阵的数值差分。对于一个数值向量或矩阵M，diff(M)计算M(2： m，： )-M(1： m-1，： )的数值差分，如下所示：

>> m=[(1： 8).^2)]  %  create a vector
M=
        1   4   9   16   25   36   49   64

>> diff(M)  %  find the differences between elements
ans=
        3   5   7   9   11   13   15

如果diff的表达式或可变参量是数值，MATLAB就非常巧妙地计算其数值差分；如果参量是符号字符串或变量，MATLAB就对其表达式进行微分。

积分

积分函数int(f)，其中f是一符号表达式，它力图求出另一符号表达式F使diff(F)=f。正如从研究微分学所了解的，积分比微分复杂得多。积分或逆求导不一定是以封闭形式存在，或许存在但软件也许找不到，或者软件可明显地求解，但超过内存或时间限制。当MATLAB不能找到逆导数时，
它将返回未经计算的命令。

>> int( ' log(x)/exp(x^2) ' )  %  attempt to integrate
ans=
        log(x)/exp(x^2)

同微分一样，积分函数有多种形式。形式int(f)相对于缺省的独立变量求逆导数；形式(f，' s ')相对于符号变量s积分；形式int(f，a，b)和int(f，' s '，a，b)，a，b是数值，求解符号表达式从a到b的定积分；形式int(f，' m ' ，' n ')和形式int(f，' s '，' m '，' n
')，其中m，n是符号变量，求解符号表达式从m到n的定积分。

>> f=' sin(s+2*x) ')  %  crate a symbolic function
f=
        sin(s+2*x)

>> int(f)  %  integrate with respect to x
ans=
        -1/2*cos(s+2*x)

>> int(f，' s ')   %  integrate with respect to s
ans=
        -cos(s+2*x)

>> int(f，pi/2，pi)  %  integrate with respect to x from  /2 to
ans=
        -cos(x)

>> int(f，' s '，pi/2，pi)  %  integrate with respect to s from  /2 to
ans=
        cos(2*x)-sin(2*x)

>> int(f，' m '，' n ')   %  integrate with respect to x from m to n
ans=
        -1/2*cos(s+2*n)+1/2*cos(s+2*m)

正如函数diff一样，积分函数int对符号数组的每一个元素进行运算。

>> F=sym( ' [a*x，b*x^2；c*x^3，d*s] ' )  %  create a symbolic array
F=
        [ a*x，b*x^2]
        [c*x^3， d*s]

>> diff(F)  %  ubtegrate the array elements with respect to x
ans=
        [1/2*a*x^2，1/3*b*x^3]
        [1/4*c*x^4，    d*s*x]

22.5  符号表达式画图

在许多的场合，将表达式可视化是有利的。MATLAB提供了函数ezplot来完成该任务。

>> y=' 16*x^2+64*x+96 '  %  expression to plot
y=
        16*x^2+64*x+96

>> ezplot(y)

图22.1  符号函数16*x^2+64*x+96  (-2 ≤x≤2 )

图22.2  符号函数16*x^2+64*x+96   0≤x≤6

正如图22.1所示，ezplot绘制了定义域为-2 ≤x≤2 的给定符号函数，并相应地调整了y轴比例，还加了网格栅和标志。在这个例子中，我们感兴趣的时间是从0到6。让我们再试一下，并指定时间范围，见图22.2

>> ezplot(y，[0 6])  %  plot y for 0≤x≤6

现在，在所感兴趣的范围内显示得略好些。一旦该图处在图形窗口内，它可以象其它图象一样作修改。

22.6  符号表达式简化和格式化

有时MATLAB返回的符号表达式难以理解，有许多工具可以使表达式变得更易读懂。第一个就是函数pretty，该命令以类似于数学课本上的形式来显示符号表达式。

>> f=taylor( ' log(x+1)/(x-5) ' )  %   6 terms is the default
f =
        -1/5*x+3/50*x^2-41/750*x^3+293/7500*x^4-1207/37500*x^5+O(x^6)

>> pretty(f)
                                   2     41   3     293   4     1207   5         6
            - 1/5 x + 3/50 x  -  -----x   + ------x   -  -------x  + O(x )
                                        750         7500        37500

符号表达式可用许多等价形式来提供。在不同的场合，某种形式可能胜于另一种。MATLAB用许多命令来简化或改变符号表达式。

>> f=sym( ' (x^2-1)*(x-2)*(x-3) ' )  %  create a function
f=
        (x^2-1)*(x-2)*(x-3)

>> collect(f)  %  collect all like terms
ans=
        x^4-5*x^3+5*x^2+5*x-6

>> hornor(ans)   %  change to Horner or nested representation
ans=
        -6+(5+(5+(-5+x)*x)*x)*x

>> factor(ans)  %  express as a product of polynomials
ans=
        (x-1)*(x-2)*(x-3)*(x+1)

>> expand(f)  %  distribute products over sums
ans=
        x^4-5*x^3+5*x^2+5*x-6

Simplify是功能强大、通用的工具。它利用各种类型代数恒等式，包括求和、积分和分数幂、三角、指数和log函数、Bessel函数、超几何函数和函数，来简化表达式。
下面例子说明函数的乘幂。

>> simplify( ' log(2*x/y) ' )
ans=
        log(2)+log(x)-log(y)

>> simplify( ' sin(x)^2+3*x+cos(x)^2-5 ' )
ans=
        3*x-4

>> simplify( ' (-a^2+1)/(1-a) ' )
ans=
        a+1

其中要讨论的最后一个函数是最有用的，但也是最不正统的。函数simple试用了几种不同的简化工具，然后选择在结果表达式中含有最少字符的那种形式。让我们看一下立方根：

>> f=' (1/x^3+6/x^2+12/x+8)^(1/3) '  %  create the expression
f=
        (1/x^3+6/x^2+12/x+8)^(1/3)

>> simple(f)  %  simplify it
simplify ：
        (2*x+1)/x
ans=
        (2*x+1)/x

>> simple(ans)  %  try it once again - another method may help
combine(trig)
        2+1/x
ans=
        2+1/x

正如所见，simple试用了几种可简化表达式的简化方式，并让看到每一个尝试的结果。有时，它帮助多次使用函数simple并对第一次的结果作不同的简化操作，如上所作。simple对于含有三角函数的表达式尤为有用。让我们试一下

>> simple( ' cos(x)+sqrt(-sin(x)^2) ' )  %  simplify a trig expression
simplify：
        cos(x)+(cos(x)^2-1)^(1/2)
radsimp：
        cos(x)+i*sin(x)
combine(trig)：
        cos(x)+(-1/2+1/2*cos(2*x))^(1/2)
factor：
        cos(x)+(-sin(x)^2)^(1/2)
expand：
        cos(x)+(-sin(x)^2)^(1/2)
convert(exp)：
        1/2*exp(i*x)+1/2/exp(i*x)+1/4*4^(1/2)*(exp(i*x)-1/exp(i*x))
convert(tan)：
        (1-tan(1/2*x)^2)/(1+tan(1/2*x)^2)+(-4*tan(1/2*x)^2/(1+tan(1/2*)^2)^2)^(1/2)
ans =
        cos(x)+i*sin(x)

--
※ 修改:·Security 於 Jun  1 22:44:02 修改本文·[FROM:   166.111.167.4]

--
☆ 来源:．哈工大紫丁香 bbs.hit.edu.cn．[FROM: bage.bbs@smth.org]

Matlab 版 (精华区)