廖川，姜麗娜，賈紅寶

（遼寧科技大學(xué)理學(xué)院，遼寧鞍山114051）

受迫倒擺系統(tǒng)動力學(xué)解的數(shù)值模擬

廖川，姜麗娜，賈紅寶

（遼寧科技大學(xué)理學(xué)院，遼寧鞍山114051）

本文分析了達(dá)芬方程的動力學(xué)特征。采用達(dá)芬方程作為受迫倒擺系統(tǒng)的模型，運用MATLAB程序模擬了倒擺系統(tǒng)的動力學(xué)解。通過調(diào)節(jié)參數(shù)并采用位移曲線、相圖、功率譜和龐加萊截面圖顯示了受迫振動倒擺系統(tǒng)的周期解和混沌解。本文對于大學(xué)物理的教學(xué)改革和控制工程應(yīng)用有一定啟發(fā)意義。

彈性系統(tǒng)；達(dá)芬方程；微擾法；倒擺

倒擺系統(tǒng)是在機(jī)械控制、航空航天器控制、精密儀器等方面擁有廣闊應(yīng)用前景的非線性系統(tǒng)［1］。在對它的研究中，多數(shù)以直線型小車驅(qū)動倒擺為基本模型［2］，目前主要的研究手段有Melnikov方法、Runge-Kutta方法等。非線性系統(tǒng)的受迫振動在理論和工程技術(shù)中都具有重要的地位［3］，計算機(jī)模擬技術(shù)使得非線性系統(tǒng)得到更加深入的研究，自上世紀(jì)70年代以來非線性問題已成為科學(xué)研究的前沿之一［4］。由于達(dá)芬方程能反映很多非線性現(xiàn)象的本質(zhì)，并在振動控制、高階諧波平衡［5］、機(jī)械動力學(xué)、控制系統(tǒng)的自激振動［6］、結(jié)構(gòu)動力學(xué)等領(lǐng)域得到了廣泛應(yīng)用［7］，而微擾法［8］則是求解達(dá)芬方程的常用方法。所以本文運用達(dá)芬方程的微擾解法得到有阻尼和無阻尼兩種情況下非線性諧振子的共振原理以及非線性受迫振動解的特性，通過MATLAB編程對受外力驅(qū)動的倒擺系統(tǒng)的振動解進(jìn)行模擬，繪制出倒擺系統(tǒng)的位移曲線、相圖、功率譜、龐加萊截面圖并進(jìn)行分析。期望對于大學(xué)物理教學(xué)改革以及控制工程的應(yīng)用有所裨益。

1 無阻尼達(dá)芬方程的諧振解

達(dá)芬方程的一般形式為

式中：α,κ,μ為常數(shù)。

采用微擾法研究無阻尼情況下的達(dá)芬方程。引入無限小參量ε（0＜ε＜＜1）并把非線性和強(qiáng)迫力都看成是對系統(tǒng)的微擾，無阻尼（α=0）時的達(dá)芬方程為

與式（1）比較可得

式中：ω可認(rèn)為是無微擾時系統(tǒng)的頻率；ω0是無微擾、無阻尼時的頻率。

作為周期解，式（4）應(yīng)滿足周期條件x(τ+2π) =x(τ)，既然φ待定，這就可以選擇較方便的初始條件

現(xiàn)在尋找方程式（4）展開成小量ε的級數(shù)解

將式（5）和式（6）代入方程（4），令左右兩邊ε的同次冪項相等，得

利用以上條件可逐級求解xi（i=0，1，2，…），結(jié)果為

式中：A0,A1,A2為常數(shù)，如要確定常數(shù)的大小就需要下一級近似解，例如：若需確定A0，則必須確定第一級近似解x1(τ)；若需確定A1，則須確定第二級近似解x2(τ)。如此下去，可得更高級近似。

2 有阻尼達(dá)芬方程的諧振解

與式（1）比較系數(shù)，可得式（3）和εζω=α。采用處理無阻尼時的相同微擾方法，可以得到x1為周期函數(shù)所應(yīng)滿足的條件

由式（11）可得有阻尼時的零級近似位相，式（10）和式（11）可改寫為

式（12）與式（13）相除，并利用εζω=α，可得

將式（3）代入式（14）得

由于滿足周期條件，有Ω=ω，于是

式（15）表明，有阻尼時振動與外力并不同步［4］。由式（10）和式（11）還可得到零級近似的振幅與頻率的關(guān)系，過程如下：

由式（12）和式（13）兩式的完全平方和得

利用式（3）并將a,b,ξ代入（16）得

由于滿足周期條件，有Ω=ω，可得

當(dāng)μ＜0時，系統(tǒng)的等價自振頻率將隨振幅的增加而減少。當(dāng)μ＞0時，系統(tǒng)的等價自振頻率將隨振幅的增加而增大。

