アトラクター

　"ストロガッツ 非線形ダイナミクスとカオス"を読み始めた。5月末までにこれの8章までを読むことを目標とする。誓います！（強気）ただしわからない部分をとばすのは可(えー

　備忘録もかねて、必要の有無にかかわらずコーディングできそうなところをjupyterにまとめていくことにした。そして先日githubを始動させたので、そこにアップロードすることでモチベーションを保っていくことにした。

https://github.com/nkknj/strogatz

　本日は2章の練習問題を解いていたのだが、その中でアトラクターの表示に手間取った話。

import numpy as np
import matplotlib.pyplot as plt

x=np.arange(-3.5, 1, 0.1)
x_d=np.exp(-x)*np.sin(x)
plt.plot(x, x_d, linestyle='-', color='b', label='x_d')

#x軸上に流れの方向へ三角形を描出
for i in range(x.shape[0]):
    if x_d[i]>0:
        plt.plot(x[i], 0, marker='>', color='dodgerblue')
    elif x_d[i]<0:
        plt.plot(x[i], 0, marker='<', color='dodgerblue')
    else:
        #流れがゼロの点
        if x_d[i+1]>0:
            #直後が正の時
            plt.plot(x[i], 0, color='red', marker='s', markersize=5)
        else:
            #直後が負の時
            plt.plot(x[i], 0, color='red', marker='o', markersize=5)

plt.legend()
plt.xlabel('x')
plt.grid()
plt.show()

　この問題の直前までは単純なxのn次式とか三角関数だったので、xを0.1おきとかpi/10おきとかで設定すればピタッとx_d=0になる点が出てきて、簡単にアトラクターだとわかった。しかし指数関数が出てくるとなかなかぴったりゼロとはいかない。今後のためにも他の方法を考えた方がよい。

　そんなわけで、2点間でゼロをまたいだ時に、その2点を結んだ線分とx軸の交点をアトラクターとしてマーカーを表示する感じにした。

import numpy as np
import matplotlib.pyplot as plt

x=np.arange(-3.5, 1, 0.1)
x_d=np.exp(-x)*np.sin(x)
plt.plot(x, x_d, linestyle='-', color='b', label='x_d')

#x軸上に流れの方向へ三角形を描出
for i in range(x.shape[0]):
    if x_d[i]>0:
        plt.plot(x[i], 0, marker='>', color='dodgerblue')
    elif x_d[i]<0:
        plt.plot(x[i], 0, marker='<', color='dodgerblue')
    else:
        #流れがゼロの点
        if x_d[i+1]>0:
            #直後が正の時
            plt.plot(x[i], 0, color='red', marker='s', markersize=5)
        else:
            #直後が負の時
            plt.plot(x[i], 0, color='red', marker='o', markersize=5)

i=0
while i+1<x.shape[0]:
    if x_d[i]*x_d[i+1]<0:
        if x_d[i]>0:
            plt.plot(x[i]+(x[i+1]-x[i])*(x_d[i]/(x_d[i]-x_d[i+1])), 0, color='red', marker='o', markersize=5)
        else:
            plt.plot(x[i]+(x[i+1]-x[i])*(-x_d[i]/(x_d[i+1]-x_d[i])), 0, color='red', marker='s', markersize=5)
   
    i+=1
   
plt.legend()
plt.xlabel('x')
plt.grid()
plt.show()

　単純にwhile文加えて、後ろから1つずつ確認する仕組み。ただしぴったりゼロになった場合の対処は、x軸に△表示するついでにお願いする形となっている。

　この方法で対処できないのは、x軸に接しているのにx_d=ジャストゼロが得られなかった場合である。現時点では練習問題にそのような問題がみられないので問題とはならないが、そのような場合どうすればいいのだろうか。