ロジスティック方程式の差分化

カオス力学系の話は詳しいサイトがたくさんあるのでそちらをどうぞ。

ロジスティック方程式のような非線形の微分方程式を不用意に差分化すると
本来の解とは全く異なる奇妙な挙動を示すことがある。
初めて目にしたときはとても面白いと思ったが、同時にこうも思ったものだ。
「差分化のしかたが悪かっただけじゃないの？」と。
気になったので自分でやってみることにした。

今回使用するロジスティック方程式は
$\begin{equation} \begin{aligned} \begin{array}{ll} \frac{dy}{dt} &= a y \left( 1-y\right) \\ &= a \left( y - y^{2} \right) \qquad (1) \end{array} \end{aligned} \end{equation}$
これを正しく差分化するのが目的になる。

まず差分で使う演算子について考えよう。
差分演算子 $\delta$ と平均演算子 $\nabla$ を以下のように定義？する。
$\begin{equation} \begin{aligned} \begin{array}{ll} \delta f(x) &= f(\nabla x + \frac{\delta x}{2}) - f(\nabla x - \frac{\delta x}{2}) \\ \nabla f(x) &= \frac{ f(\nabla x + \frac{\delta x}{2}) + f(\nabla x - \frac{\delta x}{2}) }{2} \end{array} \end{aligned} \end{equation}$
ここで、 $x$ は任意の変数、 $f(x)$ は $x$ の任意の関数である。

残念ながら、循環しているためこれだけでは定義になっていない。
そこで以下の性質をもつ変数 $n$ の存在を仮定する。
$\begin{equation} \begin{aligned} \begin{array}{ll} \delta n &= 1 \\ \nabla n &= n \end{array} \end{aligned} \end{equation}$
そして変数 $n$ をインデックスとよぶことにする。

インデックスを変数にとれば
$\begin{equation} \begin{aligned} \begin{array}{ll} \delta f(n) &= f(n + \frac{1}{2}) - f(n - \frac{1}{2}) \\ \nabla f(n) &= \frac{ f(n + \frac{1}{2}) + f(n - \frac{1}{2}) }{2} \end{array} \end{aligned} \end{equation}$
となって、中央差分になっていることが分かる。

なぜ回りくどい定義にしたか。
$y = y(x), x = x(t)$ を考えよう。
$\begin{equation} \begin{aligned} \begin{array}{ll} \delta y &= y(\nabla x + \frac{\delta x}{2}) - y(\nabla x - \frac{\delta x}{2}) \\ &= y( \frac{ x(\nabla t + \frac{\delta t}{2}) + x(\nabla t - \frac{\delta t}{2})}{2} + \frac{ x( \nabla t + \frac{\delta t}{2}) - x(\nabla t - \frac{\delta t}{2})}{2} ) \\ &\quad - y( \frac{ x(\nabla t + \frac{\delta t}{2}) + x(\nabla t - \frac{\delta t}{2})}{2} - \frac{ x( \nabla t + \frac{\delta t}{2}) - x(\nabla t - \frac{\delta t}{2})}{2} ) \\ &= y( x(\nabla t + \frac{\delta t}{2}) ) - y( x(\nabla t - \frac{\delta t}{2}) ) \end{array} \end{aligned} \end{equation}$
となって、合成関数 $y(x(t))$ の差分演算子をとったものと一致する。
つまり気軽に変数変換ができるということで、等間隔差分にはない利点である。

ちなみに差分商に関するチェーンルール
$\begin{equation} \begin{aligned} \frac{\delta y}{\delta t} = \left. \frac{\delta y}{\delta x} \right|_{x=x(t)} \frac{\delta x}{\delta t} \end{aligned} \end{equation}$
も成り立つが、ただの通分である。

基本的な性質。
$\begin{equation} \begin{aligned} \begin{array}{ll} \delta a &= 0 \qquad （aは定数）\\ \nabla a &= a \qquad （aは定数）\\ \delta \left( yz \right) &= \delta y \nabla z + \nabla y \delta z \\ \nabla \left( yz \right) &= \nabla y \nabla z + \frac{1}{4} \delta y \delta z \\ \delta \frac{1}{y} &= - \frac{\delta y}{y^{(2)}} \\ \nabla \frac{1}{y} &= \frac{\nabla y}{y^{(2)}} \\ \delta \log y &= 2 \tanh^{-1} \frac{\delta y}{2 \nabla y} \\ \nabla \log y &= \frac{1}{2} \log y^{(2)} \\ \end{array} \end{aligned} \end{equation}$
ここで
$\begin{equation} \begin{aligned} y^{(2)} = (\nabla y + \frac{\delta y}{2})(\nabla y - \frac{\delta y}{2}) \end{aligned} \end{equation}$
である。
証明は四則演算だけで簡単なのでしない。

