ここでは、生物の移動を考慮したときに避けて通れない、反応拡散方程式モデルについてその数値解法の基礎を解説する.
サンプルコードはXXX.
[1] モデルの導出(拡散項の導出)
[2] 拡散方程式の離散化アルゴリズム
[3] 離散化アルゴリズムの安定性
[4] 拡散方程式の初期条件と境界条件
[5] 簡単な反応拡散方程式(指数成長+固定境界条件:プランクトンの最小パッチ問題)
空間構造のある個体群動態をモデル化するには,局所的な個体群動態(空間について考える必要ない)と,質量保存則に基づく連続方程式とランダムな個体の移動,それぞれについて独立して定式化し,最後に組み合わせればよい.
局所的な個体群動態の部分は非常に簡単で,常微分方程式を使って出生と死亡を定式化すればよいだけ.
連続方程式については,簡単のため以下のような一次元空間を考える.河川のような横に伸びる空間を想定(川も実際には三次元であるが)し,川の上流が左側,下流が右側だとして,位置x・時刻tにおける生物量をN(x, t)と置く.連続方程式と拡散項を導くには別に生物ではなくてもよいことにも注意.なぜなら出生・死亡という生物らしさのある過程はすでに「局所的な個体群動態」として分離済みであるからである.つまり連続方程式(および拡散項)を考えるときは、まったく反応しないカラーインクの分子みたいのが川に流れていると考えてた方がよい(川ではあるが媒体である水は動いていないことを想定する.そうしないと「移流項」という別の過程も考えないといけなくなる).次図のような感じで定式化できる.
x軸が空間を表している.x = x1からx = x1 + Δxまでの領域(薄い青色になっている部分)に注目し,その領域内の全個体数(全分子数)をカウントする状況を考えよう.しかもこのカウントを短い間隔の二つの時刻t=t1とt=t1+Δtで行い,その差を取ることを考える.つまり短い時間間隔Δtの間にどれだけ個体数が変化するかを定量化すると思えばよい.この個体数変化の計算の仕方は二通りある.一つ目の計算は上の式の左辺に該当しもう一つは右辺に該当する.
まず,左辺の計算から始めよう.この領域に時刻t1に存在する全個体数は、青い部分の面積に対応するので、厳密にはx = x1からx = x1 + Δxまでの領域で積分すればよい.しかし後程Δxがゼロに近づく極限を考えるのでΔxの二乗以上の項(高次の無限小)は無視してよく,幅Δx,高さN(x1, t1)の長方形で近似してよい.同様に時刻t1+Δtに存在する全個体数も幅Δx,高さN(x1, t1+Δt)の長方形の面積として近似できる.この二つの時刻での差を取れば,時間間隔Δtの間にどれだけ個体数が変化するかの計算となる.
一方,もう一つの計算(右辺)には質量保存則 (mass conservation law, or mass balance)を用いる.今考えているのが一次元空間であることを思い出そう.かつ出生・死亡という生態学的過程はすでに切り離されているので考える必要がない.つまり新しい個体が生じることも現存する個体が消滅することも考えなくてよい.だとするとこの時間間隔Δtの間に個体数に変化があるとすればそれは個体の移動によるものしかありえない.もっといえば、考えている領域の両端の壁を行き来した個体だけがこの個体数変化に影響を及ぼすことになる.そこでそのような個体の移動量、特に単位時間当たりの移動量(=流量)を右向き(xが大きくなる方向)を正とし,かつ流量は地点xに依存するとしてJ(x)と定義する.すると時間間隔Δtの間に左側の壁を超えて流入した個体数の正味の数(=流入した個体と流出した個体の数の差)はJ1Δt = J(x1)Δtと計算できる.同様に右側の壁を越えて流出した個体数の正味の数はJ2Δt = J(x1+Δx)Δtと計算できる.さあ,この流入数と流出数の差をとれば個体数の変化が計算できる.
残りはテクニカルな計算のみであり,J(x1+Δx)をx=x1のまわりでTaylor展開すればJのxに関する偏微分係数を係数とするΔxの一次項と高次項がでてくる.Δtで両辺を割りΔt→0の極限をとればNのtに関する偏微分係数がでてくるし,右辺・左辺のΔxの一次項は両辺をΔxで割ればO(1)のオーダーになって残るのに対して右辺に残るO(Δx)はΔx→0の極限でちゃんと消える.このような計算で連続方程式(Equation of continuity)を導出できるわけである.
続いて拡散項はランダムな(方向性の無い)個体の動きから導出することができる。この導出過程は高校や大学基礎レベルの統計力学的アプローチを用いる(理想気体の気体方程式 pV = nRTの導出と同じ方法)。導出の流れは以下のような図にまとめることができる。
l(エル)・バーは、「平均自由行程」を表し、個体間(粒子間)の衝突およびそれによる個体間相互作用によって個体の移動方向が変わってしまわないで個体が一定方向に動くことができる距離の期待値である。同様にΔt・バーは個体が別個体との衝突をせずに一定方向に動くことができる時間の期待値です。最後のv・バーは一定方向に移動している時の速度の期待値です。この3つの変数のうち、二つが独立で残りの一つはそれら二つによって決まるということになります。
さて、ごく短い期間(Δt・バー)での生物個体の移動を定量化する設定として、上の左図x軸上x = x1に観測者がいるとします。この境界x=x1では、境界を左側から右側に通過する個体と右側から左側に個体がいますが、それらを差し引きして左側から右側へ通過する純個体数を正の流れとします。特に単位時間当たりの流れをJ(x1)と定義すると、今注目している期間で移動する純個体数はJ(x1)*Δt・バーとなります。