剛體運動三維可視化模擬和復(fù)擺相關(guān)問題探討

林炯輝，劉玉穎，宋 敏

(1.中國農(nóng)業(yè)大學(xué) 工學(xué)院，北京 100083；2.中國農(nóng)業(yè)大學(xué) 理學(xué)院，北京 100083)

自17世紀(jì)牛頓經(jīng)典力學(xué)建立以來，就一直得到普遍的應(yīng)用，它打開了近代自然科學(xué)的大門，在其基礎(chǔ)上發(fā)展出了彈性力學(xué)、流體力學(xué)、結(jié)構(gòu)力學(xué)等多門學(xué)科. 在牛頓力學(xué)中，如果已知物體的受力和初始狀態(tài)，對于有序系統(tǒng)其未來的運動狀態(tài)是確定的. 學(xué)者們根據(jù)經(jīng)典力學(xué)的定律和萬有引力定律曾經(jīng)精確地預(yù)言彗星和小行星等的運動，并且得到了驗證.

近年來研究者們在大學(xué)物理教學(xué)中引入Matlab計算機(jī)數(shù)值模擬，培養(yǎng)學(xué)生用計算機(jī)進(jìn)行學(xué)習(xí)和研究物理問題的能力[1-4]，將VPython引進(jìn)大學(xué)物理課堂的方式也逐漸興起[5-7]. VPython可以利用簡單的代碼模擬物體在三維空間中的運動，實時觀測物體的速度變化及其運動軌跡. 質(zhì)點平動、剛體定軸轉(zhuǎn)動和復(fù)擺振動是基本的運動形式，本文選用VPython軟件對剛體的平動、轉(zhuǎn)動、振動進(jìn)行模擬與探究，利用轉(zhuǎn)動定律、角動量定理、動量定理等基本規(guī)律，從剛體的定軸轉(zhuǎn)動出發(fā)，研究其在恒定力矩和變化力矩作用下的運動狀態(tài)以及其能量特征，進(jìn)而研究剛體的平面平行運動；通過模擬不同角度復(fù)擺的運動，觀察復(fù)擺受迫振動時的運動特征，檢驗相關(guān)結(jié)論.

1 Vpython概況

VPython是Python默認(rèn)的3D模塊，有著和Python一樣風(fēng)格的編程方式，可以快速創(chuàng)建三維場景和動畫. Vpython可創(chuàng)建多種實體類型，如球體、圓柱體、長方體、錐體等，通過代碼改變其在空間中的位置即可實現(xiàn)運動模擬. 利用VPython模擬剛體轉(zhuǎn)動的代碼框架為：① 定義常數(shù)(m，g，k等)和運動初始值(初始角位置，初始角動量等)；② 設(shè)置步長Δt，Δt越小模擬越精確；③ 進(jìn)入while循環(huán)，在循環(huán)中根據(jù)當(dāng)前狀態(tài)計算合外力矩；④ 用該力矩更新系統(tǒng)角動量；⑤ 利用角動量更新系統(tǒng)角位置. 程序中以差分代替微分，得到角速度、角位置對時間的迭代式(Euler-Cromer方法)：

ωi+1=ωi+f(θi,ωi,ti)Δt，

θi+1=θi+ωi+1Δt，

ti+1=ti+Δt

Δt越小，模擬的精度越高.

2 剛體轉(zhuǎn)動的VPython模擬

2.1 轉(zhuǎn)動定律簡介

在大學(xué)物理的學(xué)習(xí)中，我們依據(jù)牛頓第二定律導(dǎo)出了剛體定軸轉(zhuǎn)動定律M=Jα，即剛體所受的對于某定軸的合外力矩等于剛體對此定軸的轉(zhuǎn)動慣量與剛體在此合外力矩作用下所獲得的角加速度的乘積；力矩的時間累積效應(yīng)——角動量定理dL=Mdt，即剛體角動量的增量等于所受合外力矩的角沖量.

2.2 恒定力矩下勻質(zhì)圓柱體的定軸轉(zhuǎn)動

