[glossで実験] ファン・デル・ポール振動子

ファン・デル・ポール振動子は相空間¹にリミットサイクルが存在する微分方程式によって記述される振動子です。

\frac{d^2x}{dt^2} - \mu (1 - x^2)\frac{dx}{dt} + x = 0

これをHaskellで実装すると以下のようになるでしょう。

type Point = (Double, Double)

vanderpol :: Float -> Point -> Point
vanderpol dt (x, x') = let mu = 2.5
                           x'' = mu * (1 - x^2) * x' - x
                       in (x + dt * x', x' + dt * x'')

Pointは(位置, 速度)という座標を表すタプルです。微分方程式から2階微分の値を計算し、オイラー法によって次のステップの座標を計算しています。

この微分方程式の時間発展を実際にgloss²で描画して眺めてみましょう。

main :: IO ()
main = simulate inWindow white 24 initModel draw (\_ -> step)
  where
  inWindow = InWindow "Haskell Day 2018" (640, 480) (100, 100)
  initModel = [(1.0, 0.0)] -- 初期状態
  draw = scale 50 50 . line -- そのままだと小さいので50倍に拡大する
  step dt ps@(p:_) = vanderpol dt p : ps

(1.0, 0.0)から始めて時間発展をlineで描画しています。stepはvanderpolによって計算した次のステップの状態を今までの計算履歴のリストの先頭に付け加えています。これを実行すると以下のようなシミュレーションになります。

みなさんもぜひ自分で実装してみて初期状態や$\mu$の値を変えたりして遊んでみてください。

[参考] Wikipedia: ファン・デル・ポール振動子

1. x軸に位置, y軸に速度を取った空間 ↩

2. glossではじめるグラフィック描画 :: Haskell入門の次に読む記事 - Qiita ↩