拡散方程式の陰解法

お世話になります。長文で失礼します。

最近になって偏微分方程式の数値解析の勉強を始めました。

１）早速質問したいですが、二次元での拡散方程式を、
　陰解法で数値的に解く場合、途中で出てくる連立方程式の係数行列は、
　どのような構造になるのでしょうか。

　手持ちの教科書のすべてが一次元の場合を例示しており、
　係数行列は三重対角行列になる、と説明しているのですが、
　二次元になると係数行列はどうなるのでしょうか。

　試しに「陰解法　対角行列」とググってみると、
　「5重対角行列」と記載のあるリンクが出てくるのですが、
　http://www-mmc.es.hokudai.ac.jp/else/cdrom/main_ …
　実際に読んでみても、「疎な5重対角行列」と言及されるだけで、
　具体例は示されていませんでした。

　結局、2次元問題の陰解法の実装にあたって
　係数行列をどのように定義すればよいのか分からず、
　添付の画像の式３ようにインデックスを振り直し、
　計算領域の行列uをベクトル化（式４）して、
　疑似的な五重対角行列Aを無理矢理定義しました（式６、周期境界条件）。

　果たして合っているのでしょうか。
　一応、空間刻みと時間刻みを揃えて陽解法と計算結果を比べましたが、
　同じような挙動になりました。
　合っていない場合、係数行列はどのような構造になるのでしょうか。

２）上記の方法が合っているのだとすると、
　2次元問題の陰解法で用いる「5重対角行列」というのは、
　式６の係数行列のように対角成分が0の隙間が出来てしまい、
　厳密な「多重対角行列」の定義から外れると思うのですが、
　この理解は正しいでしょうか。

　むしろ逆に、対角成分が0の隙間が出来てしまっていても、
　「～重対角行列」と呼んで良いものなのでしょうか。

３）同様に三次元に拡張した陰解法スキームを考えると、
　係数行列は二次元のときよりも広い0の隙間を持った、
　七重対角行列になると思うのですが、この理解も合っているでしょうか。


どなたか詳しい方がいらっしゃいましたら、
ご回答よろしくお願い致します。

質問者からの補足コメント

• 「帯行列」の存在を教えてくださったANo.2の回答もありがたかったのですが、
  今回は情報量が多かったANo.4をベストアンサーに選ばせて頂きます。
  
  ありがとうございました。。。

    補足日時：2018/03/30 14:34

No.4 ベストアンサー 回答者： stomachman 回答日時： 2018/03/23 15:44

No.3へのコメントについてです。

どうもまともに読んでいただけていないようなんで、ま、どうでもいいんですが。

> 私が例示した行列uは4x4行列で、係数行列Aは16x16行列だったのですが、
> これを「2n+1」の係数行列に直した方が良い

　係数行列のサイズは、当然ながら、空間をいくつに分割したかによって決まります。たとえばn=2において空間を100×100に分けたとすると10000個の成分がある。そのうちの「合計2n+1個の成分に0でない係数がついた一次方程式を連立したもの」、すなわち10000×10000の行列の各行に5個だけ0でない係数が入っている。この行列の10000×5個の成分だけが0でないということです。

> なぜ共役勾配法だけ

　「行列はなるべく相手にしたくないという場合」、これなら行列を直にいじくらなくても良いからです。ことにご質問の例のように拡散係数が位置依存性を持たない（定数である）場合には、10000×(2n+1)個の成分を保持することも必要ない。

この回答へのお礼

辛抱強く何度も対応して頂き、ありがとうございます。
言い換えて頂いたおかげで、何とか理解できたと思います。

ANo.1にて、

＞空間を四角く切ることにすると、一つのセルに隣接するセルがn次元の場合2n個あるんで

という表現を、「一つのセルに隣接するセルがn次元である場合」と読んでしまい、
数値1個が入る各セルが何故か「n次元」？と意味不明なまま読み進めてしまいました。
直後の「2n個あるんで」も、何が2n個あるのか分からないまま読みました。

「n次元の空間を四角く切ることにすると」という意味だったのですね泣

はじめに頂いた文面を初見で理解できない程度に、私の到達度の初心者ぶりを自覚できました。

また、共役勾配法のメリットについても追記頂きまして、ありがとうございました。
自分で実装してみて、メモリ節約のありがたみを体感してみようと思います。


この度は、ありがとうございました。

お礼日時：2018/03/23 17:08

No.3 回答者： stomachman 回答日時： 2018/03/23 12:38

