第13回(2012年1月25日)


今日の演習

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

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


今日の目標


今回は前回までの内容の応用として、幾つかの現実の現象のモデル方程式系を取り扱います。具体的に取り扱う対象は、

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

を順次紹介する予定です。


バネ振動その1:減衰振動

下の図のように、壁に一端を固定されたバネの他端におもりがつけられている状況を考えます.

バネとおもり

このおもりを少しだけ右に引っ張ってから手をはなす.おもりの中心の時刻 t におけるつり合いの位置からのずれを y(t) とします.また、時刻 0 でのおもりの位置を y(0)=1、おもりの速度を y'(0)=0 とすると、課題2と全く同じ、つぎの単振動の問題となります.

式(1)   (13)

これではおもりに働く力はバネの復元力のみで、ある意味非現実的です.今回は復元力以外におもりの速度 y' に比例する抵抗が働く場合を考えましょう(おもりと床との間の摩擦力のようなものをイメージ).改良されたモデルは次のようになります(2ky' の項が新たに加わったもので、摩擦のようなものをあらわしている).

式(2)   (14)

ここで k,ω はそれぞれ抵抗力、復元力の強さを表す正の定数です.課題は、(14)式の問題を連立1階常微分方程式に変換し、オイラー法を用いて解き、gnuplotをつかってグラフにしてください.ω = 1 として、k < ω でのグラフと k > ω でのグラフを見比べて違いに関して考察してみて下さい. (実行例はこちら)



バネ振動その2:強制振動(リズム)

物体の性質を調べるための方法として、叩いたり揺らしたりした時の応答(挙動)を分析する方法があろいます。ここではそのような解析の基礎となる、上のような摩擦抵抗のあるバネ振動子に、周期的な外力を加えた(振動を与えた)強制振動系を考えます。

(1) 次の微分方程式

dy/dt = v

dv/dt = - 0.3*v - y + cos(w*t)

を、幾つか適当な w の場合について, y, v に適当な初期値を与え数値的に解き, y の時間変化を調べよ(時刻の刻み幅を変えずに上限を少し大きめ(100程度)にして下さい.).実行例はこちら

(2) w = 0.1, 0.2, 0.3, ... 2.0 の各 w 毎について(1)を数値的に解くプログラムを作成し、下のような、バネの振動の振幅 y の w 依存性を表すグラフを作成せよ。
( 方法:例えば w = 0.1*i ( i は整数で 0 ≦ i ≦ 20 ) とし, i についての繰り返し文を用いて w を変化させつつ, 各 w 毎について(1) を数値的に解くようにプログラムを作成する. )
下の結果は 横軸を w, 縦軸を y ( if文を使って t=70 から t=100 での 時間一定間隔毎の y の値のみを表示 )とした場合、以下のような図が描かれる (データは点でプロットしています.). 各自確認する事. 

グラフ

. この図から、バネを揺らす際には、与える揺れの周期を適度に(バネ自身が持っている固有の周期に近く)調整することで、より大きな振幅での揺れを実現できる事を意味しています。これは共振(共鳴)とをばれる現象で, ブランコで体験できたり、吊り橋が壊れたり、地震での建物の揺れで問題になってたりしています。(普通、吊り橋や建物では、共振が起こらないように設計、材料選択がなされています。)逆に物体の性質が分からない場合、この共振によって物体が本来持つ振動特性が分かることになります。

  



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

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

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

次のような化学成分 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 の時間変化をグラフにしてみてください. (実はこの解もカオス解です。)



期末レポート課題

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

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

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

    戻る