一維振子運動的數值求解研究

張定梅

摘要：建立了一維振子運動的牛頓運動學方程，利用MATLAB軟件數值研究了振子的動力學特性.通過調節參數實現了簡諧振動、阻尼振動、受迫振動及復雜非線性振動.給出了各種振動的運動時間序列，討論了參數對振子的振動周期及振幅的影響，并使用功率譜和吸引子圖分析了振子產生的混沌的特性.與傳統的解析求解方法相比，數值求解方法更加快捷、直觀，避免了復雜的數學積分運算，并且通過運動軌跡圖可快速掌握振子的運動規律.

關鍵詞：一維振子;MATLAB;Runge-Kutta算法;混沌

中圖分類號：O322? 文獻標識碼：A? 文章編號：1673-260X（2019）09-0019-03

1 引言

機械振動是自然界中普遍存在的運動形式，也是大學物理課程中一項重要的教學內容.對于簡單的簡諧振動可以利用數學積分的方法給出解析解[1]，但對于存在復雜外部擾動的非線性振動系統，解析求解將變得非常困難，因此數值求解將是一種好的替代方式.關于振動問題數值求解的研究已見有關報道[2-4]，例如，文獻[2]利用拉格朗日方法建立了復雜的六彈性振子的二維運動方程，并用maple軟件進行了數值求解，給出振子運動時間序列圖.文獻[3]研究了各種組合彈簧振子中的非線性振動問題，給出了數值解和解析解并進行了對比.文獻[4]研究了Duffing振子和Van der Pol振子耦合之后的非線性的動力學行為，利用Simulink仿真給出了不同頻率及耦合系數下耦合振子的相圖和龐加萊截面圖.在大學物理教學中，一維振子運動是一種最基本的教學內容，為了能夠系統的、直觀的研究一維振子在各種外部擾動下的動力學特性，本文建立了一維振子的牛頓運動學方程，采用四階Runge-Kutta算法并利用MATLAB軟件進行數值求解，通過合理調節參數，實現了簡諧振動、阻尼振動、受迫振動及復雜非線性振動，給出了振子的運動時間序列圖，并且得到了有趣的混沌現象.

2 四階Runge-Kutta算法

四階Runge-Kutta算法是求解常微分方程及方程組常用的數值解法，它基于Taylor展開，要求被求解的函數具有較好的光滑性，并且其截斷誤差為o（h5），具有較高的求解精度.四階Runge-Kutta公式的一般表達形式為：

其中：

3 一維振子運動學方程

在實際生活中，一維振子除了受到回復力作用以外，往往還受到阻力和外部擾動，為了全面的描述一維振子的運動情況，本文建立的一維振子運動學方程如下[5-7]：

其中，m為振子的質量，-k1x為線性回復力， 為與速度成正比的阻力，-k3x3為與位移有關的三階非線性阻力項，F1cos1t+F2sin2t為周期性的外力.為簡化運算，m取1kg，并且將二階常微分方程化為一階微分方程組求解，具體表達式如下：

4 仿真結果及結果

本文采用四階Runge-Kutta算法對方程（4）進行數值求解，求解步長為0.01s，采用MATLAB軟件進行仿真.

4.1 簡諧振動

簡諧振動是振子最簡單的運動形式，此時的振子只受到回復力的作用，方程（4）的參數取值情況如圖題所示.圖1給出了振子在不同的回復系數下的振動曲線，如圖1（a）所示，當A=1s-2時，此時振子的運動軌跡為一規則正弦曲線，振幅為0.110m，振動頻率為0.16Hz.當A=3s-2時（圖1（b）），由于回復力的增大，振幅變為0.063m，振動頻率變為0.28Hz.當A=6s-2時（圖1（c）），振幅減小為0.045m，振動頻率增大為0.39Hz.從圖中可以直觀地看出隨著回復系數的增大，振子的振動振幅慢慢縮小，同時振動頻率慢慢變大.

4.2 阻尼振動

阻尼振動是指振子受到外界阻力下的運動，本文中考慮阻力與振子的運動速度成正比的情況.在阻尼振蕩下，方程（4）的參數取值情況如圖題所示.圖2（a）給出了阻尼系數B=0.01s-1時振子振動的時間序列，可以看出，由于此時阻尼系數較小，振子的振幅衰減較慢.當繼續增加阻尼系數到0.03s-1，如圖2（b）所示，振子振動到200s時振幅幾乎為零.繼續增大阻尼系數到0.06s-1（圖2（c）），振子到150s就幾乎停止振動了.因此，從圖中可以直觀地看出，隨著阻尼系數的增大，振子的振幅衰減的越來越快.

