理論推導參考該作者：https://www.jianshu.com/p/844287eb0b34

最終運用分析力學的方法，推導出如下的動力學方程

首先，在GGB中輸入初始參數(shù)，如下：

g=9.8

懸掛點A=(0,0)

擺長l_1=10

擺長l_2=10

質(zhì)量m_1=5

質(zhì)量m_2=5

終止時間t_f=100

球1的初始角度=pi/6

球2的初始角度=pi/2

球1的初始角速度=0

球2的初始角速度=0

動畫播放速度speed=slider(0.1,1,0.1)

滑動條a=slider(0,1)

設(shè)置a的播放速度為speed

輸入函數(shù)，如下：

xa(θ_{1},θ_{2},ω_{1})=l_{1} ω_{1}^(2) sin(θ_{1}-θ_{2})-g sin(θ_{2})

xb(θ_{1},θ_{2},ω_{2})=-(l_{2} m_{2} ω_{2}^(2) sin(θ_{1}-θ_{2})+(m_{1}+m_{2}) g sin(θ_{1}))

θ_{1}'(t,θ_{1},θ_{2},ω_{1},ω_{2})=ω_{1}

θ_{2}'(t,θ_{1},θ_{2},ω_{1},ω_{2})=ω_{2}

ω_{1}'(t,θ_{1},θ_{2},ω_{1},ω_{2})=((xb(θ_{1},θ_{2},ω_{2})-m_{2} cos(θ_{1}-θ_{2}) xa(θ_{1},θ_{2},ω_{1}))/(l_{1} (m_{1}+m_{2})-m_{2} l_{1} cos^(2)(θ_{1}-θ_{2})))

ω_{2}'(t,θ_{1},θ_{2},ω_{1},ω_{2})=((xa(θ_{1},θ_{2},ω_{1})-l_{1} cos(θ_{1}-θ_{2})*((xb(θ_{1},θ_{2},ω_{2})-m_{2} cos(θ_{1}-θ_{2}) xa(θ_{1},θ_{2},ω_{1}))/(l_{1} (m_{1}+m_{2})-m_{2} l_{1} cos^(2)(θ_{1}-θ_{2}))))/(l_{2}))

這里為了偷懶，我并沒有對方程進行化簡，有興趣的可以自行化簡，消去xa與xb以及合并同類項。

接著，數(shù)值解微分方程，如下：

解常微分方程組({θ_{1}',θ_{2}',ω_{1}',ω_{2}'},0,{θ_{01},θ_{02},ω_{01},ω_{02}},t_{f})

對解出的第一和第二個方程進行描點，例如：

B=描點(numericalIntegral1,a)

C=描點(numericalIntegral2,a)

時間t=x(B)

球P_{1}=(l_{1} sin(y(B)),-l_{1} cos(y(B)))

球P_{2}=P_{1}+(l_{2} sin(y(C)),-l_{2} cos(y(C)))

線f=線段(A,P_{1})

線h=線段(P_{1},P_{2})

大致模型制作完成，剩下的美觀因素，自行優(yōu)化即可。

成品如下：