%表題   MATH2 数学処理上位パッケージ
%
%履歴
%\Drireki{
%         94/04/04 酒井敏, 塩谷雅人
%         95/04/21 塩谷雅人
%        }
%
\Dchapterhead
\label{odelib}

\section{概要}
このパッケージは連立常微分方程式の積分を行うためのものである.
ここに含まれるルーチンは大きく3つのレベルにわかれる.
最下層のレベルは,
特定のアルゴリズムに従って1ステップ積分するアルゴリズムルーチンで,
その上に精度チェックをしながら2ステップ積分して適当な刻み幅を求める
ステッパルーチンがあり,
最上位に実際に与えられた積分区間の積分を行うドライバルーチンがある.
ユーザーは通常このドライバルーチンを使って積分を行う.
この際, ユーザーは使用するアルゴリズムルーチン, または,
ステッパルーチンを目的に応じて指定することができる.
\vspace{1em}

なお, このライブラリは Numerical Recipes のアルゴリズムに
基づくものである.

\section{サブルーチンのリスト}

  アルゴリズムルーチン

  \vspace{1em}
  \begin{tabular}{p{9cm}p{6cm}}
    {\tt CALL ODRKG(N,FCN,T,DT,X,DX,XOUT,WORK)} & ルンゲ-クッタ-ギル. \\
    {\tt CALL ODRK4(N,FCN,T,DT,X,DX,XOUT,WORK)} & 4次精度のルンゲ-クッタ. \\
    {\tt CALL ODRK2(N,FCN,T,DT,X,DX,XOUT,WORK)} & 2次精度のルンゲ-クッタ. \\
    {\tt CALL ODRK1(N,FCN,T,DT,X,DX,XOUT,WORK)} & 1次精度のルンゲ-クッタ.
  \end{tabular}

  \vspace{1em}
  ステッパルーチン

  \vspace{1em}
  \begin{tabular}{p{9cm}p{6cm}}
    {\tt CALL ODRKGR(N,FCN,T,DT,EPSL,X,WORK)} & 可変幅のきざみを求める. \\
    {\tt CALL ODRKGS(N,FCN,T,DT,EPSL,X,WORK)} & 固定幅のきざみを求める. \\
    {\tt CALL ODRK4R(N,FCN,T,DT,EPSL,X,WORK)} & 可変幅のきざみを求める. \\
    {\tt CALL ODRK4S(N,FCN,T,DT,EPSL,X,WORK)} & 固定幅のきざみを求める.
  \end{tabular}

  \vspace{1em}
  ドライバルーチン

  \vspace{1em}
  \begin{tabular}{p{9cm}p{6cm}}
     {\tt CALL ODRKDU(N,ALGOR,FCN,T,TEND,ISTEP,X,WORK)} &
       固定刻幅で積分を行う. \\
     {\tt CALL ODRKDV(N,STEPER,FCN,T,TINT,DT,X,WORK)} &
       要求精度を満たすように積分を行う.
  \end{tabular}

  \vspace{1em}
  内部変数管理ルーチン

  \vspace{1em}
  \begin{tabular}{p{9cm}p{6cm}}
      {\tt CALL OD}{\it p}{\tt GET(CPARA,IPARA)} & 内部変数を参照する. \\
      {\tt CALL OD}{\it p}{\tt SET(CPARA,IPARA)} & 内部変数を設定する.
  \end{tabular}

\section{サブルーチンの説明}

\subsection{ODRKG/ODRK4/ODRK2/ODRK1}

