This can be considered as the continuation from
the first one. This time, the Gauss-Jordan elimination
is used. However, the full-fledged version of
the elimination is not needed.
This elimination comes into play due to the fact
that if the elements in lower triangle of the
matrix are all zeroes, the determinant is simply
the product of the diagonals. This also applies
in the upper triangle zeroes. It means:
a ? ? ? a 0 0 0
0 b ? ? or ? b 0 0 's determinant = a*b*c*d
0 0 c ? ? ? c 0
0 0 0 d ? ? ? d
? = don't care.
So, the elimination works to shape the matrix into
one of these forms. My solution aims for the lower
triangle zeroes. Full Gauss-Jordan elimination
requires both upper and lower triangles be zeroes.
I don't want to explain this in details as you can
easily find it in good algebra books. I just give
you a debug procedure "printarr", so that you can
inject a line at mark {1} in order to learn how the
algorithm works. If you'd like to have a more
detailed look, put it at mark {2} instead.
Try it using this matrix:
1 0 -1 0
0 3 0 2
-1 0 2 1
0 5 0 7
and then this one to understand when to "re-sort":
0 0 -1 0
0 3 0 2
-1 0 2 1
0 5 0 7
Here is my solution:
----------------------------------------------------
type { you can modify this dimension }
arr = array[1..4, 1..4] of real;
procedure printarr(new_a: arr; dim: integer);
var i, j: integer;
begin
for i := 1 to dim do
begin
for j := 1 to dim do
write(new_a[i,j]:5:3,' ');
writeln;
end;
end;
{ Gauss-Jordan Elimination }
function gauss(a: arr; dim: integer): real;
var
new_a : arr;
i, j, k : integer;
factor, temp, det : real;
begin
{ copy the array }
for i := 1 to dim do
for j := 1 to dim do
new_a[i,j] := a[i,j];
det := 1.0;
{ do the elimination }
for i := 1 to dim-1 do
begin
{ if the main diagonal value is zero }
{ re-sort the array }
if (new_a[i,i] = 0) then
begin
for j := i+1 to dim do
begin
if (new_a[j,i] <> 0) then
begin
for k := 1 to dim do
begin
temp := new_a[i,k];
new_a[i,k] := new_a[j,k];
new_a[j,k] := temp;
end;
{ For Gauss-Jordan Elimination, }
{ if we do a switch, the determinant }
{ switches sign. }
det := -det;
break;
end;
end;
end;
{ if after the resorting, the value is still zero }
{ then the determinant is definitely zero }
if (new_a[i,i] = 0) then
begin
gauss := 0; exit;
end;
{ eliminate the lower rows to achieve triangular zeroes }
for j := i+1 to dim do
begin
if (new_a[j,i] <> 0) then
begin
factor := (new_a[j,i] * 1.0) / new_a[i,i];
for k := i to dim do
begin
new_a[j,k] := new_a[j,k] - factor * new_a[i,k];
{2}
end;
{1}
end;
end;
end;
{ calculate the main diagonal }
for i := 1 to dim do
det := det * new_a[i,i];
gauss := det;
end;
----------------------------------------------------
Roby Joehanes, © 2001