Mathematica——常微分方程、拉氏變換與級數(shù)實驗

1、創(chuàng)3.5常微分方程、拉氏變換與級數(shù)實驗學習目標1. 會用Mathematica求解微分方程（組）；2. 能用Mathematica求微分方程（組）的數(shù)值解；3. 會利用Mathematica進行拉氏變換與逆變換；4. 能進行幕級數(shù)和傅里葉級數(shù)的展開。一、常微分方程（組）Mathematica能求常微分方程（組）的準確解，能求解的類型大致覆蓋了人工求解的范圍， 功能很強。但不如人靈活（例如在隱函數(shù)和隱方程的處理方面），輸出的結果與教材上的答 案可能在形式上不同。另外，Mathematica求數(shù)值解也很方便，且有利于作出解的圖形。在本 節(jié)中，使用Laplace變換解常微分方程（組）的例子也是十分成

2、功的，過去敬而遠之的方法 如今可以輕而易舉的實現(xiàn)了。求準確解的函數(shù)調(diào)用格式如下：DSolveeqn ,yx ,x求方程 eqn 的通解 y（x），其中自變量是X。求滿足初始條件y （xo） = yoDSolveeqn，yxo= =yo，yx，x 的特解y (x)。DSolveeqn1，eqn2，y 1 x，y2x，x求方程組的通解。DSolveequ1，y1x°= =y10，y 1x，y2x，x 求方程組的特解。說明：應當特別注意，方程及各項參數(shù)的表述方式很嚴格，容易出現(xiàn)輸入錯誤。微分方 程的表示法只有通過例題才能說清楚。例1解下列常微分方程（組）：5 y 吩（X 1）2，1y2八（

3、x x3）y，(3)r fy =z=Z=_y創(chuàng)3.5常微分方程、拉氏變換與級數(shù)實驗y = z（4）丿'的通解及滿足初始條件y （0） =0，z （0） =1的特解。工=-y解：In1: =DSolvey ' x= =2yx/ (x+1) + (x+1) A (5/2)，yx，xf rOut1= "yxT 2(1+x)7/2+(1+x)2c13ln2: =DSolvey ' x= = (1+yxF2 ) /(x+xA3 ) yx)，yx，xOut2= yx >In3: =DSolvey ' x= =zx , z' x= = -yx, yx

4、, zx , xOut3=yx C1Cosx+ C2Sinx,zx C2Cosx- C1Si nxIn4: =DSolvey ' x= =zx , z' x= = -yx , y0= =0 , z0= =1,yx , zx , xOut4=yx Sinx, zx Cosx提示：認真觀察上例，可以從中學習輸入格式，未知函數(shù)總帶有自變量，等號用連續(xù)鍵 入兩個等號表示，這兩點由于不習慣會出錯！導數(shù)符號用鍵盤上的撇號，連續(xù)兩撇表示二階 導數(shù)，這與習慣相同。自變量、未知量、初始值的表示法與普通變量相同。表示。說明：輸出結果總是盡量用顯式解表出，有時反而會使表達式變得復雜，這與教科書的 習

5、慣不同。當求顯式解遇到問題時，會給出提示。通解中的任意常數(shù)用C1，C2，例2求解下列微分方程：(1) y 3y 3y y=(x-5)e，(2) x2(y)2 =1，(3) . y 二 xy。解：In1: =DSolve y x +3y x +3y ' x + yx = = (x - 5) Exp-x，yx，xOut1=/ 2、Lx5x23、x_5x + + e x< 2丿123丿yx > ex22-34、-區(qū) xe»C1 exC2 ex2C3 .3 4In2: =Simplify%1Out2= yx e(-20x3 x4 24C1 24xC2 24x2C3) 24

6、In3: =DSolvexA2 + y ' xA2 = = 1，yx，xOut3= yx-Ixd -x ArcSinx C1，2 2泊知1-八。1ln4: =DSolveSqrty ' x = = x yx ，yx，x3Out4= yx 3 -x -C1說明：由以上可以看出對方程的類型并無限制，但是輸出的答案未必符合習慣，例如第 一個方程的答案需要化簡，有時即使化簡后也未必與教材上的答案一致。例3求微分方程 魚'2xy =xe"的通解。dx解：In1: =DSolvey x+2x yx= = x Ea (幟八2)，yx，x1 2 2 2°ut1=yx

