diff mbox

[Ada] Fix Forward_Eliminate routine to allow use with complex matrices

Message ID 20111013105308.GA5329@adacore.com
State New
Headers show

Commit Message

Arnaud Charlet Oct. 13, 2011, 10:53 a.m. UTC 
Use proper "abs" function returning a real for comparing magnitude of
elements. The previous local implementation using "-" only was correct
for real values. This prepares for the pure Ada reimplementation of
Generic_Complex_Arrays.

Tested on x86_64-pc-linux-gnu, committed on trunk

2011-10-13  Geert Bosch  <bosch@adacore.com>

	* s-gearop.ads (Forward_Eliminate): Add "abs" formal function
	returning a Real.
	* s-gearop.adb (Forward_Eliminate): Remove local "abs" function
	and use formal.
	* a-ngrear.adb (Forward_Eliminate): Adjust instantiation for
	new profile.
diff mbox

Patch

Index: a-ngrear.adb
===================================================================
--- a-ngrear.adb	(revision 179909)
+++ a-ngrear.adb	(working copy)
@@ -33,7 +33,7 @@ 
 --  reason for this is new Ada 2012 requirements that prohibit algorithms such
 --  as Strassen's algorithm, which may be used by some BLAS implementations. In
 --  addition, some platforms lacked suitable compilers to compile the reference
---  BLAS/LAPACK implementation. Finally, on many platforms there may be more
+--  BLAS/LAPACK implementation. Finally, on some platforms there are be more
 --  floating point types than supported by BLAS/LAPACK.
 
 with Ada.Containers.Generic_Anonymous_Array_Sort; use Ada.Containers;
@@ -59,6 +59,7 @@ 
 
    procedure Forward_Eliminate is new Ops.Forward_Eliminate
     (Scalar        => Real'Base,
+     Real          => Real'Base,
      Matrix        => Real_Matrix,
      Zero          => 0.0,
      One           => 1.0);
Index: s-gearop.adb
===================================================================
--- s-gearop.adb	(revision 179909)
+++ s-gearop.adb	(working copy)
@@ -161,9 +161,6 @@ 
       pragma Assert (M'First (1) = N'First (1) and then
                      M'Last  (1) = N'Last (1));
 
-      function "abs" (X : Scalar) return Scalar is
-        (if X < Zero then Zero - X else X);
-
       --  The following are variations of the elementary matrix row operations:
       --  row switching, row multiplication and row addition. Because in this
       --  algorithm the addition factor is always a negated value, we chose to
@@ -274,14 +271,14 @@ 
       for J in M'Range (2) loop
          declare
             Max_Row : Integer := Row;
-            Max_Abs : Scalar := Zero;
+            Max_Abs : Real'Base := 0.0;
 
          begin
             --  Find best pivot in column J, starting in row Row
 
             for K in Row .. M'Last (1) loop
                declare
-                  New_Abs : constant Scalar := abs M (K, J);
+                  New_Abs : constant Real'Base := abs M (K, J);
                begin
                   if Max_Abs < New_Abs then
                      Max_Abs := New_Abs;
@@ -290,7 +287,7 @@ 
                end;
             end loop;
 
-            if Zero < Max_Abs then
+            if Max_Abs > 0.0 then
                Switch_Row (M, N, Row, Max_Row);
                Divide_Row (M, N, Row, M (Row, J));
 
Index: s-gearop.ads
===================================================================
--- s-gearop.ads	(revision 179909)
+++ s-gearop.ads	(working copy)
@@ -65,12 +65,14 @@ 
 
    generic
       type Scalar is private;
+      type Real is digits <>;
       type Matrix is array (Integer range <>, Integer range <>) of Scalar;
+      with function "abs" (Right : Scalar) return Real'Base is <>;
       with function "-" (Left, Right : Scalar) return Scalar is <>;
       with function "*" (Left, Right : Scalar) return Scalar is <>;
       with function "/" (Left, Right : Scalar) return Scalar is <>;
-      with function "<" (Left, Right : Scalar) return Boolean is <>;
-      Zero, One : Scalar;
+      Zero : Scalar;
+      One  : Scalar;
    procedure Forward_Eliminate
      (M   : in out Matrix;
       N   : in out Matrix;