4.3 受迫振動

受迫振動指的是振子在外部周期力作用下的運動，振子在不受外力作用時，由于自身恢復力的作用會產生一個固有振蕩頻率，此頻率與回復系數有關，本文中回復系數A取2s-2，相應的本征角頻率為1.414Hz.圖3考慮了三種不同的外部擾動頻率下振子的運動情況.圖3（a）為D=0.5ms-2，ω1=1Hz情況下的受迫振動振動曲線，此時外部擾動頻率與振子的固有頻率不相同，但差距較小，振子在72s后趨于規則振蕩，振蕩幅度為0.52m.當D=0.5ms-2，ω1=1.414Hz時（圖3（b）），此時外部擾動頻率與振子的固有頻率相同，出現了共振，振子的振蕩幅度在規則振蕩后達到了4.2m.當D=0.5ms-2，ω1=3Hz時（圖3（c）），此時外部擾動頻率與振子的固有頻率差距較大，穩定后的振蕩幅度只有0.07m，未形成共振.從圖中可以清楚地看出在不同的外部擾動頻率下振子振幅及頻率的變化.

4.4 非線性振動

圖4給出了振子在非線性力作用情況下的振動曲線圖及功率譜圖.其中功率譜指的是將時域中的時間序列通過傅里葉變換轉換到頻域進行觀察，它代表了該時間序列中某一頻率間隔內的功率占總功率的比例，能夠直觀地看出時間序列中所包含的頻率成分.如圖4（a）所示，振子的運動軌跡為規則的脈沖包，時間序列（圖4（a））中雖然有多個峰值，但整體的振動是重復的，功率譜（圖4（a1））中可以看出多個頻率成分，其中最強的峰值為振動的主頻率（約為0.2GHz），主頻率前面的峰值為次諧波.圖4（b）所示，此時振子的運動出現混沌行為，時間序列（圖4（b））無明顯周期行為，并且出現了隨機的振動，功率譜（圖4（b1））出現了展寬并無明顯的峰值.

混沌是一種有趣的物理的現象，它具有自相似性、初值敏感性和內在隨機性等特性.為了進一步研究所產生的混沌的性質，圖5（a）給出了不同初始條件下的振子的運動軌跡，從圖中可以看出，即使所有參數都一樣，當初始條件不同，兩條軌跡在27s的時候出現了偏離，而且隨著時間的增加，偏離越來越大，最后軌跡完全不相同，這就證明了混沌具有初值敏感性.圖5（b）給出了該參數情況下的時間序列的吸引子圖，圖中有多個環路纏繞在一起，并且具有一定的自相似性，表現為奇怪吸引子.

5 總結

本文采用四階Runge-Kutta算法數值研究了一維振子的運動問題，利用MATLAB軟件進行了仿真模擬，通過調節參數實現了簡諧振動、阻尼振動、受迫振動及復雜非線性振動，給出了各種振動的運動時間序列圖.從圖中可以直觀地看出，對于簡諧振動，回復系數越大，振子運動的振幅越小.對于阻尼振動，隨著阻尼系數的增大，振子的振幅衰減的越來越快.

對于受迫振動，當外部擾動的頻率與振子本征頻率相同時，振子的振幅將大幅度增加.最后對于非線性振動情況，適當調節參數可實現混沌振動.本文簡單直觀地展示了振子在各種情況下的運動規律，可作為大學物理教學內容的有效補充.

——————————

參考文獻：

〔1〕王樹平.對稱非線性彈簧振子的振動周期的級數解[J].大學物理，2017，36（4）：15-16.

〔2〕包興明，袁玉全，閆安英，戚作濤對稱六彈性振子的二維非線性振動[J].大學物理，2010，35（2）：1-6.

〔3〕廖旭，任學藻.組合線性彈簧振子中的非線性振動[J].大學物理，2008，27（2）：25-28.

〔4〕王曉東，楊紹普，趙志宏.Duffing振子和Van der Pol振子耦合的動力學行為分析[J].石家莊鐵道大學學報（自然科學版），2015，28（4）：53-57.

〔5〕卓崇培.非線性物理學[M].天津：天津科學技術出版社，1996.26-61.

〔6〕劉式適，劉式達.物理學中的非線性方程[M].北京大學出版社，2000.6-66.

〔7〕倪致祥，馬濤.一種非簡諧的微振動模型[J].大學物理，2006，25（9）：14-16.