どう差分化するか。
$(1)$ 式は解けることが分かっている。
正しく差分化するということは、 $(1)$ 式の解を再現しなければならないから
差分方程式として解ける必要があるだろう。

まず $(1)$ 式を解いてみる。ベルヌーイ型だから
$\begin{equation} \begin{aligned} -\frac{1}{y^{2}} \frac{dy}{dt} = -a \left( \frac{1}{y} - 1 \right) \end{aligned} \end{equation}$
と変形する。
ここで $z = \frac{1}{y} - 1$ とおくと $\frac{dz}{dt} = -\frac{1}{y^{2}} \frac{dy}{dt}$ だから
$\begin{equation} \begin{aligned} \frac{dz}{dt} = -a z \end{aligned} \end{equation}$
$\begin{equation} \begin{aligned} \frac{1}{z} \frac{dz}{dt} = -a \end{aligned} \end{equation}$
$\begin{equation} \begin{aligned} \frac{d \log z}{dt} = -a \end{aligned} \end{equation}$
$\begin{equation} \begin{aligned} \log z = -at +const. \qquad(2) \end{aligned} \end{equation}$
$\begin{equation} \begin{aligned} z = A e^{-at} \end{aligned} \end{equation}$
$\begin{equation} \begin{aligned} \frac{1}{y} -1 = A e^{-at} \end{aligned} \end{equation}$
$\begin{equation} \begin{aligned} y = \frac{1}{1+A e^{-at}} \end{aligned} \end{equation}$
となる。

上の解法を参考にして差分化してみる。 $(2)$ 式から
$\begin{equation} \begin{aligned} \delta \log z = -a \delta t \end{aligned} \end{equation}$
$\begin{equation} \begin{aligned} 2 \tanh^{-1} \frac{\delta z}{2 \nabla z} = -a \delta t \end{aligned} \end{equation}$
$\begin{equation} \begin{aligned} \delta z = - 2 \tanh \frac{a \delta t}{2} \nabla z \qquad (3) \end{aligned} \end{equation}$
ただし $\tanh \left( -x \right) = - \tanh x$ を使った。
$z = \frac{1}{y} - 1$ だったから。
$\begin{equation} \begin{aligned} \delta z = - \frac{\delta y}{y^{(2)}} \end{aligned} \end{equation}$
$\begin{equation} \begin{aligned} \nabla z = \frac{\nabla y}{y^{(2)}} - 1 \end{aligned} \end{equation}$
これらを $(3)$ 式に代入すると
$\begin{equation} \begin{aligned} - \frac{\delta y}{y^{(2)}} = - 2 \tanh \frac{a \delta t}{2} \left( \frac{\nabla y}{y^{(2)}} - 1\right) \end{aligned} \end{equation}$
整理すると。
$\begin{equation} \begin{aligned} \frac{\delta y}{\delta t} = a^{*} \left( \nabla y - y^{(2)} \right) \qquad (4) \end{aligned} \end{equation}$
となる。ここで
$\begin{equation} \begin{aligned} a^{*} = \frac{2}{\delta t} \tanh \frac{a \delta t}{2} \end{aligned} \end{equation}$
とおいた。

$\delta t \to 0$ のとき
$\begin{equation} \begin{aligned} \frac{\delta y}{\delta t} \to \frac{dy}{dt} \end{aligned} \end{equation}$
$\begin{equation} \begin{aligned} \nabla y \to y \end{aligned} \end{equation}$
$\begin{equation} \begin{aligned} y^{(2)} \to y^{2} \end{aligned} \end{equation}$
$\begin{equation} \begin{aligned} a^{*} = \frac{2}{\delta t} \tanh \frac{a \delta t}{2} \to a \end{aligned} \end{equation}$
であるから、 $(4)$ 式は元のロジスティック方程式に戻る。

実際に $(4)$ 式を解くには $(2)$ 式に向かって変形していけばよいだけである。
ただし
$\delta x = 0 \Rightarrow x = const.$
はそれほど自明ではない。
いま $t = kn$ （ $k$ は定数）、 $x = x(t)$ として
$\delta x = x( k(n + 1/2))-x( k(n - 1/2 )) = 0$
を仮定する。
インデックスをシフトして移項すれば
$x( k(n + 1) )= x( kn )$
$x( kn + k )= x( kn )$
$x( t + k )= x( t )$
ここで $n=0$ とおけば $t=0$ であるから
$x( k )= x( 0 ) \qquad (5)$
$k$ をいろいろ変えたすべての変数 $t$ について $(5)$ 式が成り立つから $x = const.$ である。