IGOR Pro にはいろんな解析機能とそのサンプルが付属していますが、今回は微分方程式に挑戦してみます。今の教育指導要領では高校の数学で扱われなくなったそうですが、微分方程式とは、その名の通り方程式に微分が入っている方程式のことですね（説明になっていないなぁ）。
▼ [おさらい] 変数分離による微分方程式の解法
すっかりにぶってしまっている記憶に刺激を与えるつもりで、ちょっと教科書を紐解いてみました。
dy/dx = g(y)f(x)
という方程式を考えます。この方程式は、左辺に y を x で微分した項が含まれているので「微分方程式」の仲間ですね。g(y) ≠ 0 であれば両辺を g(y) で割って、両辺に dx を掛けると
dy/g(y) = f(x)dx
となり、x を右辺に、y を左辺に分離できます。あとは、初期条件と g(y) と f(x) の複雑さ（面倒さ）次第ですが、両辺を積分するときっと解ける...気がしますよね。
例えば、次のような簡単な例を考えます。
dy/dx = -ky
変形し、積分定数を C0 として積分すると、
dy/y= -kdx
ln y = -kx + C0
∴ y = e(-kx+C0)
となります。ec0 は定数ですので改めて C とおくと、
y = C e－kx
と、解くことができます。解けると、気持ちイイです。数学も、このくらいまでなら好きでいられる気がします。
▼ でも、数値積分
こうやって、解析解がいつでも簡単に求まればよいのですが、そうは行かないのが世の常。そういう時には、初期条件からどんどん計算していって、微分方程式が表現している現象を確認する方法、数値解析（数値積分）が利用されます。化学反応の速度について次の式を考えてみましょう。
dC/dt = -kC
......上の例と同じ式ですね。解析的に解けるのですが、簡単なので数値解析についてもこの例を使います。化学反応としては、よく混合された反応槽で、物質 C が大過剰に存在する他の物質と反応して減少していく1次反応のケースに相当しますね。もちろん、C は濃度、k は反応速度定数です。
さて、IGOR Pro で微分方程式の数値解を得るには IntegrateODE という操作関数を使います（やっと IGOR Pro の話になりました）。
IntegrateODE derivFunc, cwaveName, ywaveSpec
derivFunc は導関数（ユーザ定義関数）、cwaveName は定数パラメータを与えるウェーブ。ywaveSpec は計算結果を受け取るウェーブで、初期推定値はこのウェーブに与えます。
導関数はユーザ定義関数として、都度定義してあげる必要があります。今回の例では、次のような関数の定義が必要です。
Function FirstOrder(pw, xx, yw, dydx)
Wave pw //パラメータ（ウェーブ）
Wave yw //yデータ（ウェーブ）
Variable xx //xデータ（変数）
Wave dydx //導関数（ウェーブ）
dydx[0] = -pw[0]*yw[0]
End
pw はパラーメータのウェーブで pw[0] には k が含まれます。xx は導関数を計算する x を与える変数ですが、この例では関数に x の項がありませんので使用されません。ただし、導関数としては定義に含める必要があります。 yw は y データ（濃度データ）、dydx は導関数です。
この関数は、データレンジ全体に対して実行されるのではなく、1点ずつ実行されます。今回の微分方程式は式が1つだけなので、実際にはyw[0] と dydx[0] だけが利用されるわけです。（導関数の定義方法については、また別の機会でもご紹介しますので、今回はここまでで。）
さて、定義された関数で数値解を得るには、例えば次のようなコマンドを実行します。
Make/D/O/N=101 YY // 結果を受け取るウェーブを作成
YY[0] = 1 // 初期値を設定
Display YY // グラフを作成
Make/D/O PP = {0.05} // パラメータ k を指定
IntegrateODE FirstOrder, PP, YY // 実行！
最初に、結果を受け取るウェーブ（YY）を作成し、次に初期値を YY[0] に与え、グラフに表示します。PP はパラメータ（k）を保持する 1 ポイントだけのウェーブを作成し、その値として 0.05 を与えます。そして IntegrateODE を実行です。
Fig.1 IntegrateODE の実行結果
きれいな減衰曲線が得られます。k = 0.01 すなわち PP = 0.01 として再度実行するともっと緩やかな減衰、0.1 とするともっとシャープな減衰を確認できます。
Fig2. パラメータ k による変化
▼ まとめ
微分方程式の数値解析を IntegrateODE という操作関数で行いました。例に上げた微分方程式が解析的に解ける系なので、あまりメリットは見当たりませんが、IntegrateODE を使うと、こんなことができるんだということは、お伝えできたかと思います。今回のポイントは、導関数をユーザ定義関数として定義するところです。
近いうちに、もうちょっと複雑な系、連立微分方程式の数値解析を IntegrateODE で行う例を紹介することにします。
井上＠技術開発部でした。