有一個質(zhì)量為2 kg，長度為1 m，半徑為0.1 m的勻質(zhì)圓柱體，轉(zhuǎn)軸為過其質(zhì)心且與圓柱體中心軸線垂直. 勻質(zhì)圓柱體最初處于靜止?fàn)顟B(tài)，對其施加一恒定力矩，力矩方向沿轉(zhuǎn)軸正方向，大小為2 Nm,對其運動進(jìn)行分析.

對于勻質(zhì)圓柱體，轉(zhuǎn)軸(圖中所示z軸)位于圓柱體中心且垂直于中心軸線，建立如圖1所示的坐標(biāo)系.

圖1 勻質(zhì)圓柱體，其轉(zhuǎn)軸為z軸

(1)

我們在VPython中對該勻質(zhì)圓柱體的轉(zhuǎn)動進(jìn)行模擬，利用VPython顯示勻質(zhì)圓柱體尺寸、顏色以及轉(zhuǎn)軸位置(如圖2所示)，將棒的運動可視化.

圖2 勻質(zhì)圓柱體繞定軸轉(zhuǎn)動及對應(yīng)坐標(biāo)系，轉(zhuǎn)軸沿z軸正向

在VPython模擬中，可以直觀地觀察它的運動過程，程序中使用的物理公式及對應(yīng)的程序循環(huán)語句如下.

whilet<2:

rate(100)

torque=vector(0,0,2) #M=(0,0,2)

L=L+torque*dt#dL=M×dt

dtheta=omega*dt#dθ=|ω|×dt

theta=theta+dtheta

rod.rotate(angle=dtheta,axis=(0,0,1),origin= axle.pos)

t=t+dt

VPython模擬結(jié)果表明：t=2 s時，角位置θ=23.3044661359 rad，角速度ω=23.3021359223 rad/s，ω/θ=0.999900，與理論計算結(jié)果符合得很好.

本文還探究了勻質(zhì)圓柱體轉(zhuǎn)動過程中角速度、角位置隨時間的變化關(guān)系：由于恒定力矩作用，角加速度為一恒量，角速度隨時間線性增大，角位置為時間的二次函數(shù)，與時間t的平方成正比(如圖3所示).

圖3 角速度、角位置與時間關(guān)系圖像

2.3 正余弦變化的力矩下勻質(zhì)圓柱體的定軸轉(zhuǎn)動

(t=0,ω=0)

(2)

-0.699cos(5t)+0.699,(t=0,θ=0)

(3)

進(jìn)行VPython模擬，將力矩定義式更改為torque=3cos(5t).模擬結(jié)果表明：角速度ω、角位置θ隨時間t均做周期性變化，角速度ω、角位置θ隨時間的變化曲線如圖4所示. 理論與模擬結(jié)果均表明角速度、角位置的變化周期均等于力矩變化周期；且在運動的一個周期內(nèi)，圓柱體先加速后減速(圖5a)，當(dāng)角速度減為0時，轉(zhuǎn)角達(dá)到最大(圖5b).隨后勻質(zhì)圓柱體反向加速(圖5c)到最大值再減速，減速至角速度為零時，此時θ=0，棒回到初始位置(圖5d).

圖4 角速度、角位置隨時間變化關(guān)系曲線

圖5 勻質(zhì)圓柱體轉(zhuǎn)動的3D場景

力矩M和轉(zhuǎn)動動能E隨時間變化曲線如圖6所示，勻質(zhì)圓柱體動能變化周期等于力矩周期的一半.動能為

力矩大時，轉(zhuǎn)動動能不一定大.通過動能與力矩的函數(shù)關(guān)系圖像(圖7)可知：當(dāng)力矩為極值±3Nm時，動能最小為0，當(dāng)力矩為0時，動能為最大值1.048 J.

圖6 動能、力矩隨時間變化關(guān)系曲線

圖7 動能與力矩關(guān)系圖像

2.4 平動加轉(zhuǎn)動

現(xiàn)在在轉(zhuǎn)動軸上施加一恒力，勻質(zhì)圓柱體的運動變?yōu)橘|(zhì)心的平動和繞過質(zhì)心軸的轉(zhuǎn)動.由動量定理P2-P1=FΔt，得到動量迭代式P=P+Fdt，由P=mv，求得任意時刻的瞬時平動速度v和dt內(nèi)的位移ds=vdt.在while循環(huán)中更新轉(zhuǎn)軸的位置(設(shè)力沿x軸正方向，故只需改變勻質(zhì)圓柱體質(zhì)心x坐標(biāo)).

