After three irrelevant commit ( r706, r707, r708), this is the first relevant

modifications. I believe that in the complex case, when NORMX = ZERO but ALPHI
is not ZERO (i.e. ALPHA has a complex imaginary part), then TAU needs to be
such that:

            " ( 1 - TAU ) * ( ALPHA ) = ABS( ALPHA ) "

Since we have
            ALPHR = REAL( ALPHA )
            ALPHI = AIMAG( ALPHA )
            XNORM = DLAPY2( ALPHR, ALPHI )
The way to do this is to set TAU with

            TAU = CMPLX( ONE - ALPHR / XNORM, ALPHI / XNORM )

as opposed to

            TAU = CMPLX( ONE - ALPHR / XNORM, -ALPHI / XNORM )

(Note: XNORM is used as temporary variable here)

diff --git a/SRC/clarfp.f b/SRC/clarfp.f
index 22435db..7eb21ee 100644 (file)
--- a/SRC/clarfp.f
+++ b/SRC/clarfp.f
@@ -109,7 +109,7 @@
          ELSE
 *           Only "reflecting" the diagonal entry to be real and non-negative.
             XNORM = SLAPY2( ALPHR, ALPHI )
-            TAU = CMPLX( ONE - ALPHR / XNORM, -ALPHI / XNORM )
+            TAU = CMPLX( ONE - ALPHR / XNORM, ALPHI / XNORM )
             DO J = 1, N-1
                X( 1 + (J-1)*INCX ) = ZERO
             END DO

diff --git a/SRC/zlarfp.f b/SRC/zlarfp.f
index ac4e9d4..b34083b 100644 (file)
--- a/SRC/zlarfp.f
+++ b/SRC/zlarfp.f
@@ -109,7 +109,7 @@
          ELSE
 *           Only "reflecting" the diagonal entry to be real and non-negative.
             XNORM = DLAPY2( ALPHR, ALPHI )
-            TAU = DCMPLX( ONE - ALPHR / XNORM, -ALPHI / XNORM )
+            TAU = DCMPLX( ONE - ALPHR / XNORM, ALPHI / XNORM )
             DO J = 1, N-1
                X( 1 + (J-1)*INCX ) = ZERO
             END DO