11.2 岡田としもり俊宣 「円盤円柱領域における熱方程式に対する差分法」

http://nalab.mind.meiji.ac.jp/~mk/labo/report/open/2004-okada.pdf で公開中 である。

1998年度卒研の遠藤洋一・高木章裕・内藤達也「円盤領域における熱方程式の研究」 の後を継いだものである。 そこでは空間変数を極座標変換してから差分法を用いて、 陽解法と $\theta$ 方向のみ陰的に扱った「半陰解法」を扱い、 安定性について実験と「推測 (山勘とも言う)」で分析している。

この問題について、 筆者 (桂田) はつい最近まで、

半陰解法で安定性の条件がかなり緩和されたが、
ADI 法にしたら (長方形領域と同様に) 無条件安定になるのでは？

と能天気に予想していた。 ところが実際にプログラムを作って実験してみたら、 ADI 法にしても安定性条件は緩和されない、という予想外の事態に直面した (かなり驚いたが、考えてみれば根拠に乏しい予想であった)。 ここで少し冷静になって、差分法の1ステップを

$$\Vector{U}^{n+1}=A\Vector{U}^{n}$$

と書いたときの行列 $A$ のスペクトル半径 $r(A)$ を調べることを思いついた。 MATLAB を使えば比較的簡単に計算できるだろう…

それを誰かやりませんか、 と挑戦者を募ったところ名乗りをあげてくれたのが岡田君だった、 ということである。その結果は次のようなものとなった。

(i)

陽解法、半陰解法、ADI 法について、これまでの実験結果を確認した。

(ii)

陽解法、半陰解法については、 時間刻み幅 $\tau$ が大きいと $r(A)>1$ であるが、 $\tau$ を小さくしていくと $r(A)\kinji 1$ (format long での表示は $1$ ) になる。 ところが、ADI 法では、試してみたいかなる $\tau$ に対しても $r(A)\kinji 1$ (format long での表示は $1$ ) になった。

…またまた予想外の結果である。 時間刻み幅 $\tau$ が大きい場合、 (i) からは不安定と思われるのに、 (ii) からは $r(A)=1$ であるように思える。 とっさに、 $A$ の Jordan 細胞のサイズが $2$ 以上になるということか？ (長方形領域の場合とは異なり、行列 $A$ は対称でないことに注意) と考えたのだが、固有値は重根ではなさそうなので¹⁹、それはない。 (i), (ii) どちらかが間違っている、ということだろうか。 謎が深まる。

というわけで、 この課題はまた誰かの再挑戦を待つ、ということになる。 次回は長方形領域における ADI 法の安定性の証明 (結果は「常識」だが、 筆者は証明を見たことがない) も込めてやってみたい。

岡田君は、C のみならず、 桂田の 30 分レクチャーで MATLAB を使わされての奮闘。 (少々無責任な言い方になるが) よく出来たものだと思う。 彼の今回の結果を詳しくチェックすることは非常に重要である。

(念のため: 結果がほぼ確実に矛盾しているわけで、 どこかで間違えているのだろうけれど、 道の途中ということで、 責められるべきことではない。)