procedure GAUSSJ(VAR A:matrx2; N:integer; VAR B:array of real);
var
IPIV,INDXR,INDXC:array[1..50] of integer;
J,I,K,L,LL:integer;
BIG,PIVINV,DUM:real; IROW,ICOL:integer;
begin
For J:=1 To N do
IPIV[J]:=0;
For I:=1 To N do
begin
BIG:=0;
For J:=1 To N do
begin
If IPIV[J] <> 1 Then
begin
For K:=1 To N do
begin
If IPIV[K] = 0 Then
begin
If Abs(A[J, K]) >= BIG Then
begin
BIG:=Abs(A[J, K]);
IROW:=J;
ICOL:=K;
end;
end
Else if IPIV[K] > 1 Then
ShowMessage('Singular matrix.');
end;
end;
end;
IPIV[ICOL]:=IPIV[ICOL] + 1;
If IROW <> ICOL Then
begin
For L:=1 To N do
begin
DUM:=A[IROW, L];
A[IROW, L]:=A[ICOL, L];
A[ICOL, L]:=DUM;
end;
DUM:=B[IROW];
B[IROW]:=B[ICOL];
B[ICOL]:=DUM;
end;
INDXR[I]:=IROW;
INDXC[I]:=ICOL;
If A[ICOL, ICOL] = 0 Then ShowMessage('Singular matrix.');
PIVINV:=1 / A[ICOL, ICOL];
A[ICOL, ICOL]:=1;
For L:=1 To N do
A[ICOL, L]:=A[ICOL, L] * PIVINV;
B[ICOL]:=B[ICOL] * PIVINV;
For LL:=1 To N do
begin
If LL <> ICOL Then
begin
DUM:=A[LL, ICOL];
A[LL, ICOL]:=0;
For L:=1 To N do
A[LL, L]:=A[LL, L] - A[ICOL, L] * DUM;
B[LL]:=B[LL] - B[ICOL] * DUM;
end;
end;
end;
For L:=N DownTo 1 do
begin
If INDXR[L] <> INDXC[L] Then
begin
For K:=1 To N do
begin
DUM:=A[K, INDXR[L]];
A[K, INDXR[L]]:=A[K, INDXC[L]];
A[K, INDXC[L]]:=DUM;
end;
end;
end;
end;