8
6

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

More than 5 years have passed since last update.

【カオスアート】二重振り子のカオスアートを制御する♪

Last updated at Posted at 2018-08-12

二重振り子のような簡単な系がなぜカオスのような振舞になるのかは非常に興味深い。
参考の①に非常に興味深い説明があるので参照してほしいが、ここではカオスアートを制御することにより、その振舞を見ていこうと思う。

前回、微分方程式に過誤はあったが、カオスアートがある程度制御できることが分かった。
また、二重振り子には少なくとも以下の微分方程式が導出されている。
①単振り子連結
②棒連結
③単振り子連結でも見た目異なる微分方程式が導出されている
これらの系でどういうカオス領域あるいはその結果としてのカオスアートが生成されるのかを調査しようと思う。

やったこと

(1)微分方程式を無次元量で表す
(2)3つの運動方程式から得られるカオスアート
(3)2つの単振り子連結二重振り子の同値性検証
(4)ちょっと驚いたカオスアート

(1)微分方程式を無次元量で表す

そもそも単振り子の運動方程式は

ml\frac{d^2\theta}{dt^2}=-mg\sin\theta
\frac{d^2\theta}{dt^2}=-\frac{g}{l}\sin\theta

上式は触れ角Θ、角速度ω(=dΘ/dt)、そして周期Tで特徴づけられる力学系である。
すなわち微小角における単振り子の周期Tは

T=2\pi\sqrt{\frac{l}{g}}

固有振動数ωは

\omega=\frac{2\pi}{T}
      =\sqrt{\frac{g}{l}}

この固有振動数を用いて、以下のように表現できる。

\theta=C_1\sin\omega t+C_2\cos\omega t

参考のように、この周期と固有振動数は角度が大きい領域でも系を特徴づける量であることが分かる。

【参考】
振り子

単振り子と同様に、各系の微分方程式を無次元量で表わした。こうすることにより系の振舞を連結する系の固有な振動をイメージしつつ、その絡み合いとしてのカオス領域をイメージしやすくなる。
ということで、以下の関数を利用した。

単振り子連結二重振り子

\theta_1'' = \frac{-\omega^2(\sin\theta_1 - M\sin\theta_2\cos\theta_{12}) -M\sin \theta_{12}(\theta_2'^2L+\theta_1'^2\cos\theta_{12})}{(1-M\cos^2\theta_{12})}
\theta_2'' = \frac{{-\omega^2(\sin\theta_2-\cos\theta_{12}\sin\theta_1)+(\theta_1'^2+ML\theta'^2_2\cos\theta_{12})\sin\theta_{12}}}{L(1-M{\cos^2\theta_{12})}}
def doublePendulum(y, t, L, omegaJ,M):
    theta1, omega1,theta2, omega2 = y
    theta12=theta1 -theta2
    #M=m2/(m1+m2)
    #L=L2/L1
    #omegaJ=g/L1    
    dydt = [omega1,
            (omegaJ*L*(-np.sin(theta1)+M*np.cos(theta12)*np.sin(theta2))-M*L*(omega1*omega1*np.cos(theta12)+L*omega2*omega2)*np.sin(theta12))/(L-M*L*np.cos(theta12)*np.cos(theta12)),
            omega2,
            (omegaJ*(-np.sin(theta2)+np.cos(theta12)*np.sin(theta1))+(M*L*omega2*omega2*np.cos(theta12)+omega1*omega1)*np.sin(theta12))/(L-M*L*np.cos(theta12)*np.cos(theta12))]
    return dydt

棒連結二重振り子

\theta_1'' = \frac{-\omega^2L(M_1\sin\theta_1 - M\sin\theta_2\cos\theta_{12}) -M\sin \theta_{12}(\theta_2'^2L+2\theta_1'^2\cos\theta_{12})}{(1-2M\cos^2\theta_{12})}
\theta_2'' = \frac{{(2\theta_1'^2+2ML\theta'^2_2\cos\theta_{12})\sin\theta_{12}-\omega^2(\sin\theta_2-2M_1\cos\theta_{12}\sin\theta_1)}}{L(1-2M{\cos^2\theta_{12})}}
def doublePendulum2(y, t, L, omegaJ,M,M1):
    theta1, omega1,theta2, omega2 = y
    #M=2*m2/(m1+4*m2)
    #M1=(m1+2*m2)/(m1+4*m2)
    #L=L2/L1
    #omegaJ=g/L1  
    theta12=theta1 -theta2
    dydt = [omega1,
            (omegaJ*L*(-M1*np.sin(theta1)+M*np.cos(theta12)*np.sin(theta2))-M*(2*omega1*omega1*np.cos(theta12)+L*omega2*omega2)*np.sin(theta12))/(1-2*M*np.cos(theta12)*np.cos(theta12)),
            omega2,
            (omegaJ*(-np.sin(theta2)+2*M1*np.cos(theta12)*np.sin(theta1))+(2*M*L*omega2*omega2*np.cos(theta12)+2*omega1*omega1)*np.sin(theta12))/(L-2*M*L*np.cos(theta12)*np.cos(theta12))]
    return dydt    

