1.6 MATLAB的符號運算

MATLAB的優點不僅在于其強大的數值運算功能，而且也在于其強大的符號運算功能。MATLAB的符號運算是通過集成在MATLAB中的符號數學工具箱（Symbolic MathToolbox）來實現的。它可完成幾乎所有符號表達式的運算功能，如符號表達式的生成、復合和化簡，符號矩陣的求解，符號微積分的求解，符號函數的畫圖，符號代數方程的求解，以及符號微分方程的求解等。

1.6.1 符號表達式的生成

在MATLAB中的符號數學工具箱中，符號表達式是指代表數字、函數和變量的MATLAB字符串或字符串數組，它不要求變量有預先確定的值。

符號表達式可以是符號函數或符號方程。其中，符號函數沒有等號，而符號方程必須有等號。MATLAB在內部把符號表達式表示成字符串，以便與數字相區別。符號表達式可由以下三種方法生成。

1.用單引號生成符號表達式

在MATLAB中，符號表達式如同字符串一樣也可利用單引號來直接設定。例如：

      ＞＞fun=′sin（x）′

結果顯示：

      fun=
          sin（x）
      ＞＞fun=′a*x^2+b*x+c=0′

結果顯示：

      fun=
          a*x^2+b*x+c=0

2.用sym（）函數生成符號表達式

在MATLAB可自動確定變量類型的情況下，不用sym（）函數來顯式地生成符號表達式。但在某些情況下，特別是在建立符號數組時，必須要用sym（）函數來將字符串轉換成符號表達式。例如：

      ＞＞A=sym（′[sin（x）b；c d]′）

結果顯示：

      A=
        [ sin（x），b]
        [c，d]

其結果表示A為一個2 ×2維的符號矩陣。但命令

      ＞＞A=′[sin（x）b；c d]′

的結果為

      A=
        [sin（x）b；c d]

其結果表示A為一字符串。

3.用命令syms生成符號表達式

在MATLAB中，利用命令syms只能生成符號函數，而不能生成符號方程。例如

      ＞＞syms K t T；fun=K*（exp（-t/T））

結果顯示：

      fun=
          K*exp（-t/T）

另外，在MATLAB中，利用symvar（）函數可知道符號表達式中哪些變量為符號變量。同時MATLAB會自動把i，j，pi，inf，nan，eps等特殊字母不當成符號變量。例如，

      ＞＞symvar（′5*pi+j*K*（exp（-t/T）′）

結果顯示：

            ′K′
            ′T′
            ′t′

1.6.2 符號表達式的基本運算

在MATLAB的符號工具箱中，利用相關函數，可對符號表達式進行分子/分母的提取、基本代數運算、相互轉換、化簡和替換等基本運算。

1.符號表達式的提取分子/分母運算

在MATLAB中，如果符號表達式為有理分式的形式或可展開為有理分式的形式，則可通過numden（）函數來提取符號表達式中的分子與分母。其調用格式為

[num，den] =numden（f）

其中，f表示所求符號表達式；num和den表示返回所得的分子與分母。例如：

      ＞＞f=sym（′（x+d）/（a*x^2+b*x+c）′）；[num，den] =numden（f）

運行結果：

      num=
          x+d
      den=
          a*x^2+b*x+c

2.符號表達式的基本代數運算

在MATLAB中，符號表達式的加、減、乘、除四則運算及冪運算等基本的代數運算，分別由symadd（），symsub（），symmul（），symdiv（）及sympow（）函數來實現。其中求和函數symadd（）的調用格式為

h=symadd（f，g）

式中，f，g表示待運算的符號表達式；h表示結果符號表達式。其中，當f，g為符號矩陣時，以上四則運算及冪運算的命令仍然成立。以上其他函數的調用格式均與求和函數symadd（）的調用格式一致。

3.符號表達式與數值表達式的相互轉換

