伊達 >> 生体情報システム >> 課題一覧 >> 課題 2.5

課題 2.5

目的

説明を読むよりは「例題」を実践し「学習」してみてください．

例

1. octave により微分方程式の解を数値実験により求める方法

   dx/dt = -3x， x(0)=1，という微分方程式の解を求める．
   # から始まる行や，各行の # から後ろはコメントです．

   lsode-sample.m

   function dx = f(x,t)
           dx(1) = -3.0*x
   end
   
   t=linspace(0,10,101);    # 横軸の作成． 0 から 10 を 101等分
   x=lsode("f",1,t);        # 1 は x の初期値 x(0)． lsode は関数．こう書けば解いてくれる．
   plot (t,x,";x(t);");
   pause
   

   % octave lsode-sample.m 
   

2. Fitzhugh- Nagumo モデルを記述

   nagumo.m （古いバージョン　nagumo070426.m ）

   global theta gamma epsilon INPUT;
   theta= 0.01;
   gamma = 0.5;
   epsilon = 0.002;
   INPUT = 0.4;
   
   function dx = nagumo(x,t)
           global theta gamma epsilon INPUT;
           dx(1) = -x(1)*( x(1)- theta)*( x(1)-1.0 ) - x(2) + INPUT;
           dx(2) = epsilon*( x(1)  - gamma*x(2) );
   end
   
   t=linspace(0,2000,2000);
   x=lsode("nagumo",[0.5; 0.2],t); # [0.5;0.2] are initial values for x(1) and x(2)
   plot( t,x(:,1), "r;u;");  # u is plotted in red line
   xlabel ("t");
   ylabel ("u,v");
   hold on
   plot(t,x(:,2),"g;v;" );   # v is plotted in green line
   print -deps -FHelvetica:20  hoge.eps
   pause
   
   clf
   plot(x(:,1),x(:,2));   # see another picture, x-axis u, y-axis v.
   xlabel ("u");
   ylabel ("v");
   print -deps -FHelvetica:20  hoge2.eps
   pause
   

   % octave nagumo.m 
   

   
3. メモ
   プログラム上では u を x(1)， v を x(2) と記述している． v は u に比べ，非常にゆっくりと変化している場合を考える． epsilon の値が小さいのはそれを反映している． つまり，u が動いている間， v はある特定の値に固定されている と考えていい．