\begin{enumerate}
  \item 機能
      \begin{quote}
      特定のアルゴリズムに従い,
      与えられたステップ幅で1ステップの積分を行う.
      積分アルゴリズムは
      {\tt ODRKG}がルンゲ-クッタ-ギル,
      {\tt ODRK4}が4次精度のルンゲ-クッタ,
      {\tt ODRK2}が2次精度のルンゲ-クッタ,
      {\tt ODRK1}が1次精度のルンゲ-クッタ(Euler)
      である.
      \end{quote}

  \item 呼び出し方法
      \begin{quote}
      {\tt CALL ODRKG(N,FCN,T,DT,X,DX,XOUT,WORK)}\\
      {\tt CALL ODRK4(N,FCN,T,DT,X,DX,XOUT,WORK)}\\
      {\tt CALL ODRK2(N,FCN,T,DT,X,DX,XOUT,WORK)}\\
      {\tt CALL ODRK1(N,FCN,T,DT,X,DX,XOUT,WORK)}
      \end{quote}

  \item パラメーターの説明
      \begin{quote}
      \begin{tabular}{llp{10cm}}
      {\tt N}    & {\tt (I)}    & 被積分変数(方程式)の数. (i) \\
      {\tt FCN}  & 手続き名     & {\tt DX}を計算するサブルーチン名. \\
      {\tt T}    & {\tt (R)}    & 独立変数$t$の値. (i) \\
      {\tt DT}   & {\tt (R)}    & 積分ステップ幅. (i) \\
      {\tt X}    & {\tt R(N)}   & 被積分変数の$t={\tt T}$における値. (i) \\
      {\tt DX}   & {\tt R(N)}   & 被積分変数の$t={\tt T}$における微分値. (i)
                                  {\tt ODRKG}の時は{\tt DX}の値は
                                  保存されない. \\
      {\tt XOUT} & {\tt R(N)}   & 被積分変数の$t={\tt T+DT}$における値. (i)
                                  実引数として{\tt X}と同じ変数を
                                  指定してもよい. \\
      {\tt WORK} & {\tt R(N,M)} & 作業変数.
                                  M=3 ({\tt ODRK4}) または 1
                                  (その他のルーチン).
      \end{tabular}
      \end{quote}

  \item 備考
      \begin{enumerate}
      \item {\tt FCN} は サブルーチンの形でユーザーが用意する.
      その形式は以下のとおり.
      \begin{quote}
      {\tt SUBROUTINE FCN(N,T,X,DX)}
      \end{quote}

      \end{enumerate}

\end{enumerate}

%=======================================================================
%\section{ステッパルーチン}

\subsection{ODRKGR/ODRKGS/ODRK4R/ODRK4S}

\begin{enumerate}
  \item 機能
      \begin{quote}

      {\tt ODRKG} ({\tt ODRKGR}, {\tt ODRKGS}) または,
      {\tt ODRK4} ({\tt ODRK4R}, {\tt ODRK4S}) により
      2ステップの積分を行ない, 精度のチェックをして適当な{\tt DT}を求める.
      {\tt ODRKGS}と{\tt ODRK4S}は与えられた{\tt DT}で2ステップの積分をした後,
      次のステップのための{\tt DT}を求めてリターンするが,
      {\tt ODRKGR}と{\tt ODRK4R}は要求精度が満たされないとき{\tt DT}を
      小さくして再計算を行う.

      \end{quote}

  \item 呼び出し方法
      \begin{quote}
      {\tt CALL ODRKGR(N,FCN,T,DT,EPSL,X,WORK)}\\
      {\tt CALL ODRKGS(N,FCN,T,DT,EPSL,X,WORK)}\\
      {\tt CALL ODRK4R(N,FCN,T,DT,EPSL,X,WORK)}\\
      {\tt CALL ODRK4S(N,FCN,T,DT,EPSL,X,WORK)}
      \end{quote}

  \item パラメーターの説明
      \begin{quote}
      \begin{tabular}{llp{10cm}}
      {\tt N}    & {\tt (I)}    & 被積分変数(方程式)の数. (i) \\
      {\tt FCN}  & 手続き名     & {\tt DX}を計算するサブルーチン名. \\
      {\tt T}    & {\tt (R)}    & 独立変数$t$の値.
                                  2ステップ分の積分幅を加えて
                                  出力される. (i/o) \\
      {\tt DT}   & {\tt (R)}    & 積分ステップ幅.
                                  適当な値に変更されて出力される. (i/o) \\
      {\tt EPSL} & {\tt (R)}    & 要求精度. (i) \\
      {\tt X}    & {\tt R(N)}   & 被積分変数の$t={\tt T}$における値を入力し,
                                  $t={\tt T+ 2 \times DT}$ における値
                                  を出力する. \\
      {\tt WORK} & {\tt R(N,M)} & 作業変数.
                                  M=7 ({\tt ODRK4R}),
                                  M=5 ({\tt ODRK4S}, {\tt ODRKGR}),
                                  M=3 ({\tt ODRKGS}).
      \end{tabular}
      \end{quote}

  \item 備考
      \begin{enumerate}
        \item {\tt ODRKGS}, {\tt ODRK4S} は
        {\tt ODRKGR}, {\tt ODRK4R} に比べて, 作業領域が小さくてすむが,
        要求精度が満たされたくても, 再計算をしないので,
        記憶領域が許す限り, {\tt ODRKGR}, {\tt ODRK4R} を推奨する.

        なお {\tt ODRKGS}, {\tt ODRK4S} では,
        {\tt DT} を決める際に, 小さめの安全率を
        使っているので, 方程式の性質が時間的に緩やかに変わる場合には
        有効である.
        \item {\tt FCN} は サブルーチンの形でユーザーが用意する.
        その形式は以下のとおり.
        \begin{quote}
        {\tt SUBROUTINE FCN(N,T,X,DX)}
        \end{quote}
      \end{enumerate}

