☆☆☆数学の基礎概念☆☆☆
微分方程式
－微分方程式－
独立変数ｘと未知関数ｙおよびその導関数ｙ'，ｙ"，・・・を含む方程式を微分方程式という。
ｙが１変数の関数つまり独立変数が１つならば常微分方程式、多変数の関数とその偏導関数を含み独立変数が２つ以上の方程式ならば偏微分方程式という。
１つの微分方程式の中に含まれる微係数の階数の最大値は、その微分方程式の階数と呼ばれる。階数ｎの微分方程式はｎ階微分方程式となる。
微分方程式が未知関数およびその導関数の多項式である時、最高階の微係数の次数を、その微分方程式の次数という。
また、未知関数およびそのすべての導関数について、１次または０次の微分方程式は線形微分方程式と呼ばれる。
一般に、ｎ階微分方程式は、ｎ個の任意の定数を含む解を持つ。そのような解は一般解、そのｎ個の定数は任意定数という。
ｎ個の任意定数にある特定の値を代入して得られる解は特殊解と呼ばれる。
一般解の中の任意定数にいかなる値を代入しても得られない解があれば、それは特異解となる。微分方程式の一般解と特異解を求めることを、微分方程式を解くという。
１．常微分方程式
変数分離形
次の形の微分方程式は変数分離形と呼ばれている。
| dｙ |
 |
| dｘ |
|
＝ |
ｆ(ｘ)ｇ(ｙ) |
この一般解は，両辺をｇ(ｙ)[≠０と仮定]で割って、両辺をｘで積分する。
 |
| dｙ |
|
| ｇ(ｙ) |
|
|
＝ |
|
ｆ(ｘ)dｘ＋ｃ |
ここで、ｃは定数である。もしも、ｇ(ｙ0)＝０となるようなｙ0が存在するならば、ｙ＝ｙ0も解となる。
(変数分離形の事例)
| dｙ |
|
| dｘ |
|
＝ |
| ｘ2(ｙ＋１) |
|
| ｙ2(ｘ－１) |
|
この式の解を求めてみる。これを次のように変形して、変数が分離した形にする。
| dｙ |
|
| dｘ |
|
＝ |
| ｘ2 |
|
| (ｘ－１) |
|
・ |
| (ｙ＋１) |
|
| ｙ2 |
|
| ｙ2 |
|
| (ｙ＋１) |
|
dｙ＝ |
| ｘ2 |
|
| (ｘ－１) |
|
dｘ |
次に、この左辺と右辺をそれぞれ積分すると、次式を得る。
|
| ｙ2 |
|
| (ｙ＋１) |
|
|
dｙ |
＝ |
|
|
dｙ＝(１／２)ｙ2－ｙ＋log(ｙ＋１) |
|
| ｘ2 |
|
| (ｘ－１) |
|
|
dｘ |
＝ |
|
|
dｘ＝(１／２)ｘ2＋ｘ＋log(ｘ－１) |
したがって、一般解は、
(１／２)ｙ2－ｙ＋log(ｙ＋１)＝(１／２)ｘ2＋ｘ＋log(ｘ－１)＋ｃ
となる。これを両辺に２を掛けて、右辺を左辺に移項し、変形すると、
|
(ｙ＋ｘ)(ｙ－ｘ－２)＋
|
log |
| (ｙ＋１)2 |
|
| (ｘ－１)2 |
|
－２ｃ＝０ |
となる。ここで、定数項－２ｃを改めて、Ｃ＝－２ｃと置くと、
|
(ｙ＋ｘ)(ｙ－ｘ－２)＋
|
log |
| (ｙ＋１)2 |
|
| (ｘ－１)2 |
|
＋Ｃ＝０ |
となる。なお、ｙ＝－１も解となる。
同次形
次の形の微分方程式は同次形と呼ばれている。
これは、ｙ／ｘ＝ｖ，すなわち，ｙ＝ｖｘ とおいて、これを微分すれば、次式を得る。
| dｙ |
|
| dｘ |
|
＝ｖ＋ｘ |
| dｖ |
|
| dｘ |
|
これらを代入すれば、変数分離形に変形される。
| ｖ＋ｘ |
| dｖ |
|
| dｘ |
|
＝ｆ(ｖ) |
| dｘ |
|
| ｘ |
|
＝ |
| dｖ |
|
| ｆ(ｖ)－ｖ |
|
したがって、これを積分すれば、一般解が得られる。
| logｘ＝ |
|
| dｖ |
|
| ｆ(ｖ)－ｖ |
|
|
＋ｃ |
ここで、ｖ＝ｙ／ｘとすれば、一般解になる。
(同次形の解法事例)
(ｘ3＋ｙ3)dｘ－３ｘｙ2dｙ＝０
これは、両辺をｘ3で割ってみると、同次形の微分方程式であることがわかる。
ここで、ｙ＝ｘｖと置くと、dｙ＝ｖdｘ＋ｘdｖとなるので、これを代入すれば、次のように展開できる。
(１＋ｖ3)dｘ－３ｖ2(ｖdｘ＋ｘdｖ)＝０
(１－２ｖ3)dｘ－３ｘｖ2dｖ＝０
| dｘ |
|
| ｘ |
|
－ |
| ３ｖ2 |
|
| １－２ｖ3 |
|
dｖ＝０ |
logｘ＋(１／２)log(１－２ｖ3)＝ｃ'
logｘ2(１－２ｖ3)＝２ｃ'
ｘ2(１－２ｖ3)＝ｃ
但し、ｖ＝ｙ／ｘであるから、
ｘ3－２ｙ3＝ｃｘ
を得る。
－数値計算による常微分方程式の解法－
最も簡単な１階の常微分方程式は、
| dｘ |
|
| dｔ |
|
＝ |
ｆ(ｔ．ｘ) |
で与えられる。
数値計算・オイラー法
１階の常微分方程式は、近似的に、次式で与えられる。
| dｘ |
|
| dｔ |
|
≒ |
| ｘ(t+h)－ｘ(t) |
|
| ｈ |
|
＝ |
ｆ(ｔ,ｘ) |
ｘ(t+h)＝ｘ(t)＋ｈｆ(ｔ,ｘ)
オイラー法では、上式の微分方程式を数値家計算で近似的に解く。
最初、初期値ｔ0の時、ｘ(ｔ0)(＝ｘ0)を与える。
これを与えられた微分方程式に代入すれば、ｔ0における傾きとｔ0＋ｈでのｘ(t+h)の推定値が求められる。
次に、ｔ0＋ｈを新たなｔ0と考えて、同じような手続きを繰り返す。
この場合、キザミ幅ｈが大きいと、真のｘ(t0+h)とその推定値の間に大きな差が生じる。
この差はキザミ幅ｈを十分に小さくすれば、この差を小さくすることができる。
しかし、この計算を繰り返すことで、その累積誤差は次第に大きくなる。
数値計算・ルンゲ-クッタ法
オイラー法では、ｔ0での傾きを利用して、ｘ(ｔ0＋ｈ)を推定した。
ルンゲ-クッタ法では、複数の点での傾きを利用する。
ｔの値をｈ刻みに取り、
ｔn＝ｔ0＋ｎｈ（ｎ＝１，２，３，・・・）
として、ｘの値ｘnを次式で順々に求める。
ｘn+1＝ｘn＋（ｋ1＋２ｋ2＋２ｋ3＋ｋ4）／６
ここで、
| ｋ1＝ｈｆ(ｔn，ｘn) |
| ｋ2＝ｈｆ(ｔn＋ｈ／２，ｘn＋ｋ1／２) |
| ｋ3＝ｈｆ(ｔn＋ｈ／２，ｘn＋ｋ2／２) |
| ｋ4＝ｈｆ(ｔn＋ｈ，ｘn＋ｋ3) |
数値計算・ミルン法(予測子修正子法)
さらに常微分方程式の数値解の精度を上げる方法として、ミルン法(予測子修正子法)がある。
この方法は、やはり独立変数ｔをｈ刻みに分割し、
ｔn＝ｔ0＋ｎｈ（ｎ＝１，２，３，・・・）
として、ｘの値ｘnを順々に求めるが、ｘn-1まで求めた時に、ｘnの値を予測し、
ｆn＝ｆ(ｔn，ｘn)
を使って修正する。ミルン法では、ｘnの予測値ｘnPは、
ｘnP＝ｘn-4＋(４ｈ／３)(２ｆn-3－ｆn-2＋２ｆn-1)
で与える。修正されたｘnCは、
ｘnC＝ｘn-2＋(ｈ／３)(ｆn-2＋４ｆn-1＋ｆn)
で与えられる。ただし、右辺のｆnは、最初にｆ(ｔn，ｘnP)を使用するが、
その次からは修正値ｘnCを用いてｆ(ｔn，ｘnC)を求め、修正値ｘnCを反復計算する。
そして、修正値の変化が予め与えて置いた微少量εより小さくなれば反復を止める。一般的には反復回数が数回で修正値が最終値ｘnCになるようにする。
２．偏微分方程式
偏導関数を含む方程式を偏微分方程式という。例えば、独立変数がｘとｙの関数ｚに関する偏導関数を含む方程式が次のように与えられているとする。
| Ｆ |
 |
ｘ，ｙ，ｚ， |
| ∂ｚ |
|
| ∂ｘ |
|
， |
| ∂ｚ |
|
| ∂ｙ |
|
， |
| ∂2ｚ |
|
| ∂ｘ2 |
|
， |
| ∂2ｚ |
|
| ∂ｘ∂ｙ |
|
， |
| ∂2ｚ |
|
| ∂ｙ2 |
|
 |
＝０ |
これは、２階の偏導関数を含む方程式であり、２階偏微分方程式と呼ばれる。もしも、１階の偏導関数しか含まれなければ、１階偏微分方程式と呼ばれる。いま、偏導関数を次のように表示する。
| ｐ＝ |
| ∂ｚ |
|
| ∂ｘ |
|
，ｑ＝ |
| ∂ｚ |
|
| ∂ｙ |
|
，ｒ＝ |
| ∂2ｚ |
|
| ∂ｘ2 |
|
，ｓ＝ |
| ∂2ｚ |
|
| ∂ｘ∂ｙ |
|
，ｔ＝ |
| ∂2ｚ |
|
| ∂ｙ2 |
|
すると、２階偏微分方程式は次のように表示することができる。
Ｆ（ｘ，ｙ，ｚ，ｐ，ｑ，ｒ，ｓ，ｔ）＝０
また、１階偏微分方程式は次のように表示することができる。
Ｆ（ｘ，ｙ，ｚ，ｐ，ｑ）＝０
なお、偏導関数は次のような記号を用いることもある。
| ｆz＝ |
| ∂ｆ |
|
| ∂ｚ |
|
，Ｆp＝ |
| ∂Ｆ |
|
| ∂ｐ |
|
，Ｇyq＝ |
| ∂2Ｇ |
|
| ∂ｙ∂ｑ |
|
２つの独立変数ｘ，ｙを持つ偏微分方程式は、恒等的に満足する関数ｚ(ｘ，ｙ)が存在すれば、その関数が方程式の解となる。そして、解を求めることを、その方程式を解く、または積分するという。
この時、偏微分方程式には、次の３種類の解が存在する。
- 完全解 ２つの任意定数を含む解
- 一般解 １つの任意定数を含む解
- 特異解 完全解にも一般解にも含まれない解(必ずしも存在しない)
偏微分方程式には、標準形と呼ばれる次に示す形のものがある。そして、これらの形のものは、容易に完全解を求めることができる。
- Ｆ(ｐ，ｑ)＝０
- Ｆ(ｘ，ｐ，ｑ)＝０
- Ｆ(ｙ，ｐ，ｑ)＝０
- Ｆ(ｚ，ｐ，ｑ)＝０
- Ｆ(ｘ，ｐ)＝Ｇ(ｙ，ｑ)
- ｚ＝ｐｘ＋ｑｙ＋Ｆ(ｐ，ｑ)
(１階偏微分方程式の完全解の事例)
いま、ｘ，ｙ，ｚの間に次の関係があったとする。
(ｘ－ａ)2＋(ｙ－ｂ)2＋ｚ2＝１ 但し、ａ，ｂは任意定数
これをｘとｙについて、偏微分すると、次の式を得る。
２(ｘ－ａ)＋２ｚｐ＝０
２(ｙ－ｂ)＋２ｚｑ＝０
これらの３つの式から、任意定数ａ，ｂを消去する。
ｚ2(ｐ2＋ｑ2＋１)＝１
つまり、この１階偏微分方程式は、次の完全解を持つことを意味する。
(１階偏微分方程式の一般解の解法事例)
いま、次の１階偏微分方程式が与えられたとする。
(ｙ－ｚ)ｐ＋(ｚ－ｘ)ｑ＝ｚ－ｙ
この方程式は、ラグランジェの微分方程式と呼ばれ、次の形をしている。
Ｐ(ｘ，ｙ，ｚ)ｐ＋Ｑ(ｘ，ｙ，ｚ)ｑ＝Ｒ(ｘ，ｙ，ｚ)
この方程式を解くには、最初に、次の補助方程式と呼ばれる常微分方程式を解く。
この場合、その解として、次式が得られる。
ｕ(ｘ，ｙ，ｚ)＝ａ
ｖ(ｘ，ｙ，ｚ)＝ｂ
また、任意の２変数ｕとｖの関数ｆ(ｕ，ｖ)から作った次式も解となる。
ｆ(ｕ，ｖ)＝０
すなわち、これらが求める一般解となる。したがって、与えられた１階偏微分方程式の補助方程式は、次のようになる。
これから次式を得る。
dｘ＋dｙ＋dｚ＝０
ｘdｘ＋ｙdｙ＋ｚdｚ＝０
これらを積分すると、次式になる。
ｘ＋ｙ＋ｚ＝ａ
ｘ2＋ｙ2＋ｚ2＝ｂ
したがって、求める一般解は、次のようになる。
ｆ(ｘ＋ｙ＋ｚ，ｘ2＋ｙ2＋ｚ2)＝０
ここで、ｆは任意の関数である。
(２階線形偏微分方程式)
| ｐ＝ |
| ∂ｚ |
|
| ∂ｘ |
|
，ｑ＝ |
| ∂ｚ |
|
| ∂ｙ |
|
，ｒ＝ |
| ∂2ｚ |
|
| ∂ｘ2 |
|
，ｓ＝ |
| ∂2ｚ |
|
| ∂ｘ∂ｙ |
|
，ｔ＝ |
| ∂2ｚ |
|
| ∂ｙ2 |
|
Ｒｒ＋Ｓｓ＋Ｔｔ＋Ｐｐ＋Ｑｑ＋Ｚｚ＝Ｆ
(係数 Ｒ，Ｓ，・・・，Ｆはｘ，ｙの関数)
これは次の４つの形に分類される。
- ｒ＝Ｆ， ｓ＝Ｆ， ｔ＝Ｆ
- ｒ＋Ｐｐ＝Ｆ， ｓ＋Ｐｐ＝Ｆ， ｓ＋Ｑｑ＝Ｆ， ｔ＋Ｐｑ＝Ｆ
- Ｒｒ＋Ｓｓ＋Ｐｐ＝Ｆ， Ｓｓ＋Ｔｔ＋Ｑｑ＝Ｆ
- Ｒｒ＋Ｐｐ＋Ｚｚ＝Ｆ， Ｔｔ＋Ｑｑ＋Ｚｚ＝Ｆ
３．差分方程式(階差方程式)
代表的な偏微分方程式として、放物型方程式と楕円型方程式がある。放物型方程式は、熱伝導方程式や拡散方程式が代表される。
| ∂ｕ |
|
| ∂ｔ |
|
＝ｋ2 |
| ∂2ｕ |
|
| ∂ｘ2 |
|
一方、楕円型方程式は、ラプラス方程式やポアソン方程式が代表される。
| ∂2ｕ |
|
| ∂ｘ2 |
|
＋ |
| ∂2ｕ |
|
| ∂ｙ2 |
|
＝０ |
| ∂2ｕ |
|
| ∂ｘ2 |
|
＋ |
| ∂2ｕ |
|
| ∂ｙ2 |
|
＝ρ |
(放物型方程式の数値解法)
| ∂ｕ |
|
| ∂ｔ |
|
＝ｃ2 |
| ∂2ｕ |
|
| ∂ｘ2 |
|
いま、これを次の初期条件および境界条件で解くことを考える。
初期条件 ｕ(ｘ，０)＝ｆ(x)
境界条件 ｕ(０，ｔ)＝ｕ(ａ，ｔ)＝０
この時、ｘ軸の領域(0≤x≤a)をＭ等分して、その刻みをｈ＝ａ／Ｍとする。また、ｔ軸の刻みをｋで刻むことにする。そして、一般に、ｕmnは次のように表示する。
ｕmn＝ｕ(mh,nk)
この場合、放物型微分方程式は、近似的に差分式(階差式)で置き換えることができる。
| ｕm,n+1－ｕm,n |
|
| ｋ |
|
＝ｃ2 |
| ｕm-1,n－２ｕm,n＋ｕm+1,n |
|
| ｈ2 |
|
ここで、簡単のために、
ｒ＝ｋｃ2／ｈ2
置くと、
ｕm,n+1＝ｒｕm-1,n＋（１－２ｒ）ｕm,n＋ｒｕm+1,n
が得られる。これを用いれば、ｕm-1,n，ｕm,n，ｕm+1,nから、ｕm,n+1を求めることができる。
なお、ｍ＝１，２，３，・・・，Ｍ－１であり、ｕ0,n+1とｕM,n+1は、境界条件から、常に０である。
この時、ｒを小さくすると、ｋ＝ｒｈ2／ｃ2が小さくなり、時間軸の進みは遅くなる。
しかし、ｒを大きくすると、不安定になる。もしも、ｒ＞１／２にすると必ず不安定になることが証明されている。
(楕円型方程式の数値解法)
| ∂2ｕ |
|
| ∂ｘ2 |
|
＋ |
| ∂2ｕ |
|
| ∂ｙ2 |
|
＝０ |
これはラプラス方程式と呼ばれ、境界Ｄに囲まれた領域Ｃで成立する。これを下記の境界条件の下で解くことを考える。
ｕ(ｘ，ｙ)＝Φ(ｘ，ｙ)
いま、平面に間隔ｈｋ格子点を設定し、その境界はｘ軸とｙ軸の格子点で交わるとする。
領域内の格子点(ｉ，ｊ)について、番号ｉは１からＩＭまで、番号ｊはＪＬからＪＵまで与えられているとする。
格子点(ｉ，ｊ)でｕijが存在するとすれば、ラプラス方程式は、次式で近似できる。
| １ |
|
| ４ｈ2 |
|
(ｕi+1,j－２ｕi,j＋ｕi-1,j)＋ |
| １ |
|
| ４ｈ2 |
|
(ｕi,j+1－２ｕi,j＋ｕi,j-1)＝０ |
これは、整理すると、次式に変形できる。
ｕi,j＝(ｕi+1,j＋ｕi-1,j＋ｕi,j+1＋ｕi,j-1)／４
これを第ｋ近似ｕi,jkから第ｋ＋１近似ｕi,jk+1へ、順次求めて行く方法は、緩和法と呼ばれる。
ｕi,jk+1＝(ｕi+1,jk＋ｕi-1,jk＋ｕi,j+1k＋ｕi,j-1k)／４
しかし、格子点数が多くなると、収束が遅くなる。そこで、加速係数ωを用いて、次のように修正する。この方法は加速緩和法と呼ばれる。
ｕi,jk+1＝ｕi,jk＋ω｛(ｕi+1,jk＋ｕi-1,jk＋ｕi,j+1k＋ｕi,j-1k)／４－ｕi,jk｝
さらに、既に第ｋ＋１近似ｕi,jk+1が求められていれば、それを採用する方法がある。これは逐次式加速緩和法(ＳＯＲ：Successive Over-Relaxation)と呼ばれる。
ｕi,jk+1＝ｕi,jk＋ω｛(ｕi+1,jk＋ｕi-1,jk+1＋ｕi,j+1k＋ｕi,j-1k+1)／４－ｕi,jk｝
以上
（２０１３年３月２４日）
戻る
|