可求出［9］

可見式（20）與式（18）最主要的差別在于：受外力驅(qū)動的線性振動的自振頻率是常數(shù)ω0，而受外力驅(qū)動的達(dá)芬方程（非線性振動）的等價自振頻率ωe是隨振幅的變化而變化的。

3 外力驅(qū)動的倒擺系統(tǒng)

倒擺系統(tǒng)是非線性和強(qiáng)耦合的不穩(wěn)定動力學(xué)系統(tǒng)，可以采用達(dá)芬方程作為模型?？刂频箶[的過程涉及到非線性、魯棒性、以及跟蹤等問題［10］，這些都是控制領(lǐng)域的核心問題。

本文主要通過數(shù)值模擬實驗展示在適當(dāng)?shù)膮?shù)下，倒擺的強(qiáng)迫振動將出現(xiàn)混沌現(xiàn)象，并討論這種混沌現(xiàn)象。

首先，對外力驅(qū)動的倒擺運動微分方程進(jìn)行部分簡化，再進(jìn)行無量綱化可得到外力驅(qū)動的倒擺的達(dá)芬方程［11］

式中：δ=β/mΩ0，ω=Ω/Ω0，f=(A/θ0)(Ω/Ω0)2；A, Ω,Ω0,θ0為倒擺系統(tǒng)的已知常數(shù)，θ可近似表示桿對鉛垂線的偏離，β為阻尼系數(shù)。

設(shè)x=y1，dx/dt=y2，則式（21）可化為兩個一階微分方程

選擇適當(dāng)?shù)某跏贾担肕atlab繪制位移曲線如圖1所示。

由圖1可觀察出，在初始階段，由初值微擾產(chǎn)生的兩曲線振蕩的差別并不大，但是經(jīng)過較長一段時間以后兩者就有明顯差別?；蛘哒f，在初始階段內(nèi)，倒擺的振蕩的行為是可預(yù)言的，但是在長時間振蕩以后其行為就無法預(yù)言。倒擺系統(tǒng)的方程是確定性的，但是方程解的演化卻對初值敏感，這種現(xiàn)象通常稱為混沌現(xiàn)象。

圖1 混沌狀態(tài)下，初值有微擾形成的兩條曲線Fig.1 Initial perturbation formation of two curves in chaos state

下面來研究受迫倒擺系統(tǒng)周期解的形式。通過對系統(tǒng)參數(shù)的調(diào)節(jié)，采用位移曲線、相圖、功率譜和龐加萊截面圖來分析系統(tǒng)在不同參數(shù)下周期解和混沌解的特征，可以看到存在不同的吸引子，即周期l吸引子、周期2吸引子和奇怪吸引子。

取阻尼項參數(shù)δ=1.5，此時外加驅(qū)動周期力很?。磃很小），外驅(qū)力的振蕩對非線性系統(tǒng)的作用很弱，于是倒擺系統(tǒng)的運動可以近似看做為兩個獨立線性振蕩的疊加，振蕩頻率為ω=Ω，產(chǎn)生的單周期解如圖2所示。其中位移曲線只能觀察到一個振動頻率，相圖軌道只顯示一個中心，功率譜圖只存在一個頻率，對應(yīng)的龐加萊截面圖中也只有一個點，這些結(jié)論說明此時是單周期解或者叫周期解1吸引子，倒擺的周期與外加驅(qū)動力的周期相同，因此表示為主諧波運動。

取阻尼項參數(shù)δ=1.28，加大外加驅(qū)動力的振幅（即增大f），外加驅(qū)動力的周期作用使倒擺系統(tǒng)的振蕩出現(xiàn)分頻，倒擺振蕩周期τ被鎖在外驅(qū)動力周期的有理數(shù)倍數(shù)上τ=nT（n為有理數(shù)）。在這個參數(shù)δ下，倒擺出現(xiàn)雙周期解，對應(yīng)于有理倍數(shù)n=2（如圖3）：位移圖中存在兩個振動頻率，相圖軌道的吸引子有兩個周期，同樣功率譜圖顯示兩個頻率，而龐加萊截面圖中也出現(xiàn)兩個點，即此時倒擺周期解的周期等于外加驅(qū)動力周期的2倍，可稱為周期解2吸引子，代表次諧波運動。

圖2 周期解1Fig.2 Periodic solution 1

圖3 周期解2Fig.3 Periodic solution 2

