next up previous
Next: 2 次回 Up: システム数理IV 第4回 Previous: システム数理IV 第4回

1 陰的公式

陰的公式については、以下の強力な結果が知られている:

s段陰的公式の到達可能次数は 2s である。」

これのため、いくらでも次数の高い公式をつくることができる。陽的な場合に は8次よりも高次のものはほとんど知られていないのとは対照的である。

これは単なる存在証明ではなく、実際に 2s 次の公式が構築されている。こ れは実はユニークに決まり、「陰的ガウス公式」と呼ばれるものがそれである。 以下、陰的ガウス公式と、それに関係したいくつかの公式についてその導出の 原理と計算法を簡単に説明しよう。

陰的ガウス公式以外にも陰的な公式はいろいろあるが、それらについては後で 触れる。

1.1 ガウス型数値積分

ガウス型公式について述べる前に、まずそれと非常に関係の深い数値積分公式 について説明しておく。

数値積分とは、一変数であれば

と、微分方程式ではなく単に関数として与えられたものに対して、ある固定区 間 での定積分を数値的に求めることである。一般に、解析的に 積分ができないようなものを計算しないといけないことはよくあるので、そう いうものが必要になることは多い。

数値積分法として、すでに皆さんが知っているであろうものは、たとえば

といったものであろう(知らない人は、、、まあ、以下で簡単に説明します)。なお、ロンバーグ加速については後でもう少し詳しく説明する。

例えば台形公式は、等間隔に区間を分けたとして、ある部分区間での積分を

と、両端を結んだ台形の面積で近似するということに相当する。これに対して、 もっとましな方法というのを理論的に考えられないであろうか?

一つの方法は、区間の中の点をいくつかとって、それらを通る補間多項式を作 り、それを積分することである。中の点を等間隔にとったものをニュートン・ コーツ型の公式という。 シンプソン公式はこのやりかたのもっとも簡単な場 合に相当する。

一般に補間するための点をp個とった時、 p-1次の多項式を持って くればそのすべての点を通るように出来る。ニュートン・コーツ型の公式はそ のような最小次数多項式を使うものである。例えばもとの関数が実際にp-1 次以下の多項式ならば、厳密に正しいものが求まる。

この補間の原理を説明しておく。例えば補間に使う点が とあ るとすれば、以下のような多項式を考える

この多項式は、 で 1、それ以外のすべての で0になる。従って、 補間多項式は

とすればいい。

数値積分に使うには、あらかじめ の積分を求めておいて、それと を掛けて和をとるだけである。

等間隔でなく点をとることで精度をあげることはできないかというのが、ここ での問題である。

ここで、一般的に、不等間隔に点をとってなんとかするという観点から考える こともできるが、なかなか難しいので答のほうから考えていく。

答は、「ガウス・ルジャンドルの多項式の零点を点にとると、次数を最大にで きる」というものである。ガウス・ルジャンドルの多項式は、有限区間で帰納 的に定義できる直交関数系である。つまり、最初の項を定数として、k番目 の項は k次多項式であって、それまでの項すべてと直交するようにとるので ある。つまり、関数基底にグラム・シュミットの直 交化を適用したものがガウス・ルジャンドルの多項式である。

1.2 数値積分に関する話の続き

ニュートン・コーツ型の公式の一般形についてはすでに説明したが、あれでは 具体的なことが良くわからないと思うので例によって次数が低いものを構成し てみる。

を区間 (これは、すでに適当に分割したあとの狭い区間である と思って欲しい)で数値的に積分するのに、例えば関数を1度だけ計算すると すれば、以下のようにいろいろなものが考えられる

これらは、それぞれ、その区間で関数の値が定数であると近似して、その定数 を積分したというふうに解釈できる。もちろん、中点以外の任意の点をとるこ とも考えられる。

計算精度という観点からは、中点をとることには非常に特別な性質がある。つ まり、この場合だけ局所誤差が になるのである。これは、テイラー 展開の1次の項が積分すると消えることから容易に示すことができる。という わけでこれには名前がついていて、中点公式という。

さて、それでは、2度計算していいとすればどうであろうか?普通に考えつく のは先週も書いた台形公式である。つまり、関数を

で近似して、それを積分して求まる

を積分の値として使うわけである。

台形公式の局所誤差は であり、中点公式と同じである。2点も使っ て、1点しか使わない中点公式と同じでは意味がないと思う人もいるかもしれ ないが、実はそうでもない。1区間だけを考えるとたしかに台形公式は2点使う が、ある区間の左端はもちろんその左の区間の右端なので、関数の値は同じだ からである。つまり、全体をn個の区間に分けた時、中点公式では n回関 数を計算するが、台形公式ではそれより1多い n+1回計算するだけである。

区間の端ではない、一般的な位置に点をとるとすればどうであろうか?以下、 面倒なので区間 級関数

の積分をするということにする(適当な変数変換をすればいつでもこう置き換 えられる)。

積分した結果は

となって、テイラー展開の偶数項しか現れないことに注意して欲しい。

2点の x座標を としたとき、一般的には数値積分の近 似式が

 

の形になる。これで、上の積分を近似するわけだが、0次の項が消えるために , また奇数次の項が出てこないために かつ である必要があることは すぐにわかるであろう。

というわけでと決まってしまったので、あと 出来ることは の値を決めることだけである。式 11を書き直すと

となるので、

とすれば、までの項を一致させられるということがわかる。

台形公式は に相当し、 までしか消えていないので、ずいぶ ん精度が上がっていることがわかるであろう。 そのかわり、区間の端の点の再利用が出来ないので、中点公式に比べて計算量 が2倍になっているということに注意しておく。

