ガウス・ジョルダン法で連立方程式を解こうとしています。下記のコードでは係数行列は単位行列になるのですが解が正しくありません。どこがおかしいのでしょう? ちょっとしたところをミスしていると思うのですが、どうしてもわかりません。
var
A: array [1..3, 1..4] of Extended = ( ( 3, 2, 1,10),
( 1, 1, 2, 9),
( 2, 3, 1,11) );
procedure GaussJordan(N: Integer);
var
i, j, k: integer;
pivot,erase: Extended;
begin
for i := 1 to N do
begin
pivot := A[i][i]; //ピボット
for j := 1 to N+1 do
A[i][j] := A[i][j]/pivot; //ピボットで割る
for k := 1 to N do //掃き出し処理
begin
erase := A[k][i]; //消去する式の係数
for j := i to N do
if k <> i then
A[k][j] := A[k][j]-erase*A[i][j];
end;
end;
end;
正直計算式は全然分かってません。
で、力業で手元のCのソースと比べてみると
pivot := A[i][i]; //ピボット
> for j := 1 to N+1 do
for j :=i + 1 to N+1 do
こんな感じになってました。
以下の列ループがまちがっていました。やはり単純ミスでした(笑)。お騒がせいたしました。
erase := A[k][i]; //消去する式の係数
//for j := i to N do
for j := i to N+1 do
まるさんの例はソースを全部見ていないので何とも言えませんがpvotの扱いや行ループと列ループの回し方が違っていると思います。