取阻尼項參數(shù)δ=0.86，一旦增加f超過某一臨界值時，倒擺自身的振動和外驅(qū)力的耦合作用可以強(qiáng)大到使得倒擺振蕩的振幅超出原來流域的邊界，于是倒擺系統(tǒng)的解從原先的流域進(jìn)入到一個新的流域。同樣，這時倒擺的解也有可能大到超過新流域的邊界而使其又進(jìn)入另一新的流域，也有可能又回到原來的流域。倒擺的解來回在幾個不同的流域之間跳動，其相軌道既有局域不穩(wěn)定又有全局穩(wěn)定性的特點，由于漲落而帶來的這種隨機(jī)性運動稱之為混沌［12］。混沌解的特征（如圖4）：位移曲線沒有明確的周期，功率譜圖呈現(xiàn)連續(xù)譜，而相圖軌跡和龐加萊截面圖顯示隨機(jī)性吸引區(qū)域，稱之為奇怪吸引子。

圖4 混沌解Fig.4 Chaotic solution

4 結(jié)論

本文運用求解達(dá)芬方程的微擾解法。首先得到無阻尼和有阻尼兩種情況下非線性諧振子的共振原理及非線性受迫振動的一些動力學(xué)特性，然后給出倒擺系統(tǒng)的達(dá)芬方程。借助Matlab工具，通過適當(dāng)選取方程中阻尼項參數(shù)δ，并且逐漸改變外加周期性驅(qū)動力的大小分別得到了位移曲線、相圖、功率譜和龐加萊截面圖。通過分析上述四種模擬結(jié)果的圖形可以分析出倒擺的周期解1吸引子、周期解2吸引子以及混沌解對應(yīng)的奇怪吸引子，這些圖形從不同側(cè)面表征出倒擺系統(tǒng)的運動形態(tài)，非常直觀清晰。

［1］劉珊珊，周良強(qiáng)，陳芳啟.一類倒擺系統(tǒng)的混沌運動［J］.中國科技論文，2016，11（5）：494-496.

［2］MRAD F，EL-HASSAN N，MAHMOUD E H，et al.Real-time control of free-standing cart-mounted inverted pendulum using LabVIEW RT［C］//Industry Applications Conference，Conference Record of the IEEE Xplore，2000，2：1291-1298.

［3］項國波.非線性系統(tǒng)［M］.北京：知識出版社，1991：25-42.

［］［］

［4］胡靜，彭芳麟，管靖，等.理論力學(xué)中非線性問題的MATLAB數(shù)值解［J］.大學(xué)物理，2001，20（10）：39-41.

［5］劉南，白俊強(qiáng)，華俊，等.高階諧波平衡方法中非物理解來源分析及改進(jìn)方法研究［J］.力學(xué)學(xué)報，2016，48（4）：897-906.

［6］傅利峰.基于MEMS技術(shù)的非線性靜電式能量采集器的研究［D］.杭州：浙江大學(xué)機(jī)械工程學(xué)院，2016：50-70.

［7］牛江川，申永軍，楊紹普，等.基于速度反饋分?jǐn)?shù)階PID控制的達(dá)芬振子的主共振［J］.力學(xué)學(xué)報，2016，48（2）：422-429.

［8］陳予恕.非線性振動［M］.北京：高等教育出版社，2002： 259-267.

［9］湯艷芬.非線性系統(tǒng)受迫振動的周期解［J］.大學(xué)物理，2004，23（9）：25-28.

［10］張麗娟，涂亞慶，王保中等.小車-倒立擺系統(tǒng)數(shù)學(xué)模型［J］.中南大學(xué)學(xué)報（自然科學(xué)版），2005，36（8）：199-202.

［11］彭芳麟.理論力學(xué)計算機(jī)模擬［M］.北京：清華大學(xué)出版社，2002：198-203.

［12］孟艷麗，王素秋，韓晶.基于Matlab的非線性動力學(xué)系統(tǒng)分析［J］.物理實驗，2005，25（8）：42-45.

Numerical simulation of dynamic solution of inverted pendulum system

LIAO Chuan，JIANG Lina，JIA Hongbao
（School of Science，University of Science and Technology Liaoning，Anshan 114051）

The dynamic characteristics of Duffing equation is analyzed in this paper.Duffing equation is taken as the model of the inverted pendulum with external driven force，and the dynamical solution of the inverted pendulum is simulated with MATLAB program.By adjusting parameter，the periodic solutions and chaotic solutions are shown with displacement curve，orbital phase diagram，power spectrum and Poincare section.It can be taken as an example for educational reform in university physics and the application of control engineering.

elastic system；Duffing equation；perturbation method；inverted pendulum

February 22，2017）

O322

A

1674-1048（2017）02-0143-06

10.13988/j.ustl.2017.02.012

2017-02-22。

遼寧省科技廳項目(2015020230)；校基金（2014TD01）。

廖川（1996—），男，四川營山人。

姜麗娜（1962—），女，胡南長沙人，教授。