diff --git a/cmake/lapack.cmake b/cmake/lapack.cmake
index d339f0ce9..ce5d0831f 100644
--- a/cmake/lapack.cmake
+++ b/cmake/lapack.cmake
@@ -187,7 +187,7 @@ set(CLASRC
cposv.f cposvx.f cpotrf2.f cpotri.f cpstrf.f cpstf2.f
cppcon.f cppequ.f cpprfs.f cppsv.f cppsvx.f cpptrf.f cpptri.f cpptrs.f
cptcon.f cpteqr.f cptrfs.f cptsv.f cptsvx.f cpttrf.f cpttrs.f cptts2.f
- crot.f cspcon.f csprfs.f cspsv.f
+ crot.f crscl.f cspcon.f csprfs.f cspsv.f
cspsvx.f csptrf.f csptri.f csptrs.f csrscl.f cstedc.f
cstegr.f cstein.f csteqr.f csycon.f
csyrfs.f csysv.f csysvx.f csytf2.f csytrf.f csytri.f
@@ -381,7 +381,7 @@ set(ZLASRC
zposv.f zposvx.f zpotrf2.f zpotri.f zpotrs.f zpstrf.f zpstf2.f
zppcon.f zppequ.f zpprfs.f zppsv.f zppsvx.f zpptrf.f zpptri.f zpptrs.f
zptcon.f zpteqr.f zptrfs.f zptsv.f zptsvx.f zpttrf.f zpttrs.f zptts2.f
- zrot.f zspcon.f zsprfs.f zspsv.f
+ zrot.f zrscl.f zspcon.f zsprfs.f zspsv.f
zspsvx.f zsptrf.f zsptri.f zsptrs.f zdrscl.f zstedc.f
zstegr.f zstein.f zsteqr.f zsycon.f
zsyrfs.f zsysv.f zsysvx.f zsytf2.f zsytrf.f zsytri.f
diff --git a/lapack-netlib/SRC/Makefile b/lapack-netlib/SRC/Makefile
index 74db14e46..c75fd5f49 100644
--- a/lapack-netlib/SRC/Makefile
+++ b/lapack-netlib/SRC/Makefile
@@ -280,7 +280,7 @@ CLASRC_O = \
cposv.o cposvx.o cpotf2.o cpotri.o cpstrf.o cpstf2.o \
cppcon.o cppequ.o cpprfs.o cppsv.o cppsvx.o cpptrf.o cpptri.o cpptrs.o \
cptcon.o cpteqr.o cptrfs.o cptsv.o cptsvx.o cpttrf.o cpttrs.o cptts2.o \
- crot.o cspcon.o cspmv.o cspr.o csprfs.o cspsv.o \
+ crot.o crscl.o cspcon.o cspmv.o cspr.o csprfs.o cspsv.o \
cspsvx.o csptrf.o csptri.o csptrs.o csrscl.o cstedc.o \
cstegr.o cstein.o csteqr.o \
csycon.o csymv.o \
@@ -488,7 +488,7 @@ ZLASRC_O = \
zposv.o zposvx.o zpotf2.o zpotrf.o zpotri.o zpotrs.o zpstrf.o zpstf2.o \
zppcon.o zppequ.o zpprfs.o zppsv.o zppsvx.o zpptrf.o zpptri.o zpptrs.o \
zptcon.o zpteqr.o zptrfs.o zptsv.o zptsvx.o zpttrf.o zpttrs.o zptts2.o \
- zrot.o zspcon.o zspmv.o zspr.o zsprfs.o zspsv.o \
+ zrot.o zrscl.o zspcon.o zspmv.o zspr.o zsprfs.o zspsv.o \
zspsvx.o zsptrf.o zsptri.o zsptrs.o zdrscl.o zstedc.o \
zstegr.o zstein.o zsteqr.o \
zsycon.o zsymv.o \
diff --git a/lapack-netlib/SRC/cgetf2.f b/lapack-netlib/SRC/cgetf2.f
index aac989970..995ee40ec 100644
--- a/lapack-netlib/SRC/cgetf2.f
+++ b/lapack-netlib/SRC/cgetf2.f
@@ -101,7 +101,7 @@
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*
-*> \ingroup complexGEcomputational
+*> \ingroup getf2
*
* =====================================================================
SUBROUTINE CGETF2( M, N, A, LDA, IPIV, INFO )
@@ -126,16 +126,14 @@
$ ZERO = ( 0.0E+0, 0.0E+0 ) )
* ..
* .. Local Scalars ..
- REAL SFMIN
- INTEGER I, J, JP
+ INTEGER J, JP
* ..
* .. External Functions ..
- REAL SLAMCH
INTEGER ICAMAX
- EXTERNAL SLAMCH, ICAMAX
+ EXTERNAL ICAMAX
* ..
* .. External Subroutines ..
- EXTERNAL CGERU, CSCAL, CSWAP, XERBLA
+ EXTERNAL CGERU, CRSCL, CSWAP, XERBLA
* ..
* .. Intrinsic Functions ..
INTRINSIC MAX, MIN
@@ -161,10 +159,6 @@
*
IF( M.EQ.0 .OR. N.EQ.0 )
$ RETURN
-*
-* Compute machine safe minimum
-*
- SFMIN = SLAMCH('S')
*
DO 10 J = 1, MIN( M, N )
*
@@ -181,15 +175,8 @@
*
* Compute elements J+1:M of J-th column.
*
- IF( J.LT.M ) THEN
- IF( ABS(A( J, J )) .GE. SFMIN ) THEN
- CALL CSCAL( M-J, ONE / A( J, J ), A( J+1, J ), 1 )
- ELSE
- DO 20 I = 1, M-J
- A( J+I, J ) = A( J+I, J ) / A( J, J )
- 20 CONTINUE
- END IF
- END IF
+ IF( J.LT.M )
+ $ CALL CRSCL( M-J, A( J, J ), A( J+1, J ), 1 )
*
ELSE IF( INFO.EQ.0 ) THEN
*
diff --git a/lapack-netlib/SRC/crscl.f b/lapack-netlib/SRC/crscl.f
new file mode 100644
index 000000000..22919cd62
--- /dev/null
+++ b/lapack-netlib/SRC/crscl.f
@@ -0,0 +1,202 @@
+*> \brief \b CRSCL multiplies a vector by the reciprocal of a real scalar.
+*
+* =========== DOCUMENTATION ===========
+*
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
+*
+*> \htmlonly
+*> Download CRSCL + dependencies
+*>
+*> [TGZ]
+*>
+*> [ZIP]
+*>
+*> [TXT]
+*> \endhtmlonly
+*
+* Definition:
+* ===========
+*
+* SUBROUTINE CRSCL( N, A, X, INCX )
+*
+* .. Scalar Arguments ..
+* INTEGER INCX, N
+* COMPLEX A
+* ..
+* .. Array Arguments ..
+* COMPLEX X( * )
+* ..
+*
+*
+*> \par Purpose:
+* =============
+*>
+*> \verbatim
+*>
+*> CRSCL multiplies an n-element complex vector x by the complex scalar
+*> 1/a. This is done without overflow or underflow as long as
+*> the final result x/a does not overflow or underflow.
+*> \endverbatim
+*
+* Arguments:
+* ==========
+*
+*> \param[in] N
+*> \verbatim
+*> N is INTEGER
+*> The number of components of the vector x.
+*> \endverbatim
+*>
+*> \param[in] A
+*> \verbatim
+*> A is COMPLEX
+*> The scalar a which is used to divide each component of x.
+*> A must not be 0, or the subroutine will divide by zero.
+*> \endverbatim
+*>
+*> \param[in,out] X
+*> \verbatim
+*> X is COMPLEX array, dimension
+*> (1+(N-1)*abs(INCX))
+*> The n-element vector x.
+*> \endverbatim
+*>
+*> \param[in] INCX
+*> \verbatim
+*> INCX is INTEGER
+*> The increment between successive values of the vector X.
+*> > 0: X(1) = X(1) and X(1+(i-1)*INCX) = x(i), 1< i<= n
+*> \endverbatim
+*
+* Authors:
+* ========
+*
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
+*
+*> \ingroup complexOTHERauxiliary
+*
+* =====================================================================
+ SUBROUTINE CRSCL( N, A, X, INCX )
+*
+* -- LAPACK auxiliary routine --
+* -- LAPACK is a software package provided by Univ. of Tennessee, --
+* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
+*
+* .. Scalar Arguments ..
+ INTEGER INCX, N
+ COMPLEX A
+* ..
+* .. Array Arguments ..
+ COMPLEX X( * )
+* ..
+*
+* =====================================================================
+*
+* .. Parameters ..
+ REAL ZERO, ONE
+ PARAMETER ( ZERO = 0.0E+0, ONE = 1.0E+0 )
+* ..
+* .. Local Scalars ..
+ REAL SAFMAX, SAFMIN, OV, AR, AI, ABSR, ABSI, UR
+ % , UI
+* ..
+* .. External Functions ..
+ REAL SLAMCH
+ COMPLEX CLADIV
+ EXTERNAL SLAMCH, CLADIV
+* ..
+* .. External Subroutines ..
+ EXTERNAL CSCAL, CSSCAL, CSRSCL
+* ..
+* .. Intrinsic Functions ..
+ INTRINSIC ABS
+* ..
+* .. Executable Statements ..
+*
+* Quick return if possible
+*
+ IF( N.LE.0 )
+ $ RETURN
+*
+* Get machine parameters
+*
+ SAFMIN = SLAMCH( 'S' )
+ SAFMAX = ONE / SAFMIN
+ OV = SLAMCH( 'O' )
+*
+* Initialize constants related to A.
+*
+ AR = REAL( A )
+ AI = AIMAG( A )
+ ABSR = ABS( AR )
+ ABSI = ABS( AI )
+*
+ IF( AI.EQ.ZERO ) THEN
+* If alpha is real, then we can use csrscl
+ CALL CSRSCL( N, AR, X, INCX )
+*
+ ELSE IF( AR.EQ.ZERO ) THEN
+* If alpha has a zero real part, then we follow the same rules as if
+* alpha were real.
+ IF( ABSI.GT.SAFMAX ) THEN
+ CALL CSSCAL( N, SAFMIN, X, INCX )
+ CALL CSCAL( N, CMPLX( ZERO, -SAFMAX / AI ), X, INCX )
+ ELSE IF( ABSI.LT.SAFMIN ) THEN
+ CALL CSCAL( N, CMPLX( ZERO, -SAFMIN / AI ), X, INCX )
+ CALL CSSCAL( N, SAFMAX, X, INCX )
+ ELSE
+ CALL CSCAL( N, CMPLX( ZERO, -ONE / AI ), X, INCX )
+ END IF
+*
+ ELSE
+* The following numbers can be computed.
+* They are the inverse of the real and imaginary parts of 1/alpha.
+* Note that a and b are always different from zero.
+* NaNs are only possible if either:
+* 1. alphaR or alphaI is NaN.
+* 2. alphaR and alphaI are both infinite, in which case it makes sense
+* to propagate a NaN.
+ UR = AR + AI * ( AI / AR )
+ UI = AI + AR * ( AR / AI )
+*
+ IF( (ABS( UR ).LT.SAFMIN).OR.(ABS( UI ).LT.SAFMIN) ) THEN
+* This means that both alphaR and alphaI are very small.
+ CALL CSCAL( N, CMPLX( SAFMIN / UR, -SAFMIN / UI ), X, INCX )
+ CALL CSSCAL( N, SAFMAX, X, INCX )
+ ELSE IF( (ABS( UR ).GT.SAFMAX).OR.(ABS( UI ).GT.SAFMAX) ) THEN
+ IF( (ABSR.GT.OV).OR.(ABSI.GT.OV) ) THEN
+* This means that a and b are both Inf. No need for scaling.
+ CALL CSCAL( N, CMPLX( ONE / UR, -ONE / UI ), X, INCX )
+ ELSE
+ CALL CSSCAL( N, SAFMIN, X, INCX )
+ IF( (ABS( UR ).GT.OV).OR.(ABS( UI ).GT.OV) ) THEN
+* Infs were generated. We do proper scaling to avoid them.
+ IF( ABSR.GE.ABSI ) THEN
+* ABS( UR ) <= ABS( UI )
+ UR = (SAFMIN * AR) + SAFMIN * (AI * ( AI / AR ))
+ UI = (SAFMIN * AI) + AR * ( (SAFMIN * AR) / AI )
+ ELSE
+* ABS( UR ) > ABS( UI )
+ UR = (SAFMIN * AR) + AI * ( (SAFMIN * AI) / AR )
+ UI = (SAFMIN * AI) + SAFMIN * (AR * ( AR / AI ))
+ END IF
+ CALL CSCAL( N, CMPLX( ONE / UR, -ONE / UI ), X, INCX )
+ ELSE
+ CALL CSCAL( N, CMPLX( SAFMAX / UR, -SAFMAX / UI ),
+ $ X, INCX )
+ END IF
+ END IF
+ ELSE
+ CALL CSCAL( N, CMPLX( ONE / UR, -ONE / UI ), X, INCX )
+ END IF
+ END IF
+*
+ RETURN
+*
+* End of CRSCL
+*
+ END
diff --git a/lapack-netlib/SRC/zgetf2.f b/lapack-netlib/SRC/zgetf2.f
index c247f8645..7c63dbbee 100644
--- a/lapack-netlib/SRC/zgetf2.f
+++ b/lapack-netlib/SRC/zgetf2.f
@@ -101,7 +101,7 @@
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*
-*> \ingroup complex16GEcomputational
+*> \ingroup getf2
*
* =====================================================================
SUBROUTINE ZGETF2( M, N, A, LDA, IPIV, INFO )
@@ -127,7 +127,7 @@
* ..
* .. Local Scalars ..
DOUBLE PRECISION SFMIN
- INTEGER I, J, JP
+ INTEGER J, JP
* ..
* .. External Functions ..
DOUBLE PRECISION DLAMCH
@@ -135,7 +135,7 @@
EXTERNAL DLAMCH, IZAMAX
* ..
* .. External Subroutines ..
- EXTERNAL XERBLA, ZGERU, ZSCAL, ZSWAP
+ EXTERNAL XERBLA, ZGERU, ZRSCL, ZSWAP
* ..
* .. Intrinsic Functions ..
INTRINSIC MAX, MIN
@@ -181,15 +181,8 @@
*
* Compute elements J+1:M of J-th column.
*
- IF( J.LT.M ) THEN
- IF( ABS(A( J, J )) .GE. SFMIN ) THEN
- CALL ZSCAL( M-J, ONE / A( J, J ), A( J+1, J ), 1 )
- ELSE
- DO 20 I = 1, M-J
- A( J+I, J ) = A( J+I, J ) / A( J, J )
- 20 CONTINUE
- END IF
- END IF
+ IF( J.LT.M )
+ $ CALL ZRSCL( M-J, A( J, J ), A( J+1, J ), 1 )
*
ELSE IF( INFO.EQ.0 ) THEN
*
diff --git a/lapack-netlib/SRC/zrscl.f b/lapack-netlib/SRC/zrscl.f
new file mode 100644
index 000000000..970f6de75
--- /dev/null
+++ b/lapack-netlib/SRC/zrscl.f
@@ -0,0 +1,203 @@
+*> \brief \b ZDRSCL multiplies a vector by the reciprocal of a real scalar.
+*
+* =========== DOCUMENTATION ===========
+*
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
+*
+*> \htmlonly
+*> Download ZDRSCL + dependencies
+*>
+*> [TGZ]
+*>
+*> [ZIP]
+*>
+*> [TXT]
+*> \endhtmlonly
+*
+* Definition:
+* ===========
+*
+* SUBROUTINE ZRSCL( N, A, X, INCX )
+*
+* .. Scalar Arguments ..
+* INTEGER INCX, N
+* COMPLEX*16 A
+* ..
+* .. Array Arguments ..
+* COMPLEX*16 X( * )
+* ..
+*
+*
+*> \par Purpose:
+* =============
+*>
+*> \verbatim
+*>
+*> ZRSCL multiplies an n-element complex vector x by the complex scalar
+*> 1/a. This is done without overflow or underflow as long as
+*> the final result x/a does not overflow or underflow.
+*> \endverbatim
+*
+* Arguments:
+* ==========
+*
+*> \param[in] N
+*> \verbatim
+*> N is INTEGER
+*> The number of components of the vector x.
+*> \endverbatim
+*>
+*> \param[in] A
+*> \verbatim
+*> A is COMPLEX*16
+*> The scalar a which is used to divide each component of x.
+*> A must not be 0, or the subroutine will divide by zero.
+*> \endverbatim
+*>
+*> \param[in,out] X
+*> \verbatim
+*> X is COMPLEX*16 array, dimension
+*> (1+(N-1)*abs(INCX))
+*> The n-element vector x.
+*> \endverbatim
+*>
+*> \param[in] INCX
+*> \verbatim
+*> INCX is INTEGER
+*> The increment between successive values of the vector SX.
+*> > 0: SX(1) = X(1) and SX(1+(i-1)*INCX) = x(i), 1< i<= n
+*> \endverbatim
+*
+* Authors:
+* ========
+*
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
+*
+*> \ingroup complex16OTHERauxiliary
+*
+* =====================================================================
+ SUBROUTINE ZRSCL( N, A, X, INCX )
+*
+* -- LAPACK auxiliary routine --
+* -- LAPACK is a software package provided by Univ. of Tennessee, --
+* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
+*
+* .. Scalar Arguments ..
+ INTEGER INCX, N
+ COMPLEX*16 A
+* ..
+* .. Array Arguments ..
+ COMPLEX*16 X( * )
+* ..
+*
+* =====================================================================
+*
+* .. Parameters ..
+ DOUBLE PRECISION ZERO, ONE
+ PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 )
+* ..
+* .. Local Scalars ..
+ DOUBLE PRECISION SAFMAX, SAFMIN, OV, AR, AI, ABSR, ABSI, UR, UI
+* ..
+* .. External Functions ..
+ DOUBLE PRECISION DLAMCH
+ COMPLEX*16 ZLADIV
+ EXTERNAL DLAMCH, ZLADIV
+* ..
+* .. External Subroutines ..
+ EXTERNAL DSCAL, ZDSCAL, ZDRSCL
+* ..
+* .. Intrinsic Functions ..
+ INTRINSIC ABS
+* ..
+* .. Executable Statements ..
+*
+* Quick return if possible
+*
+ IF( N.LE.0 )
+ $ RETURN
+*
+* Get machine parameters
+*
+ SAFMIN = DLAMCH( 'S' )
+ SAFMAX = ONE / SAFMIN
+ OV = DLAMCH( 'O' )
+*
+* Initialize constants related to A.
+*
+ AR = DBLE( A )
+ AI = DIMAG( A )
+ ABSR = ABS( AR )
+ ABSI = ABS( AI )
+*
+ IF( AI.EQ.ZERO ) THEN
+* If alpha is real, then we can use csrscl
+ CALL ZDRSCL( N, AR, X, INCX )
+*
+ ELSE IF( AR.EQ.ZERO ) THEN
+* If alpha has a zero real part, then we follow the same rules as if
+* alpha were real.
+ IF( ABSI.GT.SAFMAX ) THEN
+ CALL ZDSCAL( N, SAFMIN, X, INCX )
+ CALL ZSCAL( N, DCMPLX( ZERO, -SAFMAX / AI ), X, INCX )
+ ELSE IF( ABSI.LT.SAFMIN ) THEN
+ CALL ZSCAL( N, DCMPLX( ZERO, -SAFMIN / AI ), X, INCX )
+ CALL ZDSCAL( N, SAFMAX, X, INCX )
+ ELSE
+ CALL ZSCAL( N, DCMPLX( ZERO, -ONE / AI ), X, INCX )
+ END IF
+*
+ ELSE
+* The following numbers can be computed.
+* They are the inverse of the real and imaginary parts of 1/alpha.
+* Note that a and b are always different from zero.
+* NaNs are only possible if either:
+* 1. alphaR or alphaI is NaN.
+* 2. alphaR and alphaI are both infinite, in which case it makes sense
+* to propagate a NaN.
+ UR = AR + AI * ( AI / AR )
+ UI = AI + AR * ( AR / AI )
+*
+ IF( (ABS( UR ).LT.SAFMIN).OR.(ABS( UI ).LT.SAFMIN) ) THEN
+* This means that both alphaR and alphaI are very small.
+ CALL ZSCAL( N, DCMPLX( SAFMIN / UR, -SAFMIN / UI ), X,
+ $ INCX )
+ CALL ZDSCAL( N, SAFMAX, X, INCX )
+ ELSE IF( (ABS( UR ).GT.SAFMAX).OR.(ABS( UI ).GT.SAFMAX) ) THEN
+ IF( (ABSR.GT.OV).OR.(ABSI.GT.OV) ) THEN
+* This means that a and b are both Inf. No need for scaling.
+ CALL ZSCAL( N, DCMPLX( ONE / UR, -ONE / UI ), X, INCX )
+ ELSE
+ CALL ZDSCAL( N, SAFMIN, X, INCX )
+ IF( (ABS( UR ).GT.OV).OR.(ABS( UI ).GT.OV) ) THEN
+* Infs were generated. We do proper scaling to avoid them.
+ IF( ABSR.GE.ABSI ) THEN
+* ABS( UR ) <= ABS( UI )
+ UR = (SAFMIN * AR) + SAFMIN * (AI * ( AI / AR ))
+ UI = (SAFMIN * AI) + AR * ( (SAFMIN * AR) / AI )
+ ELSE
+* ABS( UR ) > ABS( UI )
+ UR = (SAFMIN * AR) + AI * ( (SAFMIN * AI) / AR )
+ UI = (SAFMIN * AI) + SAFMIN * (AR * ( AR / AI ))
+ END IF
+ CALL ZSCAL( N, DCMPLX( ONE / UR, -ONE / UI ), X,
+ $ INCX )
+ ELSE
+ CALL ZSCAL( N, DCMPLX( SAFMAX / UR, -SAFMAX / UI ),
+ $ X, INCX )
+ END IF
+ END IF
+ ELSE
+ CALL ZSCAL( N, DCMPLX( ONE / UR, -ONE / UI ), X, INCX )
+ END IF
+ END IF
+*
+ RETURN
+*
+* End of ZRSCL
+*
+ END