whilet<2:

P=P+F*dt

#dPtranslate=F*dt

v=P/M

rod.pos=rod.pos+v*dt

#dx=vdt

axle.pos.x=axle.pos.x+mag(v*dt)

在相同的力矩和大小不同的外力作用下，系統(tǒng)的運動軌跡不同，如圖8所示.力矩均為M=3cos5t，外力大小分別為0.1 N、0.5 N、1 N沿x軸正向(方向向右)，模擬結(jié)果表明：勻質(zhì)圓柱體端點的軌跡呈鋸齒狀，且隨著外力的增大，“鋸齒”越來越疏松(如圖8所示).

圖8 相同的力矩和不同外力作用下勻質(zhì)圓柱體的運動軌跡，力矩均為M=3cos5t

3 勻質(zhì)棒的擺動(大角度、小角度)的VPython模擬

3.1 重力矩下的棒的擺動(大角度、小角度)

在大學(xué)物理振動學(xué)章節(jié)中，我們了解了簡諧振動是最簡單、最基本的振動，小角度復(fù)擺即是簡諧運動的一個特例. 在研究復(fù)擺的運動周期時，令物體作小振幅擺動(擺角<5°). 根據(jù)轉(zhuǎn)動定律得

(3)

式中l(wèi)為復(fù)擺質(zhì)心到轉(zhuǎn)軸的距離，當(dāng)θ較小時，可近似看作sinθ≈θ，式(3)變?yōu)?/p>

解微分方程得

θ=Acos(ωt+ψ)

振動的圓頻率和周期分別為[8]

(4)

(5)

通過積分理論計算，得大擺角運動時復(fù)擺的周期為[10]

(6)

可以看出任意角度下擺的運動方程與運動周期的理論推導(dǎo)式都與小角度近似(sinθ≈θ)存在差異，θ0為大角度時，擺的運動已非簡諧運動.在VPython中，對一根懸掛在低摩擦軸上的桿的運動(忽略阻力)進(jìn)行建模(圖9)，設(shè)桿的質(zhì)量為m，桿長為L，并設(shè)立一組遞增的初始偏置角度.模擬結(jié)果如圖10所示，虛線為小角度近似得到的余弦曲線:

θ=θ0cos(ω0t),

實線為擺桿實際角位移曲線.當(dāng)θ0為小角度時(<10°)，可見擺桿實際角位移圖線與簡諧運動角位移理論計算結(jié)果曲線相重合，擺桿的運動為簡諧運動；當(dāng)θ0>10°時，隨著初始角的增大，實際角位移曲線與簡諧運動角位移理論計算結(jié)果曲線偏差越大，擺桿不再作簡諧運動(如圖10(c)-(f)所示)；而當(dāng)θ0=180°時，擺桿達(dá)到往復(fù)擺動與圓周運動臨界點.

圖9 擺桿Vpython建模實體

圖10 不同初始偏轉(zhuǎn)角度下擺桿實際角位移曲線與簡諧運動角位移理論曲線對比圖

令縱坐標(biāo)表示角速度，橫坐標(biāo)表示角位置，圓棒的運動狀態(tài)對應(yīng)坐標(biāo)軸上的點，這樣的圖形叫做相圖.相圖沒有直觀地顯示角速度、角位置與時間的關(guān)系，卻能反映出運動的全局性.相圖法是非線性力學(xué)中最基本的研究方法[11].將不同初始條件下擺桿的運動狀態(tài)用ω-θ相圖表示如圖11所示.

圖11 不同初始條件下擺桿受重力矩的運動狀態(tài)ω-θ相圖

由圖11可知，當(dāng)θ0較小時，相圖軌線為圓形閉合曲線，表明運動為簡諧振動；θ0越大，軌線越“方”.

3.2 考慮阻尼和周期性外力矩時擺運動特性的研究

