偏微分方程式にトライ（波動方程式1D） - 理系ツールのクールな使いこなし

独立変数が複数になるので、ODEに比べ複雑さは桁違いだ。それでも適用アプリケーションが存在するので、できるだけ紹介していきたい。

方程式には1階、2階という階数があるが、自然界の法則を表現するPDEには2階が非常に多い。何故多いのかは上手く説明できないが、私的には2次関数の考え方をアナロジー的に利用しているからと思っている。これまで3次関数の延長線での応用など聞いた事が無いし、線形代数の公式も3次になると途端に面倒になるし、2から3へのステップアップは「数学の専門分野」の話になるのかもしれない。まずは整理されている2階線形からスタート。

2階線形PDEの分類

現在、様々な参考書では次のように扱われている。

2階線形PDEは、独立変数が２つの場合は、

式（1）で表せる。ここでuはxとyの関数で、A～Gはxとyの関数または定数である。最終的には表1のように大きく3種類に分類される。

判別式	根	型	代表PDE
B²-4AC > 0	2実根	双曲型	波動方程式
B²-4AC = 0	重根（実根）	放物型	熱伝導方程式
B²-4AC < 0	共役な2複素根	楕円型	ラプラス方程式

表1

これ以降、3種類のPDEについて、例を挙げていく。ODEの説明では解析解と非解析解に分けたが、PDEは具体的に初期条件などを与えて全て非解析解（級数解、数値計算解）とし、出来るだけ視覚的にも工夫を凝らして説明を進めていく。

A1 　波動方程式（一次元）パート1

一次元波動方程式を調べていくと、無限に長い弦の振動から、２つの任意関数を含む一般解「ダランベールの解」にたどり着く（式（2））。

更に工学的利用を考えて初期条件を考慮すると、以下の「ストークスの波動公式」が得られる。

初期条件は、

ここで、αとβはxを独立変数とする任意関数である。数学的な振る舞いを調べるにはこれらの形が都合良いらしい。

工学的には初期条件だけでなく境界条件も必要な場合が多い、これにアプローチするにはフーリエの方法と言われ、変数分離やフーリエ級数などを用いて進められる。ここまでは解析解と言って良い。

計算例

（1）数式処理ツールMapleの計算結果を以下に示す。PDEなのでpdsolveコマンド類を使っている。初期条件はu(x,0)=?、1回微分はD[2](u)(x,0)=?の様に入力する決まり（1番目の引数で微分する場合はD[1]、2番目の引数で微分する場合はD[2]）になっている。#を先頭に付けると、それ以降は計算実行に影響を与えない説明文になる。unionは初期条件と境界条件を連携させるために入れてある。

解plotを見ると初期条件は弦の真ん中を上に引張ったイメージである事が解るだろう。MapleにはPDE数値計算手法が色々と用意されているが、時間tの2階微分を用いる波動方程式の解法は残念ながら有限差分法のみである（それでpdsolve(,,,numeric)を使った）。それでもプログラミングを利用するまでもなく解を可視化できる所が凄い。

------------------------------------------------------------------------------------

-------------------------------------------------------------------------------------

（2）次はWolfram CloudのMathematicaだ。

今回は計算に必要なメモリーが多く、無料でクラウド計算できる範囲を超えてしまい結果が得られなかった。取り敢えず、果敢に取り組んだ事だけ見てもらおう。NDコマンドを利用した。Mathematicaを使用できる環境にいる方は是非試してほしい。

-------------------------------------------------------------------------------------

（3）次はFlexPDEという名前のツールで、PDE Solutions Inc.のHPからダウンロードできる。その名の通りPDE数値解析専門ツールだ。Liteバージョンなら、制限はあるものの試用可能だ。計算原理は有限要素法（以下FEM）なので、機械系の解析がメリット多い感じだ。FEMに乗せられれば、構成方程式（数学記号に似ている）を選ばない。その構成方程式と共に初期条件や境界条件を入力する。様々なPlotコマンドがあるので、最後に試行錯誤しながら解をPlotしていく。数値解析ツールなので、得られるのは数値データだけで数式は得られない。

初めての紹介なので、プログラムコードの補足を加える。赤い字は約束された記号である。TITLEは解く問題の名前。SELECTでは計算誤差を1/10000以下にするように指示できる。VARIABLEではその名の通り変数を指定する。FlexPDEは時間の2階は積分できないので、本問題では二つのPDEに分けて連立して解くしかない。DEFINITIONSでは構成式中の定数などを指定する。本問題ではcとLそして計算スタート時間t0だ。INITIAL VALUESでは初期値を指定する。本問題ではvはuを時間tで微分したものである。EQUATIONSでは、その名の通り構成方程式を指定する。本問題では先に説明したように時間の2階微分を1階に落とすために連立式にした。BOUDARIESは境界条件を指定する。本問題では大きさはy方向、時間軸はx方向にして二次元で表現する。FEMはこの二次元平面を要素分割して計算を進めていく。start(0,L)は計算スタート地点の座標で、このスタート地点から二次元平面を反時計回りに計算していく習わしである。最後のLine to closeで二次元平面を閉じている。座標間に初期条件や境界条件を指示する。valueはディリクレ形境界条件で、naturalはノイマン型境界条件だ。TIMEは計算する時間範囲を指定している。PLOTSは様々なplotコマンドを入力するエリアだ。最後はENDでプログラムを閉じる。