計算数学演習第14回

Chaosについて+期末レポート

今日の演習

  •  具体的な問題をといてもらう.

なぜ天気の長期予報は難しい?:大気の流れのカオス

地球上の気象は、地球を取り巻く大気の流れ、地熱や太陽の熱の流れ、水蒸気ー水の相転移といった様々な出来事が絡み合うことで、非常に複雑に変動します。
気象学者のエドワード・ローレンツは特殊な状況下の大気の流れを表す流体力学の方程式を簡略化することにより次のような方程式を導きだしました。

(1)   \begin{eqnarray*} \frac{dx}{dt} &=& P(y-x) \\ \frac{dy}{dt} &=& -y -xz + Rx \\ \frac{dz}{dt} &=& xy - bz \end{eqnarray*}

この方程式は「ローレンツ方程式」と呼ばれていて,カオス最も有名なモデルの一つです。

カオスの厳密な定義はここでは触れません(実は人によって定義が違ったりします).
ですが,カオスのもっている初期値鋭敏性という性質(ちょっとした初期値のずれが時間とともに指数関数的に増大すること)を理解してもらうために,以下の課題(特に課題2)をやってもらいます.

題1

ローレンツ方程式(以降,全てパラメータはP = 10.0, b = 8.0/3.0, R=28) を初期値(x(0)=1.0, y(0)=0.5, z(0)=0.5)のもとオイラー法で解き,どのような3次元空間上でどのような軌道を描くかgnuplotを用いて確認せよ.
この際,時間の刻み幅h=0.001としt=100まで計算するものとする.
ただし,このままだと若干出力されるデータが多すぎるので,tが0.01変化する毎に t, x, y, zを出力すること.

splot “output.txt” using 2:3:4 w l

と打ち込めば,3次元のグラフが得られる.(第一成分はtなので,using 2:3:4とした)
また,(t, x)のグラフをみていると,二つの渦(?)の間の移り変わりが周期的ではない事がわかる.

課題2(ローレンツ方程式の初期値鋭敏性)

微分方程式は初期値を与えるとその後の振る舞いが一意的に定まる.
だが,物理現象のモデルであると考えると,初期値というものを無限に精度良く知る事は一般には不可能である.
(例えば,計算の為に現在の気温が必要だったとしても,我々の身の回りの温度計ではせいぜい10分の1度程度の精度までしかわからない.
25.3℃と読み取れても,もしかしたら25.32℃かもしれないし25.29℃かもしれない)
初期値がわずかに違っていたとしても,tが小さいうちはそのずれはあまり問題にならないはずである.
では,tがどれくらい大きくなったら,その初期値の曖昧さが無視できなくなるだろうか?
初期値を少しだけずらすと,その後の振る舞いがどのように変化するかを調べる.
初期値を

x(0)=1.0+k, y(0)=0.5, z(0)=0.5

として,

k=0のときとk=0.0001のときのデータを比較してみよ.
3次元のグラフをみるより,x(t)のグラフをみるとわかりやすい.
また,そのずれを10分の1にしてみたら(k=0.00001とおいたら)どうなるだろうか?
kを10分の1にしたら,10倍長い時間精度よく計算できるといえるだろう?


発展問題(ローレンツ方程式)

ローレンツ方程式において,x(t), y(t), z(t)の時系列は明らかに周期的ではなく、一見不規則に移り変わっているように見える。
だが、実は完全に不規則な訳ではなく、ある種の規則性があること示す為に、次のことをやってもらう。

t > 100という条件のもと、tが小さい方から数えてn番目のz(t)の極大値をZ_nとおく。
(Z_n, Z_(n+1) )をn =1, 2, … , 200まで出力するプログラムを作成せよ。
また、得られた出力をもとに(Z(n), Z(n+1))の散布図(データを線で結ばない図)を作図し,データがどのような形の分布を持つのか確認すること。

例: Z(1) = 3, Z(2)= 5, Z(3) = 6, Z(4) = 4,..のときは、
二次元平面上に点(3, 5), (5, 6), (6, 4), .. がプロットされる。 

これらの点は,二次元平面状にまばらに広がっているのではなく、
ある種の構造を持っていることがわかるはずである。
これは、Z(n)からZ(n+1)が予測可能であることを意味している。


期末レポート問題

ローレンツ方程式についてレポートを作成してもらう.
紙のレポートを作成して,数学事務のレポートボックス(2月6日までに設置予定)に提出してください.
レポート内の図は,手書き,もしくは自分でgnuplotで描いたものを使うこと.

また,レポート内容に関連するプログラムの内一つをrep20180126.cという名前で保存して,bb9から提出してください.
ただし,数値計算はrunge-kutta法を用いて,時間の刻み幅hはh<0.02であること.
(オイラー法を用いた場合も受け付けますが1割程減点します)


(冒頭にコメントとして,どのようなプログラムであるか説明すること.高度な内容であるほどポイントは高いです.)

レポート
・ローレンツ方程式の場合
課題1,課題2,発展問題から,ローレンツ方程式(P = 10.0, b = 8.0/3.0, R=28)についてどのようなことがわかったかまとめてください.
(発展問題までできていない人は,できている範囲でわかったことを議論してください.)

ローレンツ方程式はカオスで有名ですが,p,b,Rの値によってはカオス以外の様々な解(一点に収束する場合や,周期解をもつ場合など)が現れるので,
余裕のある人はそれについても調べてみてください.(やってくれたら,少し加点します)

レポートは手書きでもwordでもtexを用いたものでもかまいません.
ただし,参考にした文献やWebサイトがある場合は,最後に参考文献として記載してください.

諸注意:
レポートはwordでもtexを用いたものでも手書きでもかまいません.
数値計算結果のグラフは自分で作成したgnuplotの図を用いてください.
図が思った位置に表示できない場合や手書きのレポートの場合は,
図だけ印刷し切り抜いてレポートに(糊やテープで)貼付けても結構です.

イメージを伝える為に手書きの図を使うことは許可しますが,
手書きのグラフを「数値計算の結果」としては認めませんので注意してください.

ただし,参考にした文献やWebサイトがある場合は,最後に参考文献として記載してください.

提出期限:2018年2月15日(木)午後4時半(時間厳守)
提出場所:数学科事務レポートボックス