\end{enumerate}

%=======================================================================
%\section{ドライバルーチン}

\subsection{ODRKDU}

\begin{enumerate}
  \item 機能
      \begin{quote}
      各アルゴリズムルーチンにより, 固定刻幅で積分を行う.
      \end{quote}

  \item 呼び出し方法
      \begin{quote}
      {\tt CALL ODRKDU(N,ALGOR,FCN,T,TEND,ISTEP,X,WORK)}
      \end{quote}

  \item パラメーターの説明
      \begin{quote}
      \begin{tabular}{llp{10cm}}
      {\tt N}     & {\tt (I)}    & 被積分変数(方程式)の数. (i) \\
      {\tt ALGOR} & 手続き名     & 使用するアルゴリズムルーチン名 (i) \\
      {\tt FCN}   & 手続き名     & {\tt DX}を計算するサブルーチン名 (i) \\
      {\tt T}     & {\tt (R)}    & 積分を始める独立変数$t$の値.
                       RETURN 時には最後のステップ$t$ {\tt (TEND) }
                       が出力される. (i/o) \\
      {\tt TEND}  & {\tt (R)}    & 積分を終了する独立変数$t$の値. (i) \\
      {\tt ISTEP} & {\tt (I)}    & ステップ数. (i) \\
      {\tt X}     & {\tt R(N)}   & 被積分変数の$t={\tt T}$における値を入力し,
                                   $t={\tt TEND}$ における値を出力する.
                                   (i/o) \\
      {\tt WORK}  & {\tt R(N,M)} & 作業変数.
                                  M=5 ({\tt ODRK4}), M=3 (その他).
      \end{tabular}
      \end{quote}

  \item 備考
      \begin{enumerate}
        \item 使用するアルゴリズムルーチンは, 特に支障がない限り,
        {\tt ODRKG} (Runge-Kutta-Gill) を推奨する.
        \item {\tt FCN} は サブルーチンの形でユーザーが用意する.
        その形式は以下のとおり.
        \begin{quote}
        {\tt SUBROUTINE FCN(N,T,X,DX)}
        \end{quote}
      \end{enumerate}

\end{enumerate}

%-----------------------------------------------------------------------
\subsection{ODRKDV}

