[やること]
微分方程式で表現されるカオスをPythonで解いて、軌跡をプロットします。
[カオスについて]
カオスは式で明確に定義されているのに、初期値の選び方で将来の状態が予測できない現象のことです。
参考文献[1][2]がわかりやすいです
[結果1]
ローレンツ方程式
p,r,b=10,28,8/3
初期値 x,y,z=1,1,1
初期値 x,y,z=1.0001,1,1
初期値 x,y,z=1.001,1,1
初期値 x,y,z=1.01,1,1
横軸に時間、縦軸にxを取ったプロットです。
初期値が少し変わるだけで軌跡が大きく変わることがわかります。
[結果2]
レスラー方程式
a,b,c=0.2,1.0,5.7
初期値 x,y,z=1.0,0.0,0.0
初期値 x,y,z=8.01,0.0,0.0
初期値 x,y,z=8.02,0.0,0.0
[プログラム]
微分方程式を解くライブラリodeintについては、参考文献[3]がわかりやすいです
import numpy as np from scipy.integrate import odeint import matplotlib.pyplot as plt from mpl_toolkits.mplot3d import Axes3D #レスラー方程式 def rossler(var,t,a,b,c): x_next=-var[1]-var[2] y_next=var[0]+a*var[1] z_next=b*var[0]-(c-var[0])*var[2] return np.array([x_next,y_next,z_next]) #ローレンツ方程式 def lorentz(var,t,p,r,b): x_next=-p*var[0]+p*var[1] y_next=-var[0]*var[2]+r*var[0]-var[1] z_next=var[0]*var[1]-b*var[2] return np.array([x_next,y_next,z_next]) if __name__ == "__main__": #レスラー方程式のパラメータと初期値 #abc=(0.2,1.0,5.7) #xyz_first=[1.0,0.0,0.0] #ローレンツ方程式のパラメータと初期値 prb=(10,28,8/3) xyz_first=[1.01,1,1] #ステップ数 t=np.linspace(0,1000,100000) #微分方程式を解く #var=odeint(chaos,xyz_first,t,args=abc) var=odeint(lorentz,xyz_first,t,args=prb) #グラフの枠を作っていく fig = plt.figure() ax = Axes3D(fig) #軸にラベルを付ける ax.set_xlabel("X") ax.set_ylabel("Y") ax.set_zlabel("Z") #プロットして表示 ax.plot(var[:, 0], var[:, 1], var[:, 2], marker="o", linestyle='None') plt.show()
[参考文献]
[1]カオスアトラクタとは
http://www.gifu-nct.ac.jp/elec/deguchi/sotsuron/hisaki/node9.html#eqresura
[2]カオス
http://www.aoni.waseda.jp/yuuka/Sim/Chaos.html
[3]Pythonで運動方程式を解く(odeint)
https://qiita.com/binaryneutronstar/items/ad5efa27fd626826846f