LU分解（ピボット選択なし仕様）
<![CDATA[
program LU(input,output);
const n=3;
type matrix=array[1..n,1..n] of real; {n$B9T(Bn$BNs$N9TNs$r:n@.(B}
vector=array[1..n] of real; {n$BNs$r:n@.(B}
var a : matrix;
b, M : vector;
procedure make_matrix(var a : matrix); { $B9TNs(BA$B$NMWAG(B }
begin
a[1,1] := 1; a[1,2] := -4; a[1,3] := 5;
a[2,1] := -2; a[2,2] := 10; a[2,3] := -6;
a[3,1] := -1; a[3,2] := 6; a[3,3] := -2
end; { make_matrix }
procedure make_vector(var b : vector); { b$B%Y%/%H%k$NMWAG(B }
begin
b[1] := 0;
b[2] := 0;
b[3] := 1
end; { make_vector }
procedure decomp(var a:matrix; var b:vector); { $BA0?J>C5n(B }
var i,j,k:integer;
begin
for k:=1 to n-1 do
for j:=k+1 to n do
begin
for i:=k+1 to n do
begin
M[j]:=a[j,k]/a[k,k];
a[j,i]:=a[j,i]-M[j]*a[k,i]
end;
b[j]:=b[j]-M[j]*b[k]
end
end; { decomp }
procedure solve(var a:matrix; var b:vector); { $B8eB`BeF~(B }
var i,j : integer;
B:real;
begin
for i:=n downto 1 do
begin
B:=0;
for j:=i+1 to n do
B:=B-a[i,j]*b[j];
b[i]:=(b[i]+B)/a[i,i]
end;
end; { solve }
procedure print_vector(b:vector); { $BO"N)#1<!J}Dx<0$N2r(B }
var i : integer;
begin
for i:=1 to n do
writeln(i:2, b[i]:10:5); { $BI=<($9$k@:EY$O>.?tE@0J2<#57e$^$G=PMh$k(B }
end; { print_vector }
begin
make_matrix(a);
make_vector(b);
decomp(a,b);
solve(a,b);
print_vector(b);
end.
]]>
Last modified: 2000$BG/(B1$B7n(B20$BF|(B 13:51