数値計算入門 武尾英哉. 離散数学と数値計算 数学的解法の中には理論計算では求められないものもある. 例えば, 定積分は, まずは積分 ( 被積分関数の原始関数をみつけること できなければ値を得ることはできない. また, ある関数の所定の値における微分値を得るには, まずその関数の微分ができなければならない. さらに代数方程式の解を得るためには, 解析的に代数方程式を解く必要がある. ところが, これらは必ずしも解析的に導けるとは限らない. そこで本章では, 計算機を使って, 離散的な数値 ( とびとびの数値を扱う を扱って, 近似的に解を導く手法について説明する. このことを離散数学という. 本稿では, 連続的な数値による数学 ( これまでの解析的数学 を, 離散的な数値による数学として近似的に扱うことで解を求めることができる手法について説明する. このように, 離散数学を用いて種々の手法による近似計算により解を求めることを数値計算とよぶ.. 数値代数方程式 ( ニュートン法 数値計算を使って代数方程式を解く代表的な手法であるニュートン法について説明する. 代数方程式を解く ( f( を満足する を求める ということは, グラフにおいて 関数 f( の 軸切片を求めることと等価である. すなわち, 図 における を 求めることである. 以下, ニュートン法による手順について説明する. まず初めに, 予想される真の解に近いと思われる値をひとつとる ( の初期値. の初期値を とする. 次に, グラフの接線を考え, その 軸切片を計算する. 図 において, にお f '( ( f( であるので, この接線が 軸と交差す ける接線の方程式は る点は, が f ( f '( のときである. この を新たに f ( f '( ただし, f '( とする. とする. この 軸切片の値は, 予想される真の解により近いものとなるのが一般である. 3の手法を繰り返すと, 限りなく に近い値となる. すなわち, 一般的に n
から次のステップの n を求める式は次のようになる. n f( n n ただし, f '( n とする. ( f '( n 4この繰り返し計算を希望する精度の解を得るまで行う. すなわち, 実際には, 次の収束判定条件を計算し, 満足するまで繰り返すことで, 解の近似解を得る. n n または n n ( n ここで,, は収束条件である ( 例えば,. など微小な値を設定する. なお, 収束条件の詳細に関しては, 数値計算や数値解析等の関連書籍を参照されたい. ニュートン法は, 一般的に早く収束するが, 初期値 のとり方によっては収束しない こともある. また, 解が複数ある場合には, 初期値の設定により別の解に収束する場合がある. よって, 予めおおまかな解を初期値としてうまく設定する必要がある. f '( ( f( f( 図 : ニュートン法の原理 [ 例 ] を計算で求める場合として, f '( であるので, 式 ( は n n n n と書き表せる. ここで, 初期値を例えば f おくと, n は に収束していく. ( の解をニュートン法で求める. とおくと, n は に収束し, 一方, と
3 3. 数値微分 ( 差分法 数値計算を使って関数の所定値の微分値を求める代表的な手法である差分法について説明する. 数値計算では無限小を扱うことができないので, 小さいが有限の値を用いた微分演算の近似である差分を用いる. ここで, 微分はその点での傾き f ( f( f '( lm (3 で定義されるが, 差分では有限距離離れた 点の傾き f ( f( (4 である. 差分は下記のように, 前進差分, 後退差分, 中央差分の 3 つがある ( 図 参照. f( の における微分係数は, 微小区間 について近似的に下記のように求 めることができる. ( 前進差分 f f f( f( f '( ( 後退差分 f f f ( f( f '( (3 中央差分 f '( 3 (5 (6 f f f( f( 3 (7 f( f f f 3 図 : 数値微分の原理
4 4. 数値積分 ( 台形法とシンプソン法 定積分の数値計算法としては, を n 分割し各区間の両端の f ( を直線で結んで台形 の面積の総和を求める方法 ( 台形法 と, 区間を 次式で近似して面積の総和を求める方法 ( シンプソン法 などがある. 解析的に積分が不可能な場合の定積分を求めるときに有効であり, 数値積分とよばれる. ( 台形法図 3 に示すように, 区間 [,b] を n 等分して, k, k における関数 f ( の値 f (, f( を利用した微小台形を作成する. この台形の面積は, 下式で計算され, k k 微小領域の面積を近似する. f ( k f( k b f( k f( k ( k k (8 n f f は台形の上底と下底を表している. ここで, k は台形の高さ k, ( k, ( k 一方, 微小領域の面積 S は定積分を用いて下記のように表せるので, 式 (8 と併せて次式が成り立つ. k b f( k f( k S f( d (9 k n よって, 区間 [,b] の面積 S は, n b S S f( f( f( f( f( f( n f( n f( n n k b f ( f ( f ( f ( n f ( n ( n b となり, 定積分 f ( d は次式で近似できる. これを台形法とよぶ. b n b f ( d f( f( f( f( n f( n (
5 f( f( f ( k f ( k k n等分図 3: 台形法の原理 k k n b k k n等分 n 図 4: シンプソン法の原理 b ( シンプソン法図 4 に示すように区間 [,b] を n 等分して, 微小区間 k k で関数 f ( を, f(,, f(,, f( k k k k k k とを考える. 次曲線を gt ( p qt rt 記のように表せる. の 3 点を通る 次曲線で補間するこ と定義し, 微小区間 k k [, ] における面積 S は下 k 3 S f ( d g( t dt p qt rt dt pt qt rt k 3 p q r p q r p r 6p r 3 3 3 3 3 3 3 ( ここで, f ( g( pq r k f ( g( p k f ( g( pq r から, k 6p r f( k 4 f( k f( k (3 式 (3 を式 ( に代入して, g( t dt f( k 4 f( k f( k 3
6 b であるので, n k b f( k 4 f( k f( k S f( d g( t dt k n 3 よって, b n n b S f( d S f( k 4 f( k f( k 6n k k b f( 4 f ( f ( f ( 4 f ( 3 f ( 4 6n b { f ( 4 f ( f ( 3 f ( n 6n f( f( f( f( } (4 4 n n となる. これをシンプソン法とよぶ. ただし, シンプソン法では区間を n 等分するので, 偶数分割することが前提となる. 5. 数値微分方程式 ( オイラー法 本節では, 微分方程式の数値解析の一つであるオイラー法について説明する. 次式の微分方程式を数値解析によって解くこととする. d f (, d ただし, 初期値 のとき とする ( ( (5 微分の定義から, d ( ( lm (6 d これが f (, と等しい. 一方, が十分に小さいとき, 式 (6 の右辺は ( ( に近似することから, ( ( f(, (7 すなわち, ( ( f(, (8 が成り立つ.
7 初期値を, ( とすれば,, ( f(,,, ( f(,, と順 次に ( を計算できる. 微分方程式の数値解析とは, これらの点 (, の集まり を求めることであり, さらに平面上に (, をプロットすると求める関数 のグラフが 描ける. このようにして, 微分方程式の解 ( 関数 がグラフ上に得られることになる. オイラー法をまとまると, 以下の手順になる. 初期値を与え, 漸化式 ( (, 初期値, ( を決める. k k f k k を順に計算すればよいので, ( ( ( f(, f(, を計算する. ( ( ( f(, f(, を計算する. 3 ( ( ( f(, f(, を計算する. 4 3 3 f f を計算す 5 順次, k ( k ( k ( k ( k, k k ( k, k る. d [ 例 ] 微分方程式 d, 初期値を,,.5 とするとき, オイラー法により微分方程式の解を求める. f (,, 前記オイラー法の手順より, (.5 (.5 f(,.5 f(,.5.5 (.5.5 ( f(,.5 f(,.5.5.5.5 (.5 (.5 f(,.5 f(,.5.5.5 3.375 3 となり, (, (.5.5, (.5, (.5 3.375, となる. 一方, 与式の微分方程式の解析解は, ( eより, (, (.5.649, (.78, (.5 4.48, となり, 誤差は図 5 のように変化する.
8 また, 刻み幅 によっても解析解との誤差に大きく影響を与える. しかし, 刻み幅 を細かくすることは計算時間を増加することになるので, 設定には留意する必要がある. 5 4 3 (, 数値解 (, 解析解 (, 3 3.5.5 図 5: オイラー法の誤差の変化 さらに, オイラー法は, 区間 [ k, k ] d における傾き ( 微分係数 d ( k ( k を一定 ( 直線 と仮定して計算するので, ( の微分係数が大きい場合にも誤差が大き くなりやすい. 本節で説明したオイラー法は一般的には前進オイラー法とよばれる. 他にも後退オイラー法や前述のような誤差を小さくするために改良された修正オイラー法や改良オイラー法などが存在する. また, 非常に精度が高く実用的にも最も広く使われているルンゲ クッタ法などもある. これら手法の詳細に関しては, 数値計算, 数値解析, 回路シミュレーション等の関連書籍を参照されたい.