参考③の二重振り子の運動方程式

\theta_1'' = \frac{-\omega^2(\sin\theta_1 + M\sin(\theta_1-2\theta_2)) -2M\sin \theta_{12}(\theta_2'^2L+\theta_1'^2\cos\theta_{12})}{(1-M\cos2\theta_{12})}
\theta_2'' = \frac{2\sin\theta_{12}(M_1\theta_1'^2+\omega^2M_1\cos\theta_1+\theta_2'^2LM_1\cos\theta_{12})}{L(1-M{\cos2\theta_{12})}}
def dualPendulum(y, t, L, omegaJ,M,M1):
    theta1, omega1,theta2, omega2 = y
    #M=m2/(2*m1+m2)
    #M1=(m1+m2)/(2*m1+m2)
    #omegaJ=g/L1
    #L=L2/L1
    dydt = [omega1, 
            (-omegaJ*(np.sin(theta1)+M*np.sin(theta1-2*theta2))-2*np.sin(theta1-theta2)*M*(omega2*omega2*L+omega1*omega1*np.cos(theta1-theta2)))/(1*(1-M*np.cos(2*(theta1-theta2)))), 
            omega2, 
            (2*np.sin(theta1-theta2)*(omega1*omega1*M1+omegaJ*M1*np.cos(theta1)+omega2*omega2*L*M1*np.cos(theta1-theta2)))/(L*(1 - M*np.cos(2*(theta1-theta2))))]
    return dydt

(2)3つの運動方程式から得られるカオスアート

単振り子連結二重振り子

xp1
pen_xp1_digest.jpg
xp2
pen_xp2_digest.jpg

棒連結二重振り子

xp1
pen2_xp1_digest.jpg
xp2
pen2_xp2_digest.jpg

参考③の二重振り子の運動方程式

xp1
dual_xp1_digest.jpg
xp2
dual_xp2_digest.jpg

(3)2つの単振り子連結二重振り子の同値性検証

そもそも運動方程式の右辺の振舞を関数的に求めてもできるが、ここはカオスアートの同値性から検証してみた。
以下、同じパラメータに対してのカオスアートを比べてみる。
pen_flower.gif
dual_flower.gif
flower.jpg

上記の初期角度はだいたい10^-4~10^-5程度の領域で比較しており、ほぼ一致した振舞になっている。
他のカオスアートも角度が小さい領域はほぼ同じ振舞をする。
一方、上記のダイジェスト見てもわかると思うが、いくつか特徴的な絵が出ていない場合がある。以下のように、特に初期角度が大きい場合に差が異なるカオスアートが出現した。ということで、この微分方程式はほぼ単振り子連結二重振り子の運動方程式だといえるが、何か近似等をやっているようだ。
xp1dual_L11L21m1100m21y010.15601991952892066
plot_xp1dual_L11L21m1100m21y010.15601991952892066.png
xp1pen_L11L21m1100m21y010.15601991952892066
plot_xp1pen_L11L21m1100m21y010.15601991952892066n.png

xp1dual_L11L21m1100m21y010.19056316004055146
plot_xp1dual_L11L21m1100m21y010.19056316004055146.png
xp1pen_L11L21m1100m21y010.19056316004055146
plot_xp1pen_L11L21m1100m21y010.19056316004055146n.png

xp1dual_L11L21m1100m21y010.23275436927724738
plot_xp1dual_L11L21m1100m21y010.23275436927724738.png
xp1pen_L11L21m1100m21y010.23275436927724738
plot_xp1pen_L11L21m1100m21y010.23275436927724738n.png

(4)ちょっと驚いたカオスアート

これ何なんだろうというカオスアートがたくさん出てくるが、こういうの出るんですね。
これって、DLA(Diffusion-limited aggregation)で有名な樹枝状成長に類似しています。
plot_xp1Pen2_L11L21m1100m21y011.925438770268361e-07n.png
plot_xp1Pen2_L11L21m1100m21y011.2973504384052513e-07n.png
plot_xp1Pen2_L11L21m1100m21y011.4227175088568964e-06n.png
plot_xp1Pen2_L11L21m1100m21y015.233885221031016e-07n.png
plot_xp1pen_L11L21m1100m21y011.4227175088568964e-06n.png
そして、以下なんだか不思議なカオスアートをセレクトしてみました。
plot_xp1pen_L11L21m1100m241y015.233885221031016e-07n.png
plot_xp1pen_L11L21m1100m260y012.6057980132902392e-08n.png
plot_xp2pen_L1.0M9.8omegaJ0.6227838551490004y012.1115e-05n.png
plot_xp1Pen2_L11L21m1100m240y012.8576045055410773e-05n.png
plot_xp1Pen2_L11L21m1100m241y017.767774400335005e-05n.png
plot_xp1Pen2_L11L21m1100m21y011.7557728903871836e-08n.png
plot_xp2dual_L13L21m1100m220y012.6057980132902392e-08.png
plot_xp2dual_L13L21m1100m220y017.083293388161538e-08.png
plot_xp1pen_L11L21m1100m240y012.6057980132902392e-08n.png

案外、七角形や五角形が多いし、いわゆる原子のs状態しかも3sとか高次の状態に見えるのがなんとも。。。今回5000枚以上のカオスアートを見たが、まだまだ別な領域がありそうで、一年がすぐすぎてしまいそうです。

まとめ

・単振り子連結と棒連結二重振り子のカオスアートを見た
・二つの系では、ほぼ似たようなカオスアートが得られる
・低角度領域においても魅力的なカオスアートが得られることが分かった

・強連結な振動子の系は微小振動とか微小角だからと云って、線形近似するのはダメなようだ
・エネルギー保存則の検証をしたいと思う
・運動方程式について、関数のΘやω依存性から上記結論を検証する

おまけ

【参考】
二重振り子@カオス&非線形力学入門
複雑系の物理学(下)琉球大学理学部物質地球科学科集中講義
Double Pendulum@MyphysicLab.com

棒連結の微分方程式は以下のとおり

\theta_1'' = \frac{-g(m_1 + 2m_2)\sin\theta_1 + 2m_2g\sin\theta_2\cos\theta_{12} -2m_2\sin \theta_{12}(\theta_2'^2L_2+2\theta_1'^2L_1\cos\theta_{12})}{L_1(m_1+4m_2-4m_2\cos^2\theta_{12})}
\theta_2'' = \frac{{(2L_1\theta_1'^2(m_1+4m_2)+4m_2L_2\theta'^2_2\cos\theta_{12})\sin\theta_{12}-g((m_1+4m_2)\sin\theta_2-2(m_1+2m_2)\cos\theta_{12}\sin\theta_1)}}{L_2(m_1+4m_2-4m_2{\cos^2\theta_{12})}}
\theta_1'' = \frac{-\omega^2(M_1\sin\theta_1 - M\sin\theta_2\cos\theta_{12}) -M\sin \theta_{12}(\theta_2'^2L+2\theta_1'^2\cos\theta_{12})}{(1-2M\cos^2\theta_{12})}
\theta_2'' = \frac{{(2\theta_1'^2+2ML\theta'^2_2\cos\theta_{12})\sin\theta_{12}-\omega^2(\sin\theta_2-2M_1\cos\theta_{12}\sin\theta_1)}}{L(1-2M{\cos^2\theta_{12})}}

以下の無次元量とωを導入した。
無次元量で表したことにより、上式の見通しはよくなり、無次元量で記載されている参考の運動方程式と類似したものとなる。

L=\frac{L_2}{L_1}
\omega^2=\frac{g}{L_1}
M_1=\frac{m_1+2m_2}{m_1+4m_2}
M=\frac{2m_2}{m_1+4m_2}

Lagrangianは以下のKとUからL=K-Uで上記微分方程式が得られます。

K=\frac{1}{2}(m_1+4m_2)L_1^2\theta_1'^2+\frac{1}{2}m_2L_2^2\theta_2'^2+2m_2L_1L_2\theta_1'\theta_2'\cos\theta_{12}
U=(m_1+2m_2)gL_1(1-\cos\theta_1)+m_2gL_2(1-\cos\theta_2)
8
6
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
8
6

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?