在經(jīng)典力學(xué)中，已知物體初始狀態(tài)和物體運動過程中受力情況，可以通過牛頓運動定律和動量定理預(yù)測物體的運動；但自然界中并非所有物體運動都是確定且唯一的，即使給定了初始條件，物體依然可能有多條運動軌跡，非線性系統(tǒng)就是典型例子.當(dāng)給任意角度復(fù)擺一個周期性外力矩并考慮阻力，使其受迫振動，當(dāng)驅(qū)動力矩達(dá)到一定值時擺的運動狀態(tài)無法預(yù)料，且對初始條件具有敏感性，這在物理上被稱為混沌.混沌的一種著名表述是“蝴蝶效應(yīng)”，初始條件的極小偏差，都將可能會引起結(jié)果的極大差異.

受迫振動微分方程為

(7)

C為阻力系數(shù)，M為力矩最大值，ωp為力矩周期變化的角頻率，上式化簡為

(8)

圖12 截取不同初始條件下擺桿受迫阻尼振動的ω -θ相圖

由于疊加原理、傅立葉變換等工具不適于求解非線性方程，而計算機(jī)性能的提高使得研究者一般采用計算機(jī)進(jìn)行迭代運算對其進(jìn)行數(shù)值求解.通過有阻尼受迫振動的模擬，發(fā)現(xiàn)在一定范圍內(nèi)即使初始條件f變化很小，卻會引起復(fù)擺的運動形式有很大的不同.例如當(dāng)f=1.505時，相軌跡為無規(guī)則的發(fā)散曲線；f=1.508時，相軌跡為收斂的三個極限環(huán).與當(dāng)年洛倫茲在利用計算機(jī)研究地球大氣時提出的“蝴蝶效應(yīng)”相類似，對于十分接近的兩個f的值，初始時相軌跡重疊，當(dāng)時間積累到一定程度后，相軌跡會出現(xiàn)較大偏差：設(shè)定模擬精度dt=0.001 s，f初始值相差0.00001時，程序運行200 s后，兩個復(fù)擺的運動將截然不同.隨著f值的變化，相軌跡運動呈現(xiàn)出由收斂、有序走向混沌的狀態(tài).混沌狀態(tài)下相點軌跡時而做周期運動，時而做非周期運動；時而向外延伸，時而在某一位置纏繞.對于某一時刻而言復(fù)擺運動具有不確定性(如圖11(e)—(g)).如果不斷增大f值，研究f>1.7以后的狀態(tài)(不再贅述)，收斂的極限環(huán)與混沌相圖將繼續(xù)交替出現(xiàn)，說明復(fù)擺的周期性運動與混沌狀態(tài)交替出現(xiàn).

簡單的復(fù)擺系統(tǒng)實際上的運動卻“不簡單”.混沌的發(fā)現(xiàn)對經(jīng)典力學(xué)造成了沖擊，同時也令其不斷完善和發(fā)展. 隨著人們對混沌理論的深入了解，將會給我們生活帶來有利的應(yīng)用. Vpython模擬可以檢驗非線性系統(tǒng)相關(guān)結(jié)論：1)疊加原理不再成立；2)初始條件的不同，會導(dǎo)致很大差異的運動形式；3)可能出現(xiàn)完全隨機(jī)的混沌狀態(tài).

4 結(jié)論

本文利用VPython對剛體運動和剛體擺動進(jìn)行運動模擬和數(shù)值模擬. 利用轉(zhuǎn)動定律和角動量定理，對恒定力矩和周期性變化力矩下剛體的定軸轉(zhuǎn)動運動狀態(tài)進(jìn)行了探究，并進(jìn)一步觀察了剛體做平面平行運動的軌跡. 探究了大、小角度下復(fù)擺的運動方程和周期差別. 有趣的是，通過對受迫阻尼振動的模擬，檢驗了非線性系統(tǒng)相關(guān)結(jié)論.

VPython不僅將課本知識形象化、可視化，加深了理解，還能快速解決學(xué)生學(xué)習(xí)中產(chǎn)生的疑問，激發(fā)學(xué)生對問題更深入地探索與思考. 在計算機(jī)技術(shù)迅猛發(fā)展，各種3D建模、數(shù)值分析軟件“百花齊放”的今天，培養(yǎng)學(xué)生能夠熟練運用功能強(qiáng)大的軟件輔助學(xué)習(xí)和研究,對深化大學(xué)物理課堂教學(xué)改革具有一定啟示.