with Standard_Floating_Numbers_io; use Standard_Floating_Numbers_io;
with Standard_Natural_Vectors;
with Standard_Natural_Vectors_io; use Standard_Natural_Vectors_io;
with Standard_Floating_Linear_Solvers; use Standard_Floating_Linear_Solvers;
package body Minor_Computations is
function Determinant
( mat : Matrix; rows,cols : Standard_Natural_Vectors.Vector )
return double_float is
-- DESCRIPTION :
-- Computes the determinant of the matrix selecting d rows and columns,
-- where d = rows'length = cols'length.
res : double_float := 1.0;
sqm : Matrix(rows'range,rows'range);
piv : Standard_Natural_Vectors.Vector(rows'range);
inf : natural;
begin
for i in rows'range loop
piv(i) := i;
for j in cols'range loop
sqm(i,j) := mat(rows(i),cols(j));
end loop;
end loop;
lufac(sqm,rows'last,piv,inf);
for i in rows'range loop
res := res*sqm(i,i);
end loop;
for i in piv'range loop
if piv(i) > i
then res := -res;
end if;
end loop;
return res;
end Determinant;
procedure Minors ( file : in file_type;
n,m : in natural; mat : in Matrix ) is
procedure Minors ( d : in natural ) is
rows,cols : Standard_Natural_Vectors.Vector(1..d);
procedure Select_Rows ( k,start : in natural ) is
det : double_float;
begin
if k > d
then det := Determinant(mat,rows,cols);
put(file,"Minor ");
put(file,rows); put(file," X"); put(file,cols);
put(file," equals ");
put(file,det); new_line(file);
else for j in start..n loop
rows(k) := j;
Select_Rows(k+1,j+1);
end loop;
end if;
end Select_Rows;
procedure Select_Columns ( k,start : in natural ) is
begin
if k > d
then Select_Rows(1,1);
else for j in start..m loop
cols(k) := j;
Select_Columns(k+1,j+1);
end loop;
end if;
end Select_Columns;
begin
Select_Columns(1,1);
end Minors;
begin
for d in 2..m loop
Minors(d);
end loop;
end Minors;
function Number_of_Minors ( n,m : natural ) return natural is
res : natural := 0;
procedure Count_Minors ( d : in natural ) is
procedure Select_Rows ( k,start : in natural ) is
begin
if k > d
then res := res + 1;
else for j in start..n loop
Select_Rows(k+1,j+1);
end loop;
end if;
end Select_Rows;
procedure Select_Columns ( k,start : in natural ) is
begin
if k > d
then Select_Rows(1,1);
else for j in start..m loop
Select_Columns(k+1,j+1);
end loop;
end if;
end Select_Columns;
begin
Select_Columns(1,1);
end Count_Minors;
begin
for d in 2..m loop
Count_Minors(d);
end loop;
return res;
end Number_of_Minors;
function Sign_of_Minors ( n,m : natural; mat : Matrix ) return Vector is
res : Vector(1..Number_of_Minors(n,m));
tol : constant double_float := 10.0**(-10);
cnt : natural := 0;
procedure Compute_Minors ( d : in natural ) is
rows,cols : Standard_Natural_Vectors.Vector(1..d);
procedure Select_Rows ( k,start : in natural ) is
det : double_float;
begin
if k > d
then det := Determinant(mat,rows,cols);
cnt := cnt + 1;
if abs(det) < tol
then res(cnt) := 0;
elsif det > 0.0
then res(cnt) := +1;
else res(cnt) := -1;
end if;
else for j in start..n loop
rows(k) := j;
Select_Rows(k+1,j+1);
end loop;
end if;
end Select_Rows;
procedure Select_Columns ( k,start : in natural ) is
begin
if k > d
then Select_Rows(1,1);
else for j in start..m loop
cols(k) := j;
Select_Columns(k+1,j+1);
end loop;
end if;
end Select_Columns;
begin
Select_Columns(1,1);
end Compute_Minors;
begin
for d in 2..m loop
Compute_Minors(d);
end loop;
return res;
end Sign_of_Minors;
end Minor_Computations;