もうちょっと頑張って3点にしてみよう。等間隔に点をとるなら両端と中点 であり、このときは2次曲線で近似して局所誤差は3次となる。これはいうまで もなくシンプソンの公式である。

不等間隔の時何が実現できるかだが、ここでは具体的に形を書かない。2点の 場合と同様な考察から、区間の中点を1つ、両側に1点づつとると、こんどは までの項を一致させられるということがわかる。

ここまでの話だと、

となるような気がする。等間隔の場合はまあいいとして、不等間隔の場合にそ ううまくいくのであろうか?

1.3 ガウス型公式の次数

前回の資料の最後に書いたように、ガウス型(あるいは、ガウス・ルジャンド ル型)積分公式は、ルジャンドル多項式の零点を点にとるものである。ルジャ ンドル多項式は有限区間 で与えられる直交関数系である。そのk番 目の項k次多項式であり以下の性質を持つ。

定義はいろいろあるが、

と書かれることが多い。

の零点とは、の解のことである。n次ルジャンドル多項式の場合 実区間n個必ず零点がある。これを以下 と する。

さて、以下でガウス型公式の次数について考える。n点使った時に2n次 になる、いいかえれば 2n-1次までの多項式については厳密な値を与えると いうことが予想されるので、これを示すという方針で考えていく。

任意の 2n-1次多項式 に対して、の全ての零点で値が 一致するようなn-1次多項式 を構成することができる。これを 使うと、

と書ける。ここで、 n-1次の多項式である。これを区間 で積分すると、もちろん

であるが、右辺の第一項の n-1次多項式であり、n次ルジャ ンドル多項式 n-1次以下の任意の多項式と直交するので、この積分 は0になる。つまり、

となるわけである。というわけで Fの積分を計算するには Lの積分を求めれ ばいい。ところが、 Ln-1次なので、 n 点での値が決まれば一意に 決まり、積分出来てしまう。

結局、 の零点での値を使えば、 の積分が計算できたことに なる。

さらに厳密な話としては、 2n-1次よりも高次の項がある関数を近似して積 分したときにどのような誤差項が入ってくるかという議論が必要になるが、こ れは省略する。

一般に、特異点を持たないような性質のよい関数の時には、ガウス型積分公式 は非常に強力な方法であるといえる。

なお、一般の次数の場合、ガウス型公式は

の形をとる。これらの定数の値については適当な参考書を見ること

1.4 微分方程式への積分公式の適用

さて、ここまで数値積分の話をしたわけだが、微分方程式の数値解法と数値積 分はもちろん非常に良く似たものである。が、大きな違いは、数値積分では被 積分関数が

と、独立変数の陽な関数として与えられていたのに対し、微分方程式では

と積分された結果自体に被積分関数がよっていることである。このために、話 がすこし複雑になる。

例えば、もっとも簡単な中点公式の場合を考えてみよう。数値積分の場合には、 単に

であった。微分方程式の場合、問題は での x の値を我々 はまだ知らないということである。

この困難は、数値積分が近似多項式の積分として与えられたということを思い 出せば一応解決できる。中点公式の場合、数値積分ではこれはfで近似するということに相当した。同様なアイディアを微分 方程式に適用すれば、微分方程式の場合、 fを定数で近似し、その結果解で ある x を一次式で近似する、つまり

という近似をする。これから、問題の中点での値を と書くと

と書けるわけで、これは に対する代数方程式になっている。これを 解けば数値解も求まるわけである。

一般には、n次のルジャンドル多項式から導かれる公式では、fn-1次 近似多項式を構成することになる。これを前節と同様 とし、零点を とすれば、での「解」

がなり立ち、さらに

がなり立つように L を決め、それを積分して での解が求 まるということになる。

と書くと非常に大変そうだが、係数をあらかじめ計算しておけば上の手続きが 陰的ルンゲクッタ公式として表現できる。これで、例えば 4 段で 8 次の公式 が実現出来るわけである。

1.5 代数方程式の扱い

さて、上のようなわけで陰的公式というものが構成できるということはわかっ たが、これを実際に使うには、 に対する代数方程式をとかないといけ ない。この解き方について簡単に説明しておこう。

代数方程式といっても、ここで扱うものはほとんど必ず

 

の形、すなわち、ある関数の不動点を求めるという形になっている。この時は、 以下のような反復を繰り返すことによって解くという方法が考えられる。

これを直接代入法という。

これは容易にプログラム出来て、実際に広く使われる方法である。が、理論的 にはいろいろ怪しい。つまり、

というような観点からみて使いものになるのかどうかというのは明らかではな い。

しかし、陰的ルンゲクッタから出てくる方程式の場合、これで結構うまくいく ことが多い。それがなぜかということを示す前に、一般に直接代入の収束の速 さを調べておこう。

方程式(27)の解が x であるとき、 直接代入法での古い近似 解と新しい近似解それぞれの真の解との差は

となる。下の方を変形し、fを展開すれば

となる。従って、

ということがわかる。

話を陰的ルンゲクッタに戻そう。

また、もっとも簡単な中点公式(陰的中点公式)で考えてみる。 解くべき式は、 として

 

と書ける。 例によってリプシッツ連続条件

がなり立つとしてよいので、中点公式32の右辺の x による微分が で押えられるということがわかる。つまり、刻み 幅を に比べて小さくとればいいということになる。

なお、これでは済まないような場合もある。これについては後の「硬い」方程 式のところで扱うことにする。

1.6 ガウス型以外の陰的公式

陰的ガウス公式は、段数に対して理論上可能な最大次数を達成する、ある意味 で「最良」の公式である。が、実用上はいろいろ不便な点もある。例えば、

などの問題がある。それぞれの点について改良された沢山の方法が知られてい る。半陰的法については別に述べる。



Jun Makino
Sun May 30 16:57:11 JST 1999