\begin{enumerate}
  \item 機能
      \begin{quote}
      各ステッパルーチンにより, 要求精度を満たすように積分を行う.
      \end{quote}

  \item 呼び出し方法
      \begin{quote}
      {\tt CALL ODRKDV(N,STEPER,FCN,T,TINT,DT,X,WORK)}
      \end{quote}

  \item パラメーターの説明
      \begin{quote}
      \begin{tabular}{llp{10cm}}
      {\tt N}     & {\tt (I)}    & 被積分変数(方程式)の数. (i) \\
      {\tt STEPER}& 手続き名     & 使用するスッテパルーチン名. \\
      {\tt FCN}   & 手続き名     & {\tt DX}を計算するサブルーチン名 (i) \\
      {\tt T}     & {\tt (R)}    & 積分を始める独立変数$t$の値.
                       RETURN 時には最後のステップ$t$ {\tt (TEND) }
                       が出力される. (i/o) \\
      {\tt TEND}  & {\tt (R)}    & 積分を終了する$t$の値. (i) \\
      {\tt DT}    & {\tt (I)}    & ステップ幅の初期値. (i/o) \\
      {\tt X}     & {\tt R(N)}   & 被積分変数の$t={\tt T}$における値を入力し,
                                   $t={\tt TEND}$ における値を出力する.
                                   (i/o) \\
      {\tt WORK}  & {\tt R(N,M)} & 作業変数.
                                  M=7 ({\tt ODRK4}), M=5 (その他).
      \end{tabular}
      \end{quote}

  \item 備考
      \begin{enumerate}
        \item 精度は
        {\tt OD}{\it p}{\tt GET/}{\tt OD}{\it p}{\tt SET} の
        管理する内部変数 {\tt 'EPSILON'}によって指定できる.
        最初の呼び出し時に, 引数 {\tt DT} に適当な値を代入しておくと
        要求精度を満たすように {\tt DT} が設定されて返って来るので,
        これを, そのまま次回の呼び出し時に指定すればよい.
        \item {\tt ODRKDV} は
        {\tt OD}{\it p}{\tt GET/}{\tt OD}{\it p}{\tt SET} の管理する内部変数
        {\tt 'MAXSTEP'} で指定された
        回数だけ積分しても {\tt TEND} に到達しないときは,
        エラーメッセージを出力して, 実行を中止する.
        {\tt 'MAXSTEP'} は {\tt ODISET} で変更可能である.
        \item 使用するステッパルーチンは, 特に支障がない限り,
        {\tt ODRKGR} (Runge-Kutta-Gill 再計算) を推奨する.
        \item {\tt FCN} は サブルーチンの形でユーザーが用意する.
        その形式は以下のとおり.
        \begin{quote}
        {\tt SUBROUTINE FCN(N,T,X,DX)}
        \end{quote}
      \end{enumerate}

\end{enumerate}

%-----------------------------------------------------------------------
\subsection{OD{\it p}GET/OD{\it p}SET}

\begin{enumerate}
  \item 機能
      \begin{quote}
      {\tt ODRKDV} で使用する内部変数を参照/設定する.
      \end{quote}

  \item 呼び出し方法
      \begin{quote}
      {\tt CALL OD}{\it p}{\tt GET(CPARA,IPARA)} \\
      {\tt CALL OD}{\it p}{\tt SET(CPARA,IPARA)}
      \end{quote}

  \item パラメーターの説明
      \begin{quote}
      \begin{tabular}{llp{10cm}}
      {\tt CPARA} & {\tt (C*8)} & 内部変数の名前. (i)\\
      {\tt IPARA} & {\tt (I,R)} & 内部変数の値.
      \end{tabular}

      \vspace{1em}
      以下に{\tt CP}として指定できる名前のリストを記す.

      \begin{tabular}{llp{10cm}}
      {\tt 'EPSILON'} & {\tt (R)} & 要求精度. \\
      {\tt 'MAXSTEP'} & {\tt (I)} & 最大積分回数. \\
      {\tt 'NSTEP'}   & {\tt (I)} & 実際に積分した回数.
      \end{tabular}
      \end{quote}

  \item 備考
      \begin{enumerate}
        \item {\tt IPARA}としては適切な型の定数または変数を
          指定すること.
      \end{enumerate}

\end{enumerate}