ANo.1へのコメントについてです。

> 少なくとも質問に対する回答が優先的に欲しい

3つに分けるまでもなく、まとめて2n+1だと申し上げたんですがね。

この回答へのお礼

追記ありがとうございます。

ご記載の内容はシンプルで理解できましたが、私が初学者のため、
ごまかさずに問題を細分化して個別に確認したかった、という次第です。

不満点があるとすれば、
私が例示した行列uは4x4行列で、係数行列Aは16x16行列だったのですが、
これを「2n+1」の係数行列に直した方が良い理由は、
先のご説明だけでは分からなかったので、回答意図の曖昧さを感じました。

また、(半)陰解法の中でもクランク・ニコルソン法やADI法、畳み込み法など、
いくつか戦法があると思うのですが、なぜ共役勾配法だけを例示して下さったのか、
「自分が好みでよく使うんですが」程度の理由でも構わないので、
オススメの根拠が欲しかったです。この点も回答意図の曖昧さを感じました。

入門者の学習補助というより中級者の向学心を誘うようなスタンスを感じ、
私の需要とは今回は合致しませんでした。
（回答する側にとっては判断が難しいでしょうが笑


ありがとうございました。

お礼日時：2018/03/23 13:11

No.2 回答者： nishi16 回答日時： 2018/03/23 10:54

全く回答にはなっていないのですが，知っている範囲でお答えしますと，

1 ) 解き方は正しいと思います。
2 ) 示された行列は対角行列ではなく，帯行列です。
対角行列の定義は質問者様の認識で正しいです。
ただ，自分は帯行列になると教わったんですよね。
3 ) おっしゃる通りです。

残念ながら数学が専門ではないので，
そのような方がいらっしゃれば，そちらの意見の方が正しいかと思われますり

この回答へのお礼

早速ご回答ありがとうございます。

１）合っているようで安堵しました。
２）「帯行列」という用語は耳にしたことがあったのですが、
　卒業した学部カリキュラムに線形代数の講義が無かったため、
　定義は知りませんでした。勉強になりました。

クリティカルな回答を頂きまして、大変感謝申し上げます。
ご提案通りに、もう少し他の回答も待ってみたいと思います。

お礼日時：2018/03/23 11:41

No.1 回答者： stomachman 回答日時： 2018/03/23 10:46

単なる拡散方程式（熱方程式）なら、それこそ陽解法が適しているんじゃないでしょうかね。

最初は空間を粗く切ってざくざく計算し、収束に近づいたらメッシュを小さくする、という風にやれば速くて、行列には出番なし。
　というのはさておき、ともあれ陰解法でやってみたいと。
　空間を四角く切ることにすると、一つのセルに隣接するセルがn次元の場合2n個あるんで、合計2n+1個の成分に0でない係数がついた一次方程式を連立したものを解く必要がある。ご質問にあるような、端っこに「折り返し」の成分がついた「隙間のある2n+1重対角行列」っぽいもの、ってことです。まともに逆行列を計算してたら大変。
　行列はなるべく相手にしたくないという場合、答のほうから探っていく、すなわち近似解を改良していく方法があります。近似解の出発値は陽解法の１ステップで計算すればいいでしょう。
　　Ax = y
という方程式の近似解 x*の持つ誤差εを
　　Ax* = y + ε
としますと、x*に微小変動Δxを加えて
　　E = Σ|ε|^2
が小さくなるように改良する、という、非線形最小二乗法でおなじみの手続きを繰り返す。Δxを決める手法は多々ありますが、行列の計算をしないで済ませるには最急降下法を使うことになる。Aはヤコビアン
　　A[i,j]=∂ε[j] /∂x*[i]
になっているわけですから理屈は簡単ですし、適当な加速のテクニックを組み合わせることができます。ただし、どうせこれ、Δt刻みの１ステップにすぎないのですから、xにあんまり精度を求めてもしょうがない。このへんの匙加減がコツってことでしょう。

この回答へのお礼

早速ご回答ありがとうございます。

ご説明頂いた方法は、「共役勾配法」のことでしょうか。
まだ実装したことはないですが、この方法の有効性も確認したいと思います。

本文には3個ほど質問があったのですが、
長文だとやはり読んで頂けないものなんですかね汗

建設的で発展的な情報は勿論歓迎するところなのですが、
少なくとも質問に対する回答が優先的に欲しいところでした笑

お礼日時：2018/03/23 11:37