在MATLAB中，利用numeric（）函數（僅適用于MATLAB 6.5及以前的版本）或eval（）函數可將符號表達式轉換成數值表達式。反之，sym（）函數可將數值表達式轉換成符號表達式。例如：

      ＞＞f=′abs（-1）+sqrt（1）/2′，p=eval（f），n=sym（p）

運行結果：

      f=
        abs（-1）+sqrt（1）/2
      p=
        1.5000
      n=
        3/2

若已知數值多項式系數向量，則可以通過符號運算工具箱提供的函數poly2sym（）將其轉換成多項式表達式。若已知多項式表達式，則可以由函數sym2poly（）將其轉換成系數向量形式。它們調用格式為

f=poly2sym（p）和p=sym2poly（f）

其中，p為由多項式系數按降冪排列構成的系數向量；f為多項式表達式。

      ＞＞syms x；p=sym2poly（x^2+3*x+2），f=poly2sym（p）

運行結果：

      p=
        1    3    2
      f=
        x^2+3*x+2

4.符號表達式的化簡

在MATLAB中，simple（）函數可按有關數學規則把符號表達式化簡成最簡形式，其調用格式為

y=simple（f）

其中，f表示化簡前的符號表達式；y表示按有關規則化簡后的符號表達式。例如：

      ＞＞f=sym（′2*sin（x）*cos（x）′），y=simple（f）

運行結果：

      f=
        2*sin（x）*cos（x）
      y=
        sin（2*x）

另外，在MATLAB的符號數學工具箱中提供的符號表達式化簡函數還有：

      pretty（）   %將符號表達式轉換成與公式編輯器顯示的符號表達式相類似的形式；
      collect（）  %將符號表達式的同類項進行合并；
      horner（）   %將一般的符號表達式轉換成嵌套形式的符號表達式；
      factor（）   %對符號表達式進行因式分解；
      expand（）   %對符號表達式進行展開；
      simplify（） %利用各種類型的代數恒等式對符號表達式進行化簡。

5.符號表達式的替換

在MATLAB的符號數學工具箱中，subexpr（）函數和subs（）函數可以進行符號表達式的替換。其中subexpr（）函數用于把復雜表達式中所含的多個相同子表達式用一個符號代替，使其表達簡潔，其調用格式為

g=subexpr（f，′S′）

其中，f，g分別表示置換前后的符號表達式；S表示置換復雜表達式中子表達式的符號變量。復雜表達式中被置換的子表達式是自動尋找的，只有比較長的子表達式才被置換，至于比較短的子表達式，即便多次重復出現，也不被置換。例如：

      ＞＞f=solve（′x^3+a*x+1′）；r=subexpr（f，′ss′）
      ＞＞f=solve（′a*x^2+b*x+c′）；r=subexpr（f，′ss′）

subs（）函數除具有與subexpr（）函數一樣的可以用一個符號變量替換復雜表達式中所含的多個相同子表達式的作用外，還可以求解被替換的復雜符號表達式的值，其調用格式為

      g=subs（f，old，new）  %表示用new置換符號表達式f中的old后產生g；
      g=subs（f，new）     %表示用new置換符號表達式f中的自由變量后產生g。

例如，對于符號表達式f（x）=asin（x）+5，可進行下列替代運算：

      ＞＞syms a x；f=a*sin（x）+5；g1=subs（f，′sin（x）′，′y′），g2=subs（f，a，9）
      ＞＞g3=subs（f，{a，x}，{2，sym（pi/3）}），g4=subs（f，{a，x}，{2，pi/3}）
      ＞＞g5=subs（subs（f，a，2），x，0：pi/6：pi），g6=subs（f，{a，x}，{0：6，0：pi/6：pi}）

1.6.3 符號表達式的微積分

MATLAB的符號工具箱中，符號表達式的微積分包括符號序列求和、符號極值、符號微分和符號積分等運算。

1.符號序列求和

對于求和問題，在MATLAB中可利用符號序列求和函數symsum（）來實現，其調用格式為

      y=symsum（f，′x′，a，b） %求符號表達式f在指定變量x取遍[a，b]中所有整數和y；
      y=symsum（f，′x′）       %求符號表達式f在指定變量x取遍[0，x-1]中所有整數和y；
      y=symsum（f，a，b）        %求符號表達式f對獨立變量從a到b的所有整數和y。

【例1-23】 求和的值。

解 MATLAB命令如下

      ＞＞syms t k；f1=[t，k^3]；f2=[1/（2*k-1）^2，（-1）^k/k]；
      ＞＞y1=simple（symsum（f1，′t′）），y2=simple（symsum（f2，1，inf））

運行結果：

      y1=
            [1/2*t*（t-1），k^3*t]
      y2=
            [1/8*pi^2，-log（2）]

2.符號極值

在MATLAB中，符號極值由函數limit（）來實現，其調用格式為：

      y=limit（f，x，a）            %求符號表達式f對變量x趨于a時的極值y；
      y=limit（f，a）               %求符號表達式f對獨立變量趨于a時的極值y；
      y=limit（f）                  %求符號表達式f對獨立變量趨于0時的極值y；
      y=limit（f，x，a，′right′） %求符號表達式f對變量x從右邊趨于a時的極值y；
      y=limit（f，x，a，′left′）  %求符號表達式f對變量x從左邊趨于a時的極值y。

【例1-24】 求符號表達式f（x）=（a-2x）2 x，在x趨于a/6的極值。

解 MATLAB命令如下

      ＞＞syms a x；f=（a-2*x）^2*x；y=limit（f，1/6*a）

運行結果：

      y=
        2/27*a^3

【例1-25】 求二元函數的極限值。

解MATLAB命令如下

      ＞＞syms a x y；f=exp（-1/（y^2+x^2））*sin（x）^2/x^2*（1+1/y^2）^（x+a^2*y^2）；
      ＞＞L=limit（limit（f，x，1/sqrt（y）），y，inf）

運行結果：

      L=
        exp（a^2）

3.符號微分

在MATLAB中，符號微分由diff（）函數來實現。diff（）函數可同時計算數值微分與符號微分，MATLAB能根據其輸入參數的類型（數值或符號字符串），自動對其進行數值微分或符號微分。其調用格式為

      y=diff（f）          %求符號表達式f對獨立變量的微分y；
      y=diff（f，n）       %求符號表達式f對獨立變量的n次微分y；
      y=diff（f，′x′）   %求符號表達式f對變量x的微分y；

y=diff（f，′x′，n）%求符號表達式f對變量x的n次微分y。

【例1-26】 求f（x）=3ax-x3的一階微分和二階微分。

解 MATLAB命令如下

      ＞＞f=′3*a*x-x^3′；df=diff（f，′x′），ddf=diff（f，′x′，2）

運行結果：

      df=
          3*a-3*x^2
      ddf=
          -6*x

對于多元函數的偏導數也可以嵌套使用函數diff（）來求取。

【例1-27】 已知二元函數f（x，y）=（x2 -2x）e -x2-y2-xy，試求?f/?x，?f/?y和?y/?x。

解 MATLAB命令如下

      ＞＞syms x y；f=（x^2-2*x）*exp（-x^2-y^2-x*y）；
      ＞＞dfdx=simple（diff（f，x）），dfdy=diff（f，y），dydx=simple（diff（f，x））/diff（f，y）

運行結果：

      dfdx=
            -exp（-x^2-y^2-x*y）*（-2*x+2+2*x^3+x^2*y-4*x^2-2*x*y）
      dfdy=
         （x^2-2*x）*（-2*y-x）*exp（-x^2-y^2-x*y）
      dydx=
            -（-2*x+2+2*x^3+x^2*y-4*x^2-2*x*y）/（x^2-2*x）/（-2*y-x）

除了利用以上函數求解偏微分方程外，MATLAB還提供了偏微分方程工具箱，可以比較規范地求解各種常見的二階偏微分方程。在MATLAB環境下鍵入pdetool，將啟動偏微分方程求解界面。另外也可以利用Simulink求解微分方程。由于篇幅原因，這里不再介紹，具體內容可參考有關文獻。

4.符號積分

在MATLAB中，符號積分由int（）函數來實現。因為積分比微分要復雜得多，故在很多情況下，積分不一定能成功。當MATLAB進行符號積分而找不到原函數時，它將返回未經計算的函數。int（）函數的調用格式為

      y=int（f）          %求符號表達式f對獨立變量的不定積分y；
      y=int（f，′x′）   %求符號表達式f對變量x的不定積分y；
      y=int（f，a，b）    %求符號表達式f對獨立變量從a到b的定積分y；
      y=int（f，′x′，a，b）%求符號表達式f對變量x從a到b的定積分y。

【例1-28】 試求以下函數的定積分：

解 通過下面的MATLAB語句可求出所需函數的定積分：

      ＞＞syms a；f=′a/sqrt（2*pi）*exp（-x^2/2）′；y=int（f，′x′，-inf，inf）
      ＞＞syms a；y=int（′a/sqrt（2*pi）*exp（-x^2/2）′，′x′，-inf，inf）

或

結果顯示：

      y=
        a

【例1-29】 求積分方程。

解 MATLAB命令如下

      ＞＞syms x y z；
      ＞＞f=int（int（int（′x^2+y^2+z^2′，′z′，sqrt（x*y），x^2*y），′y′，sqrt（x），x^2），′x′，1，2）
      ＞＞p=eval（f）

結果顯示：

      f=
        1610027357/6563700-6072064/348075*2^（1/2）+14912/4641*2^（1/4）+64/225*2^（3/4）

      p=
        224.9215

1.6.4 符號表達式的變換

MATLAB的符號工具箱中，可以進行以下幾種符號表達式的變換。

1.Laplace變換及其反變換

在MATLAB中，給出了求解Laplace變換及其反變換的函數laplace（）和ilaplace（），其調用格式分別為

F=laplace（f，t，s）和f=ilaplace（F，s，t）

其中，f表示時域函數f（t）；t表示時間變量；F表示頻域函數F（s）；s表示頻域變量。

      ＞＞syms k t s；f=k*t^0；F=laplace（f，t，s），f1=ilaplace（F）

結果顯示：

      F=
          k/s
      f1=
          k
      ＞＞syms k T t s；f=k*exp（-t/T）；F=laplace（f，t，s），f1=ilaplace（F）

結果顯示：

      F=
          k/（s+1/T）
      f1=
          k*exp（-t/T）

2.Z變換及其反變換

在MATLAB中，給出了求解Z變換及其反變換的函數ztrans（）和iztrans（），其調用格式分別為

F=ztrans（f，n，z）和f=iztrans（F，z，n）

其中，f表示時域序列f（n）或時間函數f（t）；n表示時間序列；F表示Z域函數F（z）；z表示Z域變量。

      ＞＞syms k t z；f=k*t^1；F=ztrans（f，t，z），f1=iztrans（F）

結果顯示：

      F=
          k*z/（z-1）^2
      f1=
          k*n

3.Fourier變換及其反變換

在MATLAB中，給出了求解Fourier變換及其反變換的函數fourier（）和ifourier（），其調用格式分別為

F=fourier（f，t，ω）和f=ifourier（F，ω，t）

其中，f表示時間函數f（t）；t表示時間變量；F表示頻域函數F（ω）；ω表示頻域變量。

1.6.5 符號表達式的求解

MATLAB的符號工具箱中，符號表達式的求解包括符號代數方程的求解和符號微分方程的求解等。

1.符號代數方程求解

在MATLAB中，符號代數線性方程、符號代數非線性方程，以及符號超越方程均可利用solve（）函數對其進行求解。solve（）函數的調用格式為

                [x，y，z，…] =solve（′eq1′，′eq2′，′eq3′，…，′a′，′b′，′c′，…）

或

                [x，y，z，…] =solve（′eq1，eq2，eq3，…′，′a，b，c，…′）
                [x，y，z，…] =solve（exp1，exp2，exp3，…，′a′，′b′，′c′，…）

其中，eq1，eq2，eq3，…表示符號方程，或不含“等號”的符號表達式（或稱為函數）（此時函數是對eq1=0，eq2=0，eq3=0，…求解）；exp1，exp2，exp3，…僅表示符號表達式；a，b，c，…是符號方程的求解變量名；x，y，z，…是符號方程的解賦值的變量名。

【例1-30】 求方程3x2 -3a2 =0的解。

解 MATLAB命令如下：

      ＞＞x=solve（′3*x^2-3*a^2=0′，′x′）
或
      ＞＞x=solve（′3*x^2-3*a^2′，′x′）
      ＞＞f=′3*x^2-3*a^2=0′；x=solve（f，′x′）
      ＞＞f=′3*x^2-3*a^2′；x=solve（f，′x′）

結果顯示：

      x=
          -a
          a

【例1-31】 求方程組的解。

解 MATLAB命令如下：

      ＞＞[x，y] =solve（′3*x+y=a′，′x-y=a′，′x′，′y′）
或
      ＞＞f=′3*x+y=a′；g=′x-y=a′；[x，y] =solve（f，g，′x′，′y′）

結果顯示：

      x=
        1/2*a
      y=
        -1/2*a

【例1-32】 求非線性方程組的解。

解 可利用以下命令：

      ＞＞f=′sin（x）+y^2+ln（z）-7=0′；g=′3*x+2^y-z^3+1=0′；h=′x+y+z-5=0′；
      ＞＞[x1，y1，z1] =solve（f，g，h，′x′，′y′，′z′）

結果顯示：

      x1=
          0.5991
      y1=
          2.3959
      z1=
          2.0050

2.符號微分方程求解

在MATLAB中，可利用dsolve（）函數對符號微分方程進行求解，其調用格式為

                [y1，y2，…] =dsolve（′eq1′，′eq2′，…，′cond1′，′cond2′，…，′x′）

或

                [y1，y2，…] =dsolve（′eq1，eq2，…′，′cond1，cond2，…′，′x′）
                [y1，y2，…] =dsolve（exp1，exp2，…，′cond1，cond2，…′，′x′）

其中，eq1，eq2，…表示所求符號微分方程，或不含“等號”的符號微分表達式（此時函數是對eq1=0，eq2=0，…求解的）；exp1，exp2，…僅表示符號微分表達式；cond1，cond2，…表示初始條件或邊界條件；x表示獨立變量，當x省略時，表示獨立變量為t；y1，y2，…表示輸出量。微分方程的記述規定：當y是因變量時，用Dny表示y對x的n階導數。例如，Dny表示形如的導數。

【例1-33】 求微分方程，在邊界條件x（0）=0，x（1）=1的解。

解 MATLAB命令如下

      ＞＞f=′-2*D2x+12*a*t=0′；x=dsolve（f，′x（0）=0，x（1）=1′，′t′）

運行結果：

      x=
        a*t^3+（-a+1）*t

【例1-34】 求一階非線性微分方程的解。

解 MATLAB命令如下

      ＞＞x=dsolve（′Dx=x*（1-x^2）′）

運行結果：

      x=
        [  1/（1+exp（-2*t）*C1）^（1/2）]
        [ -1/（1+exp（-2*t）*C1）^（1/2）]

即該微分方程的解為

【例1-35】 求微分方程組在邊界條件x1（0）=0，x1（π/2）=1，x2（0）=0，x2（π/2）=-1的解。

解 MATLAB命令如下

      ＞＞f=′D2x1-x2=0，D2x2-x1=0′；
      ＞＞[x1，x2] =dsolve（f，′x1（0）=0，x1（pi/2）=1，x2（0）=0，x2（pi/2）=-1′，′t′）

運行結果：

      x1=
            sin（t）
      x2=
          -sin（t）