LU分解(ピボット選択仕様)
program LU_Pbot(input,output);{ $B%T%\%C%HA*Br;HMQ(B }
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] := -2; a[1,3] := -1;
a[2,1] := 2; a[2,2] := -4; a[2,3] := 0;
a[3,1] := -1; a[3,2] := 1; a[3,3] := 5
end; { make_matrix }
procedure make_vector(var b : vector); { b$B%Y%/%H%k$NMWAG(B }
begin
b[1] := -6;
b[2] := -6;
b[3] := 16
end; { make_vector }
procedure Pbot_c(var a:matrix; var b:vector; j:integer); { $B%T%\%C%HA*Br(B }
var i,I:integer;
L,box_a,box_b:real;
begin
i:=j;
L:=abs(a[j,j]);
while (i<=n) do begin
if (L<=abs(a[i,j])) then begin { $B@dBPCM$,:GBg$N$b$N$r8+$D$1$k(B }
L:=abs(a[i,j]);
I:=i;
end;
i:=i+1;
end;
for i:=j to n do begin { step0.5 }
box_a:=a[j,i];
a[j,i]:=a[I,i];
a[I,i]:=box_a
end;
box_b:=b[j];
b[j]:=b[I];
b[I]:=box_b;
end; { Pbot_c }
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
begin
Pbot_c(a,b,k);
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
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.?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