第14・15回(2012年2月1日, 2月8日)


今日の演習

オイラー法等を用いて、幾つか現実の現象をモデル化した方程式系を解いて、その現象を解析してみます。

今日は以下の内容について演習を行います.


今日の目標


今回取り扱う対象は、前回の続きである

 。機械や建造物、材料の変形による弾性的性質(伸び縮みとその反動)を考える基礎となる、バネ振動
 。複数の化学物質を混ぜ合わせる事で起こる、化学反応によるリズム現象(体内リズムのモデル)
 。「食う」[食われる」関係で構成される生態系に生じるリズム現象
 。気象等、長期的予測が困難な現象に現れるカオス現象

と、未直な所でもしばしば現れる「フラクタル図形」について、順次紹介する予定です。


化学反応が引き起こすリズム現象(体内リズムの数理モデル)

我々の体は、いろいろなリズムを刻んでいます(いわゆる体内時計の活動です)。有名なのは朝昼夜に合わせて体の状態を変えていく概日リズムですが、体内にはもっと長い周期、短い周期のリズムが多数あり、それは人間だけでなくあらゆる生物に存在しています。そのようなリズムは、生物が取り込んだ外部からの栄養が、相互に様々な化学反応を起こす過程で生じる、化学成分の量の時間変動として生じます。(生体内の状態は非常にダイナミックに変化し続けているのです。)

ここでは、以下のような簡単な反応過程が、体内リズムと似た性質のリズム現象を引き起こす事を見ていきます。

次のような化学成分 X(activator)と Y(inhibitor)が起こす(モデル)化学反応を考えてみます。

分子 X は外と流入(増加)流出(減少)し、Y は反応速度 b で X から生成され増加する(その分 X は増加する)一方、2つの X と Y が出会うと Y は X に戻される(Y は減少、その分 X が増加する)とします。この場合、分子 X, Y の個数を x, y とすると、x、y の増減(時間変化) dx/dt, dy/dt, は x, y, を用いて次のような関係を満たします。

 dx/dt = 1 - (1 + b)*x + x2y
 dy/dt = b*x - x2y 

この方程式は「ブリュッセレーター(ブラッセレーター)モデル(方程式)」と呼ばれています。(各自調べてみてください)

これについて、(1) b < 2 の場合、(2) b > 2 の場合を、適当な初期条件 (x,y)(幾つか試す事)から数値計算し、解の様子(各化学成分量 x, y, の時間変動)を比較してみてください。

この方程式は、ベロソフ・ジャボチンスキー反応(BZ反応)、と呼ばれる化学反応現象の本質的な部分を、数理的に表現した物になっています。BZ反応について、調べてみてください。細かい所には違いは幾つもあるのですが、体内時計の仕組みも、このような化学成分の個数の変動の仕方で制御されています。

生態系におけるリズム現象(捕食者-被食者の個体数変動)

生態系は、捕食者-被食者関係から成る食物連鎖の微妙なバランスによって、多様性を実現維持しています。そのような生態系のある部分を切り取ってみると、そこにもダイナミックな現象が潜んでいます。

次のような 被食者 X (草、プランクトン、等) と 被食者 Y(草食動物、魚、等)の関係を考えてみます。

被食者 X は周りの栄養を基に増殖(増加)し、捕食者 Y と出会うと食べられ減少する。逆に Y は X と出会うとそれを食べた分だけ増殖(増加)し、一方で一定の割合で寿命で減少するとします。この場合、生物種 X, Y の個体数を x, y とすると、x、y の増減(時間変化) dx/dt, dy/dt, は x, y, を用いて次のような関係を満たします。

 dx/dt = 3x - xy
 dy/dt = -y + xy

この方程式は「ロトカ・ボルテラモデル(方程式)」と呼ばれています。(これも各自調べてみてください)

これについて、適当な初期条件 (x,y)(幾つか試す事)から数値計算し、解の様子を比較してみてください。(ロトカ・ボルテラモデルの計算では、刻み幅を十分小さくする事。)

また「ブリュッセレーター(ブラッセレーター)モデル(方程式)」での解の様相との違い(特に初期条件依存性)も考察してください。

バネ振動その3:強制振動(カオス)

先程調べたバネの強制振動ですが、先程とはちょっと性質の違うバネ(小さい伸び縮みでは弾性力が弱いけど、やや大きく伸び縮みすると、急激に弾性力が強くなるバネ。)ではどのような振動が起こるのか、解析してみましょう。

そのようなバネの振動は例えば次の微分方程式

dy/dt = v

dv/dt = - 0.3*v - y3 + B*cos(t)

でモデル化できます。これを数値的に解いて、B = 0, 1, 6, 12 の場合の y の時間変化を調べてください.また、y と v の関係(gnuplot で横軸を y, 縦軸を v として)も plot してみて下さい. これは、バネが大きく伸び縮みすると急激に硬くなるような場合の、強制振動系です. 特に y, v の初期値が少し違う場合の y の時間変化を比較してください. 先程の強制振動とどのような特徴の違いがあるでしょうか.

*実はこの方程式は「カオス解」を持ちます。ハードディスクから原子炉まで、様々な機械の「振動」が、このような方程式を応用しながら解析されています. ところでカオスって何でしょう?調べてみてください。


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

地球上の気象は、地球を取り巻く大気の流れ、地熱や太陽の熱の流れ、水蒸気ー水の相転移といった様々な出来事が絡み合うことで、非常に複雑に変動します。その肝と成るダイナミックスを記述する方程式を、気象学者のエドワード・ローレンツが次のように導きだしています。

dx/dt = P*(y - x)

dy/dt = -y - x*z + R*x

dz/dt = x*y - b*z

この方程式は「ローレンツモデル(方程式)」と呼ばれています。(これもどういう物なのか各自調べてみてください)

この方程式をある初期条件で数値計算し、次にほんの少しだけ違う初期条件で巣値計算をするとした時、先にした結果からつぎの計算の結果(変化の仕方)がどの程度、予測できるでしょうか?
例えば P = 16, b = 4、R = 35, などとし、また初期値を x(0) = y(0) = z(0) = 5 等と、その少しずれたものとし、 時刻 0≦ t ≦100 で方程式を解いて x, y, z の時間変化をグラフにしてみてください. (実はこの解もカオス解です。)


写像力学系:生態系のカオス

次の漸化式で表される数列 (0 < Xn < 1) を求めるプログラムを作成し、n と Xn の関係を見てみてください.

Xn+1 = A * Xn * (1.0 - Xn)

A < 1 なら Xn → ∞ → 0 となるのは、すぐ分ると思います. では A > 1 なら...

A = 0.5, 2.0, 3.5, 3.7, 4.0 等(A : 0 〜 4)それぞれの場合について、初期値 X0 を 0.3, 0.6, 0.61 等(0 〜 1)と変化させ、実行、描画して見比べてください. X0, X1, X2 あたりの様子から X9, X49 あたりの値が予測可能なら予測せよ.

これはロジスティック写像と呼ばれる、生物の個体数変動を表すモデルの一つです.

例えば Xn をある年(n 年)のある地域でのとある生物の(昆虫など)の数密度とすると、Xn が多いと卵を産む可能性のある個体数が多い反面、餌が競合するため成熟して実際に卵を産める数は限られます。その双方の効果を考慮したのが、右辺の

   Xn * (1.0 - Xn)

という項で、翌年 (n+1 年) の数密度 Xn+1 はそれに比例する事になります。餌の量が多い方が、より産まれる卵が多くなると考えられるので、右辺の係数 A は餌の供給量等と考えられます。実際この餌の供給に応じて、どのように個体数変化が起こるのか、見てみてください。(詳しい事は、各自調べてみてください。)


シェルピンスキー・ギャスケット:貝殻のフラクタル

自己相似的な階層構造を持つ図形は、フラクラル図形と呼ばれます(もう少し厳密な定義が気になる人は調べてください。)フラクタルな図形(模様)も、自然界でいろいろな所で顔を出しています。今回はその中の典型例である、シェルピンスキー・ギャスケットを、これまでの演習で紹介した知見を用いて描いてみます。

この条件で x[i][j] = 1 となる ( i , j ) を( printf() で )出力し、 gnuplot を用いて表示すると、以下のような自己相似的な図形が描かれます。

例えばこのような図形は、巻貝の殻の模様等でも現れます。(j の増加が殻の成長と見なすと、上記のルールに従って色素が ON (1) - OFF (0) されることで、生化学的にこのような模様が形成される事が分かります。)



期末レポート課題

計算数学演習では,期末テストの代わりに期末レポート課題を出題します.

課題:生命現象とリズムやカオスに関連する現象について調べて、自分の考えも交えてまとめてください。

  • 提出期限は 2/13 正午, 提出場所は数学事務(B709)のポストとします。
  • 単位が必要な人は必ず提出してください.
  • 必ず手書きでまとめてください。(字は丁寧に!)
  • 他の人とは違うオリジナリティを出すようにしてください。
  • レポートはA4サイズとし,表紙には「計算数学演習 期末レポート」と明記し,学籍番号,名前を書くこと.また,レポートはホチキスなどでばらばらにならないようきちんととめること.
  • 表紙合わせて4ページ程度、簡潔かつ過不足なく、読みやすくまとめて下さい(適宜、絵や図など入れてもかまいません.).
  • レポート作成において,何か参考にした本、webページなどがある場合には,参考にした本の名前と著者名、webページのアドレスなどをレポート中に記載すること.
  • レポートは返却しないので,手元に残しておきたい人は各人前もってコピーをとること.
  • 粟津は通常、理学部 A113 にいます.質問等は御遠慮なく.

    戻る