- ベストアンサー
不完全ＬＵ分解前処理つき双共役勾配法についておしえてください。
連立方程式を解くために不完全ＬＵ分解前処理つき双共役勾配法 について勉強しています。 前処理の際に、行列Ａを不完全ＬＵ分解しその逆行列(LU)^(-1)というのを使用します。LU分解まではできたのですが、この逆行列は普通にＬＵ分解＋直接法という形でもとめるのでしょうか。だとしたら、直接法をつかっていてあまり高速化が期待できない様な気がしました。 不完全コレスキー分解つき共役勾配法(ICCG)のときは、不完全コレスキー分解後、間接的にAの逆行列をもとめて使用する方法がありましたのでなにかいい方法があるのかと思い質問しました。 はじめてのプログラミングで見当違いなことをいっているかもしれませんがよろしくおねがいします。
- nebuta123
- お礼率40% (6/15)
- 数学・算数
- 回答数5
- ありがとう数2
- みんなの回答 （5）
- 専門家の回答
質問者が選んだベストアンサー
あと、p.115の 6.2代入計算 の方が参考になるかもしれません。
その他の回答 (4)
- DreaMMaster
- ベストアンサー率38% (63/164)
小国力先生の『行列計算ソフトウェア―WS、スーパーコン、並列計算機』が手元にあれば、p.88-89の 4.3.3不完全これスキー分解と不完全LU分解 を見てください。
- DreaMMaster
- ベストアンサー率38% (63/164)
#2です。 (L U)^(-1)を左からかける、と言うことは、L^(-1)を左からかけ、更にU^(-1)を左からかける。ということはお分かりでしょうか？ Lは下三角行列ですから、逆行列をかけるためにいちいち逆行列の成分を求めなくても、良いことはおわかりでしょうか？（本質的にガウスの消去法の後半でもあり、LU分解の本質でもあるのですが…。未知ベクトルの第一成分はLの第11成分で定数項の第1成分を割れば良いです。未知ベクトルの第二成分は第一成分が分かれば計算できます。以下同様） Uの逆行列も同様です。 更に言うと、(L U)^(-1)Aをいかに簡略化するかがプログラミング上のポイントの一つですから、丁寧な本であれば、書かれているはずです。 小国力編著、『行列計算ソフトウェア―WS、スーパーコン、並列計算機』、丸善 渡部力,小国力, 名取亮 (監修)、 『Fortran77による数値計算ソフトウェア 』、丸善 などでは、きちんと書かれていて、しかもFortranとはいえソースコードも載っていますから参考になると思います。 古い本しか紹介できず申し訳ありません。最近の動向は追っかけていないので…
補足
ご丁寧な解説ありがとうございます。 まさに、私の質問は、 「更に言うと、(L U)^(-1)Aをいかに簡略化するかがプログラミング上のポイントの一つ」の部分です。 小国力先生の『行列計算ソフトウェア―WS、スーパーコン、並列計算機』を借りてきました。たとえば 10,5,2 ILUBCGでも 10.5.4ILUCGSでも手順１、手順２で 不完全LU分解して(LU)^(-1) をかけているような記述しかみあたりません。 私の持っているテキストと同じ書き方しかないように思えます。 （きっと私が見落としているだけだと思います） すみませんが、どこに記述されているか教えてください。 よろしくお願いします。
- DreaMMaster
- ベストアンサー率38% (63/164)
ICCGと同様にすれば、反復法に持ち込めます。 と言うより、反復法に持ち込めるように、完全に分解しないので不完全分解と言うのですが…。
補足
すみません。私の舌足らずな（間違い）の文のせいでどうやら 質問の真意がお伝えできていないみたいです。 元の文で「LU分解まではできたのですが」というのは「不完全LU 分解はできたと」いう意味で、DreaMMaster様のおっしゃる 「反復法に持ち込めるように、完全に分解しないので不完全分解」 の意味は理解できております。 私がＩＣＣＧと同じようにできないと思っているのは ICCGでは( L D L^(T) )^(-1) r というものを求める必要に対して わざわざ( L D L^(T) )^(-1)を求めないで( L D L^(T) )^(-1) r を一気にもとめているのに対し 不完全ＬＵ分解つきの双共役勾配法や自乗共役勾配法では (L U)^(-1)(b-Ax)や (L U)^(-1)A など (L U)^(-1)を つかう式が違った形で出てきていて、これには(L U)^(-1) 自体を求める必要があるのかという点です。 しかし、(L U)^(-1)をもとめてから(L U)^(-1)(b-Ax)や (L U)^(-1)A を求める演算は結構時間ががかりそうなので疑問に思い質問しました。 よろしくお願いします。
- stomachman
- ベストアンサー率57% (1014/1775)
「不完全LU分解」は知りませんけど、ともかくLU分解をした、つまり下三角行列Lと上三角行列Uとの積の形に分解したのだとします。面倒なので逆行列を*で表すことにすると、 (LU)* = U* L* であり、L*は下三角行列、U*は上三角行列。 LとL*がどっちも下三角行列であることを利用するとL*は簡単に計算できます。というのは、 L L* = E (単位行列) だから、 E[j,j]=L[j,j] L*[j,j] = 1 他の項は全部0なんで、これだけになっちゃう。で、L*[j,j] が決まった。 さて、 E[j,j-1] = L[j,j-1] L*[j-1,j-1] + L[j,j] L*[j,j-1]=0 であるから、L*[j,j-1]が決まった。 という調子です。 U*も同様。
補足
ご丁寧な解説どうもありがとうございます。 逆行列の求め方はよくわかりました。 でも、前処理としてわざわざこれで逆行列を求めた後、 この逆行列を何度もかける演算をして真の逆行列を求めるのでは かえって遅くなってしまうようなきがします。 実際やってみると速いんでしょうか
関連するQ&A
- ラプラス方程式の解法（情報処理として）と共役勾配法
情報処理の数学はカテゴリ違いかもしれませんが、適当なところを見出せませんでしたので。 （以下、長文で申し訳ありません。） ラプラス方程式やポアソン方程式の境界値問題を解く方法にＳＯＲ法が あります。これは連立一次方程式の反復解法の１つだとも言えると思い ます。連立一次方程式の解法は通常AX=BとなってマトリックスAと既知 ベクトルBに対して未知ベクトルＸを解くわけですが、ラプラス方程式 をＳＯＲ法で解く場合、AX=Bという形にはしません（そういう風に表現 はできますが）。例えば2次元を考えると100×100の2次元の領域を解く 場合、未知数は10000個です。これをAX=Bという風に考えたら、10000元 連立一次方程式を解くことになり、A(10000,10000)のようにメモリ確保 も大変ですが、実際にはそのように解きませんね。マトリックスＡとい うメモリの取り方をしないのです。 100×100の格子領域でSOR法で解くとプログラムステップとしては数十 行ぐらいとなり、大したことはありませんね。 しかし、情報処理の連立一次方程式の解法の本を見ると、それはあくま でも連立１次方程式があった場合ということを前提としているのでその ようなメモリの確保から解説がスタートするようです。 そこで質問ですが、ラプラスやポアソン方程式を（前処理付）共役勾配 法などで解く場合、メモリ確保の前提条件としては連立一次方程式の形 ではないはずです。本を見るとマトリックス解法としての共役勾配法が 出ていますが、ラプラス方程式、ポアソン方程式の特化して解説してい る本などないでしょうか。また、簡単に説明できるのであれば教えて頂 きたいのですが。
- ベストアンサー
- 数学・算数
- 共役勾配法の解法について
Ax = b の形の式を共役勾配法で解く方法は，インターネットでいくつか見つけることができました． しかしながら，Ax = b の解が１つでない場合などに， ||Ax-b||^2+α||Cx||^2 （ノルムは２ノルム）を最小にする問題を，共役勾配法によって解く方法がわかりません． 最小自乗法に関する知識が必要なのかと思うのですが，どうにもならないのが現状です． どなたかご存知ありませんでしょうか？
- ベストアンサー
- 数学・算数
- LU分解法（ドゥーリトル法）の問題
7x+7y+6z=13 7x+6y+5z=10 6x+5y+4z=8 この方程式をLU分解法で解き、係数行列をLU分解して得られるL行列とＵ行列の0でない成分及び解ベクトルを求めよ。マトリクスの成分は全て整数か仮分数で表現せよ、また連立方程式の順序を入れ替えずに解くこと。という問題です。 自分で考えてみたものの、講義で習ったのは定義だけで、理解することができなく全く解くことができません。 もう少し自分で考えてみろと思われると思いますが、どなたか解いてくださらないでしょうか・・・
- ベストアンサー
- 数学・算数
- すみません。行列式のLU分解について質問なんですが。
すみません。行列式のLU分解について質問なんですが。 ある行列式AをLU分解するまではいいのですが、その後、どうやって固有値を求めたらいいのでしょうか？ QR分解で求める方法の方が効率が良いのはわかっているのですが、気になるので、教えてもらえたら幸いです。 あと、あまり数学が得意でないので、噛み砕いて説明していただけるとほんと助かります。
- ベストアンサー
- 数学・算数
- 共役勾配法を用いた最適化について
弾性体の変形を見るために，停留ポテンシャル法を用いようとしているのですが，その際に目的関数にひずみエネルギーをとり 独立な24変数の関数の極小値を共役勾配法で求めようとしています． Hesse行列を手計算でだし，最適化を行うプログラムを書いたのですが どうも上手くいきません． 下にプログラムの一部分をのせます． 共役勾配法について勘違いしている可能性もあるので 問題点を教えて頂ければ幸いです． x_vectorは24変数ベクトルで，dは探索方向，HはHesse行列です． for(i=0;i<10;i++) { lambda = (x_vector*d)/(d*H*d); x_save_vector = x_vector; x_vector=x_vector - lambda*H*d;//変数べクトル再設定 if(i==0) { // 初期値 gamma = 1; } else { gamma =x_vector*x_vector/(x_save_vector*x_save_vector); } d = x_vector + gamma*d;//探索方向再設定 }
- ベストアンサー
- 数学・算数
- LU分解法を使った連立方程式の解き方
2x＋３ｙ＋３ｚ＝５ ２ｘ＋２ｙ－ｚ＝－５ ５ｘ＋４y＋２ｚ＝３ この連立方程式をLU分解法を使ってときたいのですがなんどやっても答えが合いません。 大変困っているので、どなたか力を貸していただけませんか？
- 締切済み
- 数学・算数
- 行列のLU分解について教えてください
この行列ＡのＬＵ分解を求めたいのですがたくさん資料などを見ていてもわかりません。 求め方をぜひ教えていただきたいです。 Ａの第１行： ２，４，０ Ａの第２行： ２，６，４ Ａの第３行： ３，９，８ わかる方計算過程を教えてください。
- ベストアンサー
- 数学・算数
- 逆行列(LU分解)を求める数値的なテクニック
ある行列の逆行列をLU分解で求めるプログラムを使用しています。その行列の成分の大きさの最大値と最小値の差が10の30乗ほどあります。コンピュータで計算する場合、極端に大きな数字や小さな数字のまま計算すると正しく計算できないことがあります。AA-1=Eなので、行列の各成分をX倍すれば、求めたい逆行列A-1のX倍の逆行列が求まります。これをX分の1にすれば、求めたい逆行列を求めることができます。 すなわち、各成分の値がコンピュータで処理するのに適していない場合、ある種の補正によって、計算が可能となり、出てきた逆行列はその補正とは反対の操作をすることで求めたい逆行列を求めることができます。 そこで質問です。 コンピュータで逆行列を計算するのに適した行列の成分の値の大きさはいくつ程度でしょうか。 またその適した値にするテクニックとしてはどのようなものがあるでしょうか。
- ベストアンサー
- 数学・算数
- LU分解について
こんにちは。 行列のＬＵ分解のｃによるプログラムを勉強している者です。 ＬＵ分解のｃによるプログラムについてお尋ねしたいことがございます。 ニューメリカルレシピに、ＬＵ分解の方法として 以下のようなコードが記載せれておりました。 void ludcmp(double **a, int n, int *indx, double *d) { int i,imax,j,k; double big,dum,sum,temp; double *vv; // 各行の暗黙のスケーリングを記録する． vv=vector(n); *d=1.0; // まだ行交換していない． for (i=0;i<n;i++) { // 行についてループし，暗黙のスケーリングの情報を得る． big=0.0; for (j=0;j<n;j++) if ((temp=fabs(a[i][j])) > big) big=temp; if (big == 0.0) printf("Singular matrix in routine ludcmp\n"); // 最大要素が０なら特異行列である． vv[i]=1.0/big; // スケーリングを記録する． } for (j=0;j<n;j++) { // Crout法，列についてのループ for (i=0;i<j;i++) { // 方程式(2.3.12)のi=j以外 sum=a[i][j]; for (k=0;k<i;k++) sum -= a[i][k]*a[k][j]; a[i][j]=sum; } big=0.0; for (i=j;i<n;i++) { sum=a[i][j]; for (k=0;k<j;k++) sum -= a[i][k]*a[k][j]; a[i][j]=sum; if ( (dum=vv[i]*fabs(sum)) >= big) { big=dum; imax=i; } } if (j != imax) { for (k=0;k<n;k++) { dum=a[imax][k]; a[imax][k]=a[j][k]; a[j][k]=dum; } *d = -(*d); vv[imax]=vv[j]; } indx[j]=imax; if (a[j][j] == 0.0) a[j][j]=TINY; if (j != n) { dum=1.0/(a[j][j]); for (i=j+1;i<n;i++) a[i][j] *= dum; } } free_vector(vv); } 例えば、ＬＵ分解させたい行列３×３の（｛4,7,6｝｛2,5,9｝｛3,1,8｝）だとしたら、 メイン関数にどのように書いたらＬＵ分解を実行できるでしょうか？ ニューメリカルレシピを買って読んでみたものの、よく理解できず、質問させていただきます。 稚拙な質問かもしれませんが、どうぞよろしくお願いいたします。
- ベストアンサー
- C・C++・C#
お礼
DreaMMaster様 本当にご親切な解説ありがとうございました。 とてもよくわかりました。 結局のところＬＵ分解する「意味」が わたしにわかっていなかったようです。 どうもありがとうございました。