7、-尹 x 宀這就是所給微分方程的通解。式中的C1是通解中的任意常數(shù)。上述命令也可以輸入為：DSolveDyx + 2x yx= =x EA（ - xA2）, yx , x。例4求微分方程xy + y - ex = 0在初始條件y|x=i = 2e下的特解。解：In1: =DSolvex*y x+yx-EAx= =0 ，y1= =2E，yx，xe +exOut1= yx-x二、常微分方程（組）的數(shù)值解函數(shù)NDSolve用于求給定初值條件或邊界條件的常微分方程（組）的近似解，其調(diào)用格 式如下：NDSolveeqns, y 1, y2, x , xmin , xmax求常微分方程（組）的近似解。其中

8、微分方程和初值條件的表示法如同DSolve，未知函數(shù)仍有帶自變量和不帶自變量兩種形式，通常使用后一種更方便。初值點xo可以取在區(qū)間xmin , xmax上的任何一點處，得到插值函數(shù)InterpolatingFunctiondomain, table類型的近似解，近似解的定義域 domain般 為domain, table，也有可能縮小。例5求常微分方程y' = x2 + y2，滿足初始條件y （0） = 0的數(shù)值解。解：In1: =s1=NDSolvey ' x= =xA2+yxA2 , y0= =0,y, x , -2, 2Out1=y InterpolatingFuncti

9、on-2. , 2. , < >In2: = y=y / . s11Out2=InterpolatingFunction-2. , 2., < >ln3: =Plotyx , x , -2 , 2 , AspectRatio-Automatic ,PlotRange-1.5 , 1.5圖13-43微分方程的解曲線Out3= -Graphics-上例中包含許多值得學習的實用內(nèi)容，其中第二項參數(shù)使用y而不是yx，比用yx好。如果求解區(qū)間改為x , -3 , 3,就會出現(xiàn)警告提示,實際得不到-3 , 3上的解。Out1表明返 回的解放在一個表中，不便使用，實際的解就是插值函數(shù)

10、：InterpolatingFunction-2. , 2., < >ln2的結果是用y表示解函數(shù)的名字，因此In3順利畫出解曲線如圖13-43所示例6 求常微分方程組：13x y x x3y " = x滿足初始條件x（0）=0，y（0）=1的數(shù)值解。解：In1: =$仁NDSolvex ' t= = yt - （xtA3/3 - xt）, y' t= = - xt , x0= =0 , y0= =1, x , y, t, -15, 15Out1=x f InterpolatingFunction-15. , 15., < >, y Inter

11、polatingFunction-15. , 15. , < >In2: = x=x / . s11 , 1y=y / . s11 , 2Out2=InterpolatingFunction-15. , 15., < >Out3=InterpolatingFunction-15. , 15., < > ln4: =ParametricPlotxt , yt , t, -15 , 15, AspectRatio Automatic68Out3= -Graphics-說明：上例是求一個著名方程組的近似解，其中 In2也可以改用一個賦值式x , y=x , y /

12、 . Flattens1, 一次得到兩個函數(shù)。通過求數(shù)值解容易得到它的相圖,In4繪制了解的相軌 線如圖13-44所示，圖中表明原點是奇點，極限環(huán)的形狀也已經(jīng)得到。為了應付復雜的情況，需要設置可選參數(shù):Worki ngPrecisio n參見數(shù)值積分部分的介紹AccuracyGoal計算結果的絕對誤差Precisi on Goal計算結果的相對誤差#最大步數(shù)。MaxSteps#Starti ngStepSize初始步長。以上可選參數(shù)的默認值都為 Automatic，其中AccuracyGoal和PrecisionGoal的默認值比 WorkingPrecision小10，當解趨于0時應將Acc

13、uracyGoal取成Infinity。對于常微分方程，最 大步長默認值為1000。這個函數(shù)也可以解偏微分方程，最大步長默認值為 200。例7解下列微分方程（組）：1（1）y' i，滿足初始條件y （0）=1的特解；4yX x = _3x 3y（2）収=xz+26.5x-y，滿足初始條件 x（0）=z（0）=0, y（0）=1的特解。z = xy z解：In1: =NDSolvey ' x= =l/4yx ，y0= =1，y，x， 1，AccuracyGoal20， PrecisionGoa220， WorkingPrecision25Out1=y In terpolat in

14、gFun cti on0，1.000000000000000000000000000, < >In2: =y1 / . %Out2=0.968912424710644784118519 +0.24740395925452292962341091ln3: =NDSolvex ' t= = -3 （xt -yt），y' t = = -xt zt+36.5xt -yt，z' t = = xt yt- zt，x0 = = z0 = = 0，y0= =1，x， y，z，t，0，20， MaxSteps3000Out3=x InterpolatingFunction0.

15、 ，20.，< >，y InterpolatingFunction0.，20.， < >，zInterpolatingFunction0.，20.， < >，ln4: =ParametricPlot3DEvaluatext，yt，zt / . %，t，0，20，PlotPoints 1000圖13-45 3維相軌線Out3= -Graphics3D-說明：以上范例中l(wèi)n1取高精度，而且是復系數(shù)方程。ln2是求解在x=1時的近似值,1i1i求精確解能得到準確值e4，讀者可以求e4的近似值與Out2的結果比較，驗證近似解的精確 度確實很高。ln3在求解時增大步數(shù)

16、，成功地得到了由In 4繪制的如圖13-45所示的解的相 軌線。ln4所示的繪圖語句與前面例子中的不同，現(xiàn)在只要會模仿使用它們就行了，要想弄 清原理請參閱相關Mathematica書籍。三、拉氏變換Mathematica可以進行拉普拉斯變換，其變換使用的函數(shù)調(diào)用格式如下：Lap laceTra nsformf，t， s求函數(shù)f (t)的Lap lace變換，返回自變量為 s的函數(shù)。InverseLaplaceTransformF, s，t求函數(shù) F(s)的 Laplace逆變換，返回自變量為t的函數(shù)。其中函數(shù)f (t)和F (s)也可以是函數(shù)表，這樣可一次變換多個函數(shù) 例8求函數(shù)t 4和et

17、sint的拉氏變換。解：In1: =LaplaceTransformL4，t，s24Out1=飛sln2: =LaplaceTransformExpt Sint，t，sOut2=12 - 2s s2ln3: =lnverseLaplaceTransform%1，s，t4Out3=t4ln4: =lnverseLaplaceTransform%2, s，tOut4= 1 iie(1_i)t(_1 - e2it)2ln5: =FullSimplify%Out 5=et Si nt例9 求函數(shù)f (t)= t3 eat的拉氏變換解：ln1: = LaplaceTransformQ3 Expa t，t

18、，sOut1=6(-as)4以上只是直接進行拉氏變換和逆變換的例子。以下使用拉氏變換解常微分方程，解法原 理見本書理論篇，這里完全實現(xiàn)了計算機求解。例10用拉氏變換解微分方程：x + 3x + 3x' + x = 1 滿足條件 x (0) = x ' (0) = x(0) = 0 的解。解:ln1: =f仁LaplaceTransform xt+3x t+ 3x ' t+xt,t, sOut1=LaplaceTransformxt , t, s +3s LaplaceTransformxt, t, s +23 (sLaplaceTransformxt, t, s - x

19、0) - s x0 +2 f3 (s LaplaceTransformxt, t, s - s x0 - x ' 0) s x' 0- x 0ln2: =s仁LaplaceTransform1, t, s1 sln3:=x 0= x ' 0= x0=0 ;Solvef1= =s1, LaplaceTransformxt, t, sOut2=Out4=LaplaceTra nsform xt, t, s >13s(1 s)In5:1=I nv erseLaplaceTra nsform-s(1 + s)s, t1 _L2Out5= 1 ex(2 2t t2)說明：上

20、例中的LaplaceTransformxt, t, s就是教材中的X (s), In3解出X (s), 其余過程與教科書完全相同?，F(xiàn)在可以將一切計算留給計算機，學生只要弄清解法原理及過 程。技巧：充分利用復制、粘貼功能，可以加快輸入速度，避免鍵入錯誤。上例中 ln5就可 以從Out4中將表達式復制過來。例11求微分方程組：X”-2x"-y+2y = 0i x" + y "_2x = _2e丄滿足條件 x(0)=3, x'( 0) =2, y (0) =0 的特解。解: ln1: =f仁LaplaceTransformx t - 2x' t - y&

21、#39; t + 2yt , x' t + y' t - 2xt,t, s;In2: =s1=LaplaceTransform0, - 2Exp-t , t, s;ln3: = x0=3 ; x' 0= 2 ; y0=0 ;Solvef1= =s1 , LaplaceTransformxt , t, s,LaplaceTransformyt, t, s;ln5: =lnverseLaplaceTransformFlatte nLaplaceTransformxt , t, s,LaplaceTransformyt, t, s / . % , s, tOut5=5 - e

22、-t - 3et + 2e2t, e-t (- 1 + et) 2 (1 + 2 et) In6: =Simplify%-t t2t -t t 2tOut6= 5 - e - 3e + 2e , e - 3e + 2e 說明：在上例中，不顯示任何中間結果，語句比較簡練。其中， ln1和ln2分別對方程 組的左邊和右邊進行拉氏變換，ln3解出X(s)和丫(s)0 In5比較難懂，可以參看前面 的例題，這里是從Out3中自動將解X (s)和丫(s)提取出來，再進行拉氏逆變換。Out5是x (t), y (t) , Out6將答案化簡。本例已經(jīng)將求解過程一般化，只需改變方程組和初 值的數(shù)據(jù)，就可以解

23、其它方程組了。四、級數(shù)1. 求和與求積求有限或無窮和、積的函數(shù)是：i maxSumf, i , imin , imax求二 f (i)，其中 imin 可以是-x, minimax可以是x(即+x)，但是必須滿足imin < imax。基本輸入模板中也有求和專用的符號， 使用模板輸入更方便。求多重和，也可以使用基本輸入模i max求I丨f (i),基本輸入模板中也有i=imin求多重積，也可以使用基本輸入Sumf, i , imin , imax , j , jmin, jmax， 板連續(xù)多次輸入求和符號得到。Productf , i , imin, imax求積符號Productf ,

24、 i , imin , imax, j , jmin , jmax,模板連續(xù)多次輸入求積符號得到。例12求下列級數(shù)的和與積:n(1) k2 ,(2)k =1二 1J,(3)' 丄,(4)k=1 kO0 1II ek2。k=173#解：In1: =SumkA2 , k , 1 , nOut1= 1 n(1n)(1 2n)Out2=6In2:oO='1你八2k =1712QOln3: = ' 1/ kk-1Sum： div : Sum does not conv erge.:1Out3=' -k 4 kcdln4: =j【Exp1/kA2kJOut4= e莎說明：上

25、例中第三個級數(shù)發(fā)散，Mathematica給出提示，并在不能給出結果時將輸入的式 子作為輸出。NSum和NProduct得到數(shù)值解。2. 將函數(shù)展開為幕級數(shù)將函數(shù)展開為幕級數(shù)的函數(shù)調(diào)用格式如下:Seriesf, x , xo, n為止。Seriesf, x , xo, n, y , yo, m例13展開下列函數(shù)為幕級數(shù)：sin x(1) y=tgx,( 2) y 匸x將函數(shù)f (x)在xo處展成幕級數(shù)直到n次項將函數(shù)f (x, y)先對y后對x展開。 y = f (x),( 4) y = eTyo74#解：ln1:=SeriesTanx,x ,0 , 935x2x17x762x9Out1=x3

26、153152835ox10ln2: =SeriesSinx /x , x , 0 , 92468亠x x xx10Out2= 1ox61205040362880ln3: =Seriesfx , x , 1 , 71(x-1)31 1OUt3= f1 )2f 1)(1)2 61241720f 1(x-1)6f 1(x-1)75 0 4 0ox -18#ln4: =SeriesExpx y, x , 0 , 3 , y , 0 , 2#/ 2、Out4= 1 (y oy3)x - iy oy3 x2 oy3x3 ox4l2丿說明：上例中l(wèi)n3表明也可以展開抽象的函數(shù)。對已經(jīng)展開的幕級數(shù)進行操作的兩

27、個函數(shù)是：Normalexpr將幕級數(shù)expr去掉余項轉換成多項式SeriesCoefficientexpr, n找出幕級數(shù) expr 的 n 次項系數(shù)。例14將y = arcsinx展開為幕級數(shù)，只取前9項并去掉余項。10ox解：In1: =SeriesArcSinx，x，0，9In2:=Normal%Out2=x33x55x7 35x9=x亠亠亠640112 1152Out1=3579L x ± 3x ± 5x i 35x x6401121152In3: =SeriesCoefficient%1，53Out3= 403. 傅里葉級數(shù)求傅里葉級數(shù)就是求出傅里葉系數(shù)，傅里葉系數(shù)是一個積分表達式，所以利用積分函數(shù)In tegrate就可以實現(xiàn)。例如，設周期矩形脈沖信號的脈沖寬度為t,脈沖幅度為E,周期