diff --git a/lapack-netlib/SRC/slaqps.f b/lapack-netlib/SRC/slaqps.f index 9c62ec8b6..3f8af304f 100644 --- a/lapack-netlib/SRC/slaqps.f +++ b/lapack-netlib/SRC/slaqps.f @@ -127,7 +127,7 @@ *> \param[in,out] AUXV *> \verbatim *> AUXV is REAL array, dimension (NB) -*> Auxiliar vector. +*> Auxiliary vector. *> \endverbatim *> *> \param[in,out] F diff --git a/lapack-netlib/SRC/slaqr0.f b/lapack-netlib/SRC/slaqr0.f index 1dcd3d176..318b46943 100644 --- a/lapack-netlib/SRC/slaqr0.f +++ b/lapack-netlib/SRC/slaqr0.f @@ -67,7 +67,7 @@ *> \param[in] N *> \verbatim *> N is INTEGER -*> The order of the matrix H. N .GE. 0. +*> The order of the matrix H. N >= 0. *> \endverbatim *> *> \param[in] ILO @@ -79,12 +79,12 @@ *> \verbatim *> IHI is INTEGER *> It is assumed that H is already upper triangular in rows -*> and columns 1:ILO-1 and IHI+1:N and, if ILO.GT.1, +*> and columns 1:ILO-1 and IHI+1:N and, if ILO > 1, *> H(ILO,ILO-1) is zero. ILO and IHI are normally set by a *> previous call to SGEBAL, and then passed to SGEHRD when the *> matrix output by SGEBAL is reduced to Hessenberg form. *> Otherwise, ILO and IHI should be set to 1 and N, -*> respectively. If N.GT.0, then 1.LE.ILO.LE.IHI.LE.N. +*> respectively. If N > 0, then 1 <= ILO <= IHI <= N. *> If N = 0, then ILO = 1 and IHI = 0. *> \endverbatim *> @@ -97,19 +97,19 @@ *> decomposition (the Schur form); 2-by-2 diagonal blocks *> (corresponding to complex conjugate pairs of eigenvalues) *> are returned in standard form, with H(i,i) = H(i+1,i+1) -*> and H(i+1,i)*H(i,i+1).LT.0. If INFO = 0 and WANTT is +*> and H(i+1,i)*H(i,i+1) < 0. If INFO = 0 and WANTT is *> .FALSE., then the contents of H are unspecified on exit. -*> (The output value of H when INFO.GT.0 is given under the +*> (The output value of H when INFO > 0 is given under the *> description of INFO below.) *> -*> This subroutine may explicitly set H(i,j) = 0 for i.GT.j and +*> This subroutine may explicitly set H(i,j) = 0 for i > j and *> j = 1, 2, ... ILO-1 or j = IHI+1, IHI+2, ... N. *> \endverbatim *> *> \param[in] LDH *> \verbatim *> LDH is INTEGER -*> The leading dimension of the array H. LDH .GE. max(1,N). +*> The leading dimension of the array H. LDH >= max(1,N). *> \endverbatim *> *> \param[out] WR @@ -125,7 +125,7 @@ *> and WI(ILO:IHI). If two eigenvalues are computed as a *> complex conjugate pair, they are stored in consecutive *> elements of WR and WI, say the i-th and (i+1)th, with -*> WI(i) .GT. 0 and WI(i+1) .LT. 0. If WANTT is .TRUE., then +*> WI(i) > 0 and WI(i+1) < 0. If WANTT is .TRUE., then *> the eigenvalues are stored in the same order as on the *> diagonal of the Schur form returned in H, with *> WR(i) = H(i,i) and, if H(i:i+1,i:i+1) is a 2-by-2 diagonal @@ -143,7 +143,7 @@ *> IHIZ is INTEGER *> Specify the rows of Z to which transformations must be *> applied if WANTZ is .TRUE.. -*> 1 .LE. ILOZ .LE. ILO; IHI .LE. IHIZ .LE. N. +*> 1 <= ILOZ <= ILO; IHI <= IHIZ <= N. *> \endverbatim *> *> \param[in,out] Z @@ -153,7 +153,7 @@ *> If WANTZ is .TRUE., then Z(ILO:IHI,ILOZ:IHIZ) is *> replaced by Z(ILO:IHI,ILOZ:IHIZ)*U where U is the *> orthogonal Schur factor of H(ILO:IHI,ILO:IHI). -*> (The output value of Z when INFO.GT.0 is given under +*> (The output value of Z when INFO > 0 is given under *> the description of INFO below.) *> \endverbatim *> @@ -161,7 +161,7 @@ *> \verbatim *> LDZ is INTEGER *> The leading dimension of the array Z. if WANTZ is .TRUE. -*> then LDZ.GE.MAX(1,IHIZ). Otherwize, LDZ.GE.1. +*> then LDZ >= MAX(1,IHIZ). Otherwise, LDZ >= 1. *> \endverbatim *> *> \param[out] WORK @@ -174,7 +174,7 @@ *> \param[in] LWORK *> \verbatim *> LWORK is INTEGER -*> The dimension of the array WORK. LWORK .GE. max(1,N) +*> The dimension of the array WORK. LWORK >= max(1,N) *> is sufficient, but LWORK typically as large as 6*N may *> be required for optimal performance. A workspace query *> to determine the optimal workspace size is recommended. @@ -190,19 +190,19 @@ *> \param[out] INFO *> \verbatim *> INFO is INTEGER -*> = 0: successful exit -*> .GT. 0: if INFO = i, SLAQR0 failed to compute all of +*> = 0: successful exit +*> > 0: if INFO = i, SLAQR0 failed to compute all of *> the eigenvalues. Elements 1:ilo-1 and i+1:n of WR *> and WI contain those eigenvalues which have been *> successfully computed. (Failures are rare.) *> -*> If INFO .GT. 0 and WANT is .FALSE., then on exit, +*> If INFO > 0 and WANT is .FALSE., then on exit, *> the remaining unconverged eigenvalues are the eigen- *> values of the upper Hessenberg matrix rows and *> columns ILO through INFO of the final, output *> value of H. *> -*> If INFO .GT. 0 and WANTT is .TRUE., then on exit +*> If INFO > 0 and WANTT is .TRUE., then on exit *> *> (*) (initial value of H)*U = U*(final value of H) *> @@ -210,7 +210,7 @@ *> value of H is upper Hessenberg and quasi-triangular *> in rows and columns INFO+1 through IHI. *> -*> If INFO .GT. 0 and WANTZ is .TRUE., then on exit +*> If INFO > 0 and WANTZ is .TRUE., then on exit *> *> (final value of Z(ILO:IHI,ILOZ:IHIZ) *> = (initial value of Z(ILO:IHI,ILOZ:IHIZ)*U @@ -218,7 +218,7 @@ *> where U is the orthogonal matrix in (*) (regard- *> less of the value of WANTT.) *> -*> If INFO .GT. 0 and WANTZ is .FALSE., then Z is not +*> If INFO > 0 and WANTZ is .FALSE., then Z is not *> accessed. *> \endverbatim * @@ -677,7 +677,7 @@ END IF END IF * -* ==== Use up to NS of the the smallest magnatiude +* ==== Use up to NS of the the smallest magnitude * . shifts. If there aren't NS shifts available, * . then use them all, possibly dropping one to * . make the number of shifts even. ==== diff --git a/lapack-netlib/SRC/slaqr1.f b/lapack-netlib/SRC/slaqr1.f index 2de33849d..6bb88c794 100644 --- a/lapack-netlib/SRC/slaqr1.f +++ b/lapack-netlib/SRC/slaqr1.f @@ -69,7 +69,7 @@ *> \verbatim *> LDH is INTEGER *> The leading dimension of H as declared in -*> the calling procedure. LDH.GE.N +*> the calling procedure. LDH >= N *> \endverbatim *> *> \param[in] SR1 diff --git a/lapack-netlib/SRC/slaqr2.f b/lapack-netlib/SRC/slaqr2.f index 8e1f34910..f4f8ca7f2 100644 --- a/lapack-netlib/SRC/slaqr2.f +++ b/lapack-netlib/SRC/slaqr2.f @@ -103,7 +103,7 @@ *> \param[in] NW *> \verbatim *> NW is INTEGER -*> Deflation window size. 1 .LE. NW .LE. (KBOT-KTOP+1). +*> Deflation window size. 1 <= NW <= (KBOT-KTOP+1). *> \endverbatim *> *> \param[in,out] H @@ -121,7 +121,7 @@ *> \verbatim *> LDH is INTEGER *> Leading dimension of H just as declared in the calling -*> subroutine. N .LE. LDH +*> subroutine. N <= LDH *> \endverbatim *> *> \param[in] ILOZ @@ -133,7 +133,7 @@ *> \verbatim *> IHIZ is INTEGER *> Specify the rows of Z to which transformations must be -*> applied if WANTZ is .TRUE.. 1 .LE. ILOZ .LE. IHIZ .LE. N. +*> applied if WANTZ is .TRUE.. 1 <= ILOZ <= IHIZ <= N. *> \endverbatim *> *> \param[in,out] Z @@ -149,7 +149,7 @@ *> \verbatim *> LDZ is INTEGER *> The leading dimension of Z just as declared in the -*> calling subroutine. 1 .LE. LDZ. +*> calling subroutine. 1 <= LDZ. *> \endverbatim *> *> \param[out] NS @@ -194,13 +194,13 @@ *> \verbatim *> LDV is INTEGER *> The leading dimension of V just as declared in the -*> calling subroutine. NW .LE. LDV +*> calling subroutine. NW <= LDV *> \endverbatim *> *> \param[in] NH *> \verbatim *> NH is INTEGER -*> The number of columns of T. NH.GE.NW. +*> The number of columns of T. NH >= NW. *> \endverbatim *> *> \param[out] T @@ -212,14 +212,14 @@ *> \verbatim *> LDT is INTEGER *> The leading dimension of T just as declared in the -*> calling subroutine. NW .LE. LDT +*> calling subroutine. NW <= LDT *> \endverbatim *> *> \param[in] NV *> \verbatim *> NV is INTEGER *> The number of rows of work array WV available for -*> workspace. NV.GE.NW. +*> workspace. NV >= NW. *> \endverbatim *> *> \param[out] WV @@ -231,7 +231,7 @@ *> \verbatim *> LDWV is INTEGER *> The leading dimension of W just as declared in the -*> calling subroutine. NW .LE. LDV +*> calling subroutine. NW <= LDV *> \endverbatim *> *> \param[out] WORK diff --git a/lapack-netlib/SRC/slaqr3.f b/lapack-netlib/SRC/slaqr3.f index 534e2c489..ccad338b9 100644 --- a/lapack-netlib/SRC/slaqr3.f +++ b/lapack-netlib/SRC/slaqr3.f @@ -100,7 +100,7 @@ *> \param[in] NW *> \verbatim *> NW is INTEGER -*> Deflation window size. 1 .LE. NW .LE. (KBOT-KTOP+1). +*> Deflation window size. 1 <= NW <= (KBOT-KTOP+1). *> \endverbatim *> *> \param[in,out] H @@ -118,7 +118,7 @@ *> \verbatim *> LDH is INTEGER *> Leading dimension of H just as declared in the calling -*> subroutine. N .LE. LDH +*> subroutine. N <= LDH *> \endverbatim *> *> \param[in] ILOZ @@ -130,7 +130,7 @@ *> \verbatim *> IHIZ is INTEGER *> Specify the rows of Z to which transformations must be -*> applied if WANTZ is .TRUE.. 1 .LE. ILOZ .LE. IHIZ .LE. N. +*> applied if WANTZ is .TRUE.. 1 <= ILOZ <= IHIZ <= N. *> \endverbatim *> *> \param[in,out] Z @@ -146,7 +146,7 @@ *> \verbatim *> LDZ is INTEGER *> The leading dimension of Z just as declared in the -*> calling subroutine. 1 .LE. LDZ. +*> calling subroutine. 1 <= LDZ. *> \endverbatim *> *> \param[out] NS @@ -191,13 +191,13 @@ *> \verbatim *> LDV is INTEGER *> The leading dimension of V just as declared in the -*> calling subroutine. NW .LE. LDV +*> calling subroutine. NW <= LDV *> \endverbatim *> *> \param[in] NH *> \verbatim *> NH is INTEGER -*> The number of columns of T. NH.GE.NW. +*> The number of columns of T. NH >= NW. *> \endverbatim *> *> \param[out] T @@ -209,14 +209,14 @@ *> \verbatim *> LDT is INTEGER *> The leading dimension of T just as declared in the -*> calling subroutine. NW .LE. LDT +*> calling subroutine. NW <= LDT *> \endverbatim *> *> \param[in] NV *> \verbatim *> NV is INTEGER *> The number of rows of work array WV available for -*> workspace. NV.GE.NW. +*> workspace. NV >= NW. *> \endverbatim *> *> \param[out] WV @@ -228,7 +228,7 @@ *> \verbatim *> LDWV is INTEGER *> The leading dimension of W just as declared in the -*> calling subroutine. NW .LE. LDV +*> calling subroutine. NW <= LDV *> \endverbatim *> *> \param[out] WORK diff --git a/lapack-netlib/SRC/slaqr4.f b/lapack-netlib/SRC/slaqr4.f index 12b6b2fb1..cd642e07f 100644 --- a/lapack-netlib/SRC/slaqr4.f +++ b/lapack-netlib/SRC/slaqr4.f @@ -74,7 +74,7 @@ *> \param[in] N *> \verbatim *> N is INTEGER -*> The order of the matrix H. N .GE. 0. +*> The order of the matrix H. N >= 0. *> \endverbatim *> *> \param[in] ILO @@ -86,12 +86,12 @@ *> \verbatim *> IHI is INTEGER *> It is assumed that H is already upper triangular in rows -*> and columns 1:ILO-1 and IHI+1:N and, if ILO.GT.1, +*> and columns 1:ILO-1 and IHI+1:N and, if ILO > 1, *> H(ILO,ILO-1) is zero. ILO and IHI are normally set by a *> previous call to SGEBAL, and then passed to SGEHRD when the *> matrix output by SGEBAL is reduced to Hessenberg form. *> Otherwise, ILO and IHI should be set to 1 and N, -*> respectively. If N.GT.0, then 1.LE.ILO.LE.IHI.LE.N. +*> respectively. If N > 0, then 1 <= ILO <= IHI <= N. *> If N = 0, then ILO = 1 and IHI = 0. *> \endverbatim *> @@ -104,19 +104,19 @@ *> decomposition (the Schur form); 2-by-2 diagonal blocks *> (corresponding to complex conjugate pairs of eigenvalues) *> are returned in standard form, with H(i,i) = H(i+1,i+1) -*> and H(i+1,i)*H(i,i+1).LT.0. If INFO = 0 and WANTT is +*> and H(i+1,i)*H(i,i+1) < 0. If INFO = 0 and WANTT is *> .FALSE., then the contents of H are unspecified on exit. -*> (The output value of H when INFO.GT.0 is given under the +*> (The output value of H when INFO > 0 is given under the *> description of INFO below.) *> -*> This subroutine may explicitly set H(i,j) = 0 for i.GT.j and +*> This subroutine may explicitly set H(i,j) = 0 for i > j and *> j = 1, 2, ... ILO-1 or j = IHI+1, IHI+2, ... N. *> \endverbatim *> *> \param[in] LDH *> \verbatim *> LDH is INTEGER -*> The leading dimension of the array H. LDH .GE. max(1,N). +*> The leading dimension of the array H. LDH >= max(1,N). *> \endverbatim *> *> \param[out] WR @@ -132,7 +132,7 @@ *> and WI(ILO:IHI). If two eigenvalues are computed as a *> complex conjugate pair, they are stored in consecutive *> elements of WR and WI, say the i-th and (i+1)th, with -*> WI(i) .GT. 0 and WI(i+1) .LT. 0. If WANTT is .TRUE., then +*> WI(i) > 0 and WI(i+1) < 0. If WANTT is .TRUE., then *> the eigenvalues are stored in the same order as on the *> diagonal of the Schur form returned in H, with *> WR(i) = H(i,i) and, if H(i:i+1,i:i+1) is a 2-by-2 diagonal @@ -150,7 +150,7 @@ *> IHIZ is INTEGER *> Specify the rows of Z to which transformations must be *> applied if WANTZ is .TRUE.. -*> 1 .LE. ILOZ .LE. ILO; IHI .LE. IHIZ .LE. N. +*> 1 <= ILOZ <= ILO; IHI <= IHIZ <= N. *> \endverbatim *> *> \param[in,out] Z @@ -160,7 +160,7 @@ *> If WANTZ is .TRUE., then Z(ILO:IHI,ILOZ:IHIZ) is *> replaced by Z(ILO:IHI,ILOZ:IHIZ)*U where U is the *> orthogonal Schur factor of H(ILO:IHI,ILO:IHI). -*> (The output value of Z when INFO.GT.0 is given under +*> (The output value of Z when INFO > 0 is given under *> the description of INFO below.) *> \endverbatim *> @@ -168,7 +168,7 @@ *> \verbatim *> LDZ is INTEGER *> The leading dimension of the array Z. if WANTZ is .TRUE. -*> then LDZ.GE.MAX(1,IHIZ). Otherwize, LDZ.GE.1. +*> then LDZ >= MAX(1,IHIZ). Otherwise, LDZ >= 1. *> \endverbatim *> *> \param[out] WORK @@ -181,7 +181,7 @@ *> \param[in] LWORK *> \verbatim *> LWORK is INTEGER -*> The dimension of the array WORK. LWORK .GE. max(1,N) +*> The dimension of the array WORK. LWORK >= max(1,N) *> is sufficient, but LWORK typically as large as 6*N may *> be required for optimal performance. A workspace query *> to determine the optimal workspace size is recommended. @@ -199,19 +199,19 @@ *> INFO is INTEGER *> \verbatim *> INFO is INTEGER -*> = 0: successful exit -*> .GT. 0: if INFO = i, SLAQR4 failed to compute all of +*> = 0: successful exit +*> > 0: if INFO = i, SLAQR4 failed to compute all of *> the eigenvalues. Elements 1:ilo-1 and i+1:n of WR *> and WI contain those eigenvalues which have been *> successfully computed. (Failures are rare.) *> -*> If INFO .GT. 0 and WANT is .FALSE., then on exit, +*> If INFO > 0 and WANT is .FALSE., then on exit, *> the remaining unconverged eigenvalues are the eigen- *> values of the upper Hessenberg matrix rows and *> columns ILO through INFO of the final, output *> value of H. *> -*> If INFO .GT. 0 and WANTT is .TRUE., then on exit +*> If INFO > 0 and WANTT is .TRUE., then on exit *> *> (*) (initial value of H)*U = U*(final value of H) *> @@ -219,7 +219,7 @@ *> value of H is upper Hessenberg and triangular in *> rows and columns INFO+1 through IHI. *> -*> If INFO .GT. 0 and WANTZ is .TRUE., then on exit +*> If INFO > 0 and WANTZ is .TRUE., then on exit *> *> (final value of Z(ILO:IHI,ILOZ:IHIZ) *> = (initial value of Z(ILO:IHI,ILOZ:IHIZ)*U @@ -227,7 +227,7 @@ *> where U is the orthogonal matrix in (*) (regard- *> less of the value of WANTT.) *> -*> If INFO .GT. 0 and WANTZ is .FALSE., then Z is not +*> If INFO > 0 and WANTZ is .FALSE., then Z is not *> accessed. *> \endverbatim * @@ -680,7 +680,7 @@ END IF END IF * -* ==== Use up to NS of the the smallest magnatiude +* ==== Use up to NS of the the smallest magnitude * . shifts. If there aren't NS shifts available, * . then use them all, possibly dropping one to * . make the number of shifts even. ==== diff --git a/lapack-netlib/SRC/slaqr5.f b/lapack-netlib/SRC/slaqr5.f index 65278e355..f04ee577e 100644 --- a/lapack-netlib/SRC/slaqr5.f +++ b/lapack-netlib/SRC/slaqr5.f @@ -133,7 +133,7 @@ *> \verbatim *> LDH is INTEGER *> LDH is the leading dimension of H just as declared in the -*> calling procedure. LDH.GE.MAX(1,N). +*> calling procedure. LDH >= MAX(1,N). *> \endverbatim *> *> \param[in] ILOZ @@ -145,7 +145,7 @@ *> \verbatim *> IHIZ is INTEGER *> Specify the rows of Z to which transformations must be -*> applied if WANTZ is .TRUE.. 1 .LE. ILOZ .LE. IHIZ .LE. N +*> applied if WANTZ is .TRUE.. 1 <= ILOZ <= IHIZ <= N *> \endverbatim *> *> \param[in,out] Z @@ -161,7 +161,7 @@ *> \verbatim *> LDZ is INTEGER *> LDA is the leading dimension of Z just as declared in -*> the calling procedure. LDZ.GE.N. +*> the calling procedure. LDZ >= N. *> \endverbatim *> *> \param[out] V @@ -173,7 +173,7 @@ *> \verbatim *> LDV is INTEGER *> LDV is the leading dimension of V as declared in the -*> calling procedure. LDV.GE.3. +*> calling procedure. LDV >= 3. *> \endverbatim *> *> \param[out] U @@ -185,33 +185,14 @@ *> \verbatim *> LDU is INTEGER *> LDU is the leading dimension of U just as declared in the -*> in the calling subroutine. LDU.GE.3*NSHFTS-3. -*> \endverbatim -*> -*> \param[in] NH -*> \verbatim -*> NH is INTEGER -*> NH is the number of columns in array WH available for -*> workspace. NH.GE.1. -*> \endverbatim -*> -*> \param[out] WH -*> \verbatim -*> WH is REAL array, dimension (LDWH,NH) -*> \endverbatim -*> -*> \param[in] LDWH -*> \verbatim -*> LDWH is INTEGER -*> Leading dimension of WH just as declared in the -*> calling procedure. LDWH.GE.3*NSHFTS-3. +*> in the calling subroutine. LDU >= 3*NSHFTS-3. *> \endverbatim *> *> \param[in] NV *> \verbatim *> NV is INTEGER *> NV is the number of rows in WV agailable for workspace. -*> NV.GE.1. +*> NV >= 1. *> \endverbatim *> *> \param[out] WV @@ -223,9 +204,28 @@ *> \verbatim *> LDWV is INTEGER *> LDWV is the leading dimension of WV as declared in the -*> in the calling subroutine. LDWV.GE.NV. +*> in the calling subroutine. LDWV >= NV. *> \endverbatim * +*> \param[in] NH +*> \verbatim +*> NH is INTEGER +*> NH is the number of columns in array WH available for +*> workspace. NH >= 1. +*> \endverbatim +*> +*> \param[out] WH +*> \verbatim +*> WH is REAL array, dimension (LDWH,NH) +*> \endverbatim +*> +*> \param[in] LDWH +*> \verbatim +*> LDWH is INTEGER +*> Leading dimension of WH just as declared in the +*> calling procedure. LDWH >= 3*NSHFTS-3. +*> \endverbatim +*> * Authors: * ======== * diff --git a/lapack-netlib/SRC/slarfb.f b/lapack-netlib/SRC/slarfb.f index c51f69534..d853a54ec 100644 --- a/lapack-netlib/SRC/slarfb.f +++ b/lapack-netlib/SRC/slarfb.f @@ -92,6 +92,8 @@ *> K is INTEGER *> The order of the matrix T (= the number of elementary *> reflectors whose product defines the block reflector). +*> If SIDE = 'L', M >= K >= 0; +*> if SIDE = 'R', N >= K >= 0. *> \endverbatim *> *> \param[in] V diff --git a/lapack-netlib/SRC/slarfx.f b/lapack-netlib/SRC/slarfx.f index 590e99e70..3175068b8 100644 --- a/lapack-netlib/SRC/slarfx.f +++ b/lapack-netlib/SRC/slarfx.f @@ -94,7 +94,7 @@ *> \param[in] LDC *> \verbatim *> LDC is INTEGER -*> The leading dimension of the array C. LDA >= (1,M). +*> The leading dimension of the array C. LDC >= (1,M). *> \endverbatim *> *> \param[out] WORK diff --git a/lapack-netlib/SRC/slarfy.f b/lapack-netlib/SRC/slarfy.f index 340c54413..f9ba011a2 100644 --- a/lapack-netlib/SRC/slarfy.f +++ b/lapack-netlib/SRC/slarfy.f @@ -103,7 +103,7 @@ * *> \date December 2016 * -*> \ingroup single_eig +*> \ingroup realOTHERauxiliary * * ===================================================================== SUBROUTINE SLARFY( UPLO, N, V, INCV, TAU, C, LDC, WORK ) diff --git a/lapack-netlib/SRC/slarrb.f b/lapack-netlib/SRC/slarrb.f index 988e25ff0..ac9d7bc8c 100644 --- a/lapack-netlib/SRC/slarrb.f +++ b/lapack-netlib/SRC/slarrb.f @@ -91,7 +91,7 @@ *> RTOL2 is REAL *> Tolerance for the convergence of the bisection intervals. *> An interval [LEFT,RIGHT] has converged if -*> RIGHT-LEFT.LT.MAX( RTOL1*GAP, RTOL2*MAX(|LEFT|,|RIGHT|) ) +*> RIGHT-LEFT < MAX( RTOL1*GAP, RTOL2*MAX(|LEFT|,|RIGHT|) ) *> where GAP is the (estimated) distance to the nearest *> eigenvalue. *> \endverbatim @@ -117,7 +117,7 @@ *> WGAP is REAL array, dimension (N-1) *> On input, the (estimated) gaps between consecutive *> eigenvalues of L D L^T, i.e., WGAP(I-OFFSET) is the gap between -*> eigenvalues I and I+1. Note that if IFIRST.EQ.ILAST +*> eigenvalues I and I+1. Note that if IFIRST = ILAST *> then WGAP(IFIRST-OFFSET) must be set to ZERO. *> On output, these gaps are refined. *> \endverbatim diff --git a/lapack-netlib/SRC/slarre.f b/lapack-netlib/SRC/slarre.f index ea9b8fcbc..6636235d0 100644 --- a/lapack-netlib/SRC/slarre.f +++ b/lapack-netlib/SRC/slarre.f @@ -150,7 +150,7 @@ *> RTOL2 is REAL *> Parameters for bisection. *> An interval [LEFT,RIGHT] has converged if -*> RIGHT-LEFT.LT.MAX( RTOL1*GAP, RTOL2*MAX(|LEFT|,|RIGHT|) ) +*> RIGHT-LEFT < MAX( RTOL1*GAP, RTOL2*MAX(|LEFT|,|RIGHT|) ) *> \endverbatim *> *> \param[in] SPLTOL diff --git a/lapack-netlib/SRC/slarrj.f b/lapack-netlib/SRC/slarrj.f index a721d0751..fb867595c 100644 --- a/lapack-netlib/SRC/slarrj.f +++ b/lapack-netlib/SRC/slarrj.f @@ -85,7 +85,7 @@ *> RTOL is REAL *> Tolerance for the convergence of the bisection intervals. *> An interval [LEFT,RIGHT] has converged if -*> RIGHT-LEFT.LT.RTOL*MAX(|LEFT|,|RIGHT|). +*> RIGHT-LEFT < RTOL*MAX(|LEFT|,|RIGHT|). *> \endverbatim *> *> \param[in] OFFSET diff --git a/lapack-netlib/SRC/slarrv.f b/lapack-netlib/SRC/slarrv.f index f9e3cf2b9..04519fde8 100644 --- a/lapack-netlib/SRC/slarrv.f +++ b/lapack-netlib/SRC/slarrv.f @@ -149,7 +149,7 @@ *> RTOL2 is REAL *> Parameters for bisection. *> An interval [LEFT,RIGHT] has converged if -*> RIGHT-LEFT.LT.MAX( RTOL1*GAP, RTOL2*MAX(|LEFT|,|RIGHT|) ) +*> RIGHT-LEFT < MAX( RTOL1*GAP, RTOL2*MAX(|LEFT|,|RIGHT|) ) *> \endverbatim *> *> \param[in,out] W diff --git a/lapack-netlib/SRC/slasd7.f b/lapack-netlib/SRC/slasd7.f index 2adaa5ee7..d8775b6c6 100644 --- a/lapack-netlib/SRC/slasd7.f +++ b/lapack-netlib/SRC/slasd7.f @@ -400,7 +400,7 @@ VL( I ) = VLW( IDXI ) 50 CONTINUE * -* Calculate the allowable deflation tolerence +* Calculate the allowable deflation tolerance * EPS = SLAMCH( 'Epsilon' ) TOL = MAX( ABS( ALPHA ), ABS( BETA ) ) diff --git a/lapack-netlib/SRC/slassq.f b/lapack-netlib/SRC/slassq.f index 35b40f07f..d9930a597 100644 --- a/lapack-netlib/SRC/slassq.f +++ b/lapack-netlib/SRC/slassq.f @@ -60,7 +60,7 @@ *> *> \param[in] X *> \verbatim -*> X is REAL array, dimension (N) +*> X is REAL array, dimension (1+(N-1)*INCX) *> The vector for which a scaled sum of squares is computed. *> x( i ) = X( 1 + ( i - 1 )*INCX ), 1 <= i <= n. *> \endverbatim diff --git a/lapack-netlib/SRC/slaswlq.f b/lapack-netlib/SRC/slaswlq.f index 27b5b8067..5eb9cc5f9 100644 --- a/lapack-netlib/SRC/slaswlq.f +++ b/lapack-netlib/SRC/slaswlq.f @@ -1,3 +1,4 @@ +*> \brief \b SLASWLQ * * Definition: * =========== @@ -18,9 +19,20 @@ *> *> \verbatim *> -*> SLASWLQ computes a blocked Short-Wide LQ factorization of a -*> M-by-N matrix A, where N >= M: -*> A = L * Q +*> SLASWLQ computes a blocked Tall-Skinny LQ factorization of +*> a real M-by-N matrix A for M <= N: +*> +*> A = ( L 0 ) * Q, +*> +*> where: +*> +*> Q is a n-by-N orthogonal matrix, stored on exit in an implicit +*> form in the elements above the digonal of the array A and in +*> the elemenst of the array T; +*> L is an lower-triangular M-by-M matrix stored on exit in +*> the elements on and below the diagonal of the array A. +*> 0 is a M-by-(N-M) zero matrix, if M < N, and is not stored. +*> *> \endverbatim * * Arguments: @@ -150,10 +162,10 @@ SUBROUTINE SLASWLQ( M, N, MB, NB, A, LDA, T, LDT, WORK, LWORK, $ INFO) * -* -- LAPACK computational routine (version 3.8.0) -- +* -- LAPACK computational routine (version 3.9.0) -- * -- LAPACK is a software package provided by Univ. of Tennessee, -- * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. -- -* November 2017 +* November 2019 * * .. Scalar Arguments .. INTEGER INFO, LDA, M, N, MB, NB, LWORK, LDT diff --git a/lapack-netlib/SRC/slasyf_aa.f b/lapack-netlib/SRC/slasyf_aa.f index ed4ef6291..76f632602 100644 --- a/lapack-netlib/SRC/slasyf_aa.f +++ b/lapack-netlib/SRC/slasyf_aa.f @@ -284,8 +284,9 @@ * * Swap A(I1, I2+1:M) with A(I2, I2+1:M) * - CALL SSWAP( M-I2, A( J1+I1-1, I2+1 ), LDA, - $ A( J1+I2-1, I2+1 ), LDA ) + IF( I2.LT.M ) + $ CALL SSWAP( M-I2, A( J1+I1-1, I2+1 ), LDA, + $ A( J1+I2-1, I2+1 ), LDA ) * * Swap A(I1, I1) with A(I2,I2) * @@ -325,13 +326,15 @@ * Compute L(J+2, J+1) = WORK( 3:M ) / T(J, J+1), * where A(J, J+1) = T(J, J+1) and A(J+2:M, J) = L(J+2:M, J+1) * - IF( A( K, J+1 ).NE.ZERO ) THEN - ALPHA = ONE / A( K, J+1 ) - CALL SCOPY( M-J-1, WORK( 3 ), 1, A( K, J+2 ), LDA ) - CALL SSCAL( M-J-1, ALPHA, A( K, J+2 ), LDA ) - ELSE - CALL SLASET( 'Full', 1, M-J-1, ZERO, ZERO, - $ A( K, J+2 ), LDA) + IF( J.LT.(M-1) ) THEN + IF( A( K, J+1 ).NE.ZERO ) THEN + ALPHA = ONE / A( K, J+1 ) + CALL SCOPY( M-J-1, WORK( 3 ), 1, A( K, J+2 ), LDA ) + CALL SSCAL( M-J-1, ALPHA, A( K, J+2 ), LDA ) + ELSE + CALL SLASET( 'Full', 1, M-J-1, ZERO, ZERO, + $ A( K, J+2 ), LDA) + END IF END IF END IF J = J + 1 @@ -432,8 +435,9 @@ * * Swap A(I2+1:M, I1) with A(I2+1:M, I2) * - CALL SSWAP( M-I2, A( I2+1, J1+I1-1 ), 1, - $ A( I2+1, J1+I2-1 ), 1 ) + IF( I2.LT.M ) + $ CALL SSWAP( M-I2, A( I2+1, J1+I1-1 ), 1, + $ A( I2+1, J1+I2-1 ), 1 ) * * Swap A(I1, I1) with A(I2, I2) * @@ -473,13 +477,15 @@ * Compute L(J+2, J+1) = WORK( 3:M ) / T(J, J+1), * where A(J, J+1) = T(J, J+1) and A(J+2:M, J) = L(J+2:M, J+1) * - IF( A( J+1, K ).NE.ZERO ) THEN - ALPHA = ONE / A( J+1, K ) - CALL SCOPY( M-J-1, WORK( 3 ), 1, A( J+2, K ), 1 ) - CALL SSCAL( M-J-1, ALPHA, A( J+2, K ), 1 ) - ELSE - CALL SLASET( 'Full', M-J-1, 1, ZERO, ZERO, - $ A( J+2, K ), LDA ) + IF( J.LT.(M-1) ) THEN + IF( A( J+1, K ).NE.ZERO ) THEN + ALPHA = ONE / A( J+1, K ) + CALL SCOPY( M-J-1, WORK( 3 ), 1, A( J+2, K ), 1 ) + CALL SSCAL( M-J-1, ALPHA, A( J+2, K ), 1 ) + ELSE + CALL SLASET( 'Full', M-J-1, 1, ZERO, ZERO, + $ A( J+2, K ), LDA ) + END IF END IF END IF J = J + 1 diff --git a/lapack-netlib/SRC/slasyf_rk.f b/lapack-netlib/SRC/slasyf_rk.f index b1b37177f..c16708365 100644 --- a/lapack-netlib/SRC/slasyf_rk.f +++ b/lapack-netlib/SRC/slasyf_rk.f @@ -321,7 +321,7 @@ * of A and working backwards, and compute the matrix W = U12*D * for use in updating A11 * -* Initilize the first entry of array E, where superdiagonal +* Initialize the first entry of array E, where superdiagonal * elements of D are stored * E( 1 ) = ZERO @@ -649,7 +649,7 @@ * of A and working forwards, and compute the matrix W = L21*D * for use in updating A22 * -* Initilize the unused last entry of the subdiagonal array E. +* Initialize the unused last entry of the subdiagonal array E. * E( N ) = ZERO * diff --git a/lapack-netlib/SRC/slatdf.f b/lapack-netlib/SRC/slatdf.f index 5496f9db4..495d32502 100644 --- a/lapack-netlib/SRC/slatdf.f +++ b/lapack-netlib/SRC/slatdf.f @@ -85,7 +85,7 @@ *> RHS is REAL array, dimension N. *> On entry, RHS contains contributions from other subsystems. *> On exit, RHS contains the solution of the subsystem with -*> entries acoording to the value of IJOB (see above). +*> entries according to the value of IJOB (see above). *> \endverbatim *> *> \param[in,out] RDSUM @@ -260,7 +260,7 @@ * * Solve for U-part, look-ahead for RHS(N) = +-1. This is not done * in BSOLVE and will hopefully give us a better estimate because -* any ill-conditioning of the original matrix is transfered to U +* any ill-conditioning of the original matrix is transferred to U * and not to L. U(N, N) is an approximation to sigma_min(LU). * CALL SCOPY( N-1, RHS, 1, XP, 1 ) diff --git a/lapack-netlib/SRC/slatsqr.f b/lapack-netlib/SRC/slatsqr.f index d6d682799..b56b0d41e 100644 --- a/lapack-netlib/SRC/slatsqr.f +++ b/lapack-netlib/SRC/slatsqr.f @@ -1,3 +1,4 @@ +*> \brief \b SLATSQR * * Definition: * =========== @@ -19,8 +20,22 @@ *> \verbatim *> *> SLATSQR computes a blocked Tall-Skinny QR factorization of -*> an M-by-N matrix A, where M >= N: -*> A = Q * R . +*> a real M-by-N matrix A for M >= N: +*> +*> A = Q * ( R ), +*> ( 0 ) +*> +*> where: +*> +*> Q is a M-by-M orthogonal matrix, stored on exit in an implicit +*> form in the elements below the digonal of the array A and in +*> the elemenst of the array T; +*> +*> R is an upper-triangular N-by-N matrix, stored on exit in +*> the elements on and above the diagonal of the array A. +*> +*> 0 is a (M-N)-by-N zero matrix, and is not stored. +*> *> \endverbatim * * Arguments: @@ -149,10 +164,10 @@ SUBROUTINE SLATSQR( M, N, MB, NB, A, LDA, T, LDT, WORK, $ LWORK, INFO) * -* -- LAPACK computational routine (version 3.7.0) -- +* -- LAPACK computational routine (version 3.9.0) -- * -- LAPACK is a software package provided by Univ. of Tennessee, -- * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. -- -* December 2016 +* November 2019 * * .. Scalar Arguments .. INTEGER INFO, LDA, M, N, MB, NB, LDT, LWORK diff --git a/lapack-netlib/SRC/sorgtsqr.f b/lapack-netlib/SRC/sorgtsqr.f new file mode 100644 index 000000000..748760d63 --- /dev/null +++ b/lapack-netlib/SRC/sorgtsqr.f @@ -0,0 +1,306 @@ +*> \brief \b SORGTSQR +* +* =========== DOCUMENTATION =========== +* +* Online html documentation available at +* http://www.netlib.org/lapack/explore-html/ +* +*> \htmlonly +*> Download SORGTSQR + dependencies +*> +*> [TGZ] +*> +*> [ZIP] +*> +*> [TXT] +*> +* Definition: +* =========== +* +* SUBROUTINE SORGTSQR( M, N, MB, NB, A, LDA, T, LDT, WORK, LWORK, +* $ INFO ) +* +* .. Scalar Arguments .. +* INTEGER INFO, LDA, LDT, LWORK, M, N, MB, NB +* .. +* .. Array Arguments .. +* REAL A( LDA, * ), T( LDT, * ), WORK( * ) +* .. +* +*> \par Purpose: +* ============= +*> +*> \verbatim +*> +*> SORGTSQR generates an M-by-N real matrix Q_out with orthonormal columns, +*> which are the first N columns of a product of real orthogonal +*> matrices of order M which are returned by SLATSQR +*> +*> Q_out = first_N_columns_of( Q(1)_in * Q(2)_in * ... * Q(k)_in ). +*> +*> See the documentation for SLATSQR. +*> \endverbatim +* +* Arguments: +* ========== +* +*> \param[in] M +*> \verbatim +*> M is INTEGER +*> The number of rows of the matrix A. M >= 0. +*> \endverbatim +*> +*> \param[in] N +*> \verbatim +*> N is INTEGER +*> The number of columns of the matrix A. M >= N >= 0. +*> \endverbatim +*> +*> \param[in] MB +*> \verbatim +*> MB is INTEGER +*> The row block size used by SLATSQR to return +*> arrays A and T. MB > N. +*> (Note that if MB > M, then M is used instead of MB +*> as the row block size). +*> \endverbatim +*> +*> \param[in] NB +*> \verbatim +*> NB is INTEGER +*> The column block size used by SLATSQR to return +*> arrays A and T. NB >= 1. +*> (Note that if NB > N, then N is used instead of NB +*> as the column block size). +*> \endverbatim +*> +*> \param[in,out] A +*> \verbatim +*> A is REAL array, dimension (LDA,N) +*> +*> On entry: +*> +*> The elements on and above the diagonal are not accessed. +*> The elements below the diagonal represent the unit +*> lower-trapezoidal blocked matrix V computed by SLATSQR +*> that defines the input matrices Q_in(k) (ones on the +*> diagonal are not stored) (same format as the output A +*> below the diagonal in SLATSQR). +*> +*> On exit: +*> +*> The array A contains an M-by-N orthonormal matrix Q_out, +*> i.e the columns of A are orthogonal unit vectors. +*> \endverbatim +*> +*> \param[in] LDA +*> \verbatim +*> LDA is INTEGER +*> The leading dimension of the array A. LDA >= max(1,M). +*> \endverbatim +*> +*> \param[in] T +*> \verbatim +*> T is REAL array, +*> dimension (LDT, N * NIRB) +*> where NIRB = Number_of_input_row_blocks +*> = MAX( 1, CEIL((M-N)/(MB-N)) ) +*> Let NICB = Number_of_input_col_blocks +*> = CEIL(N/NB) +*> +*> The upper-triangular block reflectors used to define the +*> input matrices Q_in(k), k=(1:NIRB*NICB). The block +*> reflectors are stored in compact form in NIRB block +*> reflector sequences. Each of NIRB block reflector sequences +*> is stored in a larger NB-by-N column block of T and consists +*> of NICB smaller NB-by-NB upper-triangular column blocks. +*> (same format as the output T in SLATSQR). +*> \endverbatim +*> +*> \param[in] LDT +*> \verbatim +*> LDT is INTEGER +*> The leading dimension of the array T. +*> LDT >= max(1,min(NB1,N)). +*> \endverbatim +*> +*> \param[out] WORK +*> \verbatim +*> (workspace) REAL array, dimension (MAX(2,LWORK)) +*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK. +*> \endverbatim +*> +*> \param[in] LWORK +*> \verbatim +*> The dimension of the array WORK. LWORK >= (M+NB)*N. +*> If LWORK = -1, then a workspace query is assumed. +*> The routine only calculates the optimal size of the WORK +*> array, returns this value as the first entry of the WORK +*> array, and no error message related to LWORK is issued +*> by XERBLA. +*> \endverbatim +*> +*> \param[out] INFO +*> \verbatim +*> INFO is INTEGER +*> = 0: successful exit +*> < 0: if INFO = -i, the i-th argument had an illegal value +*> \endverbatim +*> +* Authors: +* ======== +* +*> \author Univ. of Tennessee +*> \author Univ. of California Berkeley +*> \author Univ. of Colorado Denver +*> \author NAG Ltd. +* +*> \date November 2019 +* +*> \ingroup singleOTHERcomputational +* +*> \par Contributors: +* ================== +*> +*> \verbatim +*> +*> November 2019, Igor Kozachenko, +*> Computer Science Division, +*> University of California, Berkeley +*> +*> \endverbatim +* +* ===================================================================== + SUBROUTINE SORGTSQR( M, N, MB, NB, A, LDA, T, LDT, WORK, LWORK, + $ INFO ) + IMPLICIT NONE +* +* -- LAPACK computational routine (version 3.9.0) -- +* -- LAPACK is a software package provided by Univ. of Tennessee, -- +* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- +* November 2019 +* +* .. Scalar Arguments .. + INTEGER INFO, LDA, LDT, LWORK, M, N, MB, NB +* .. +* .. Array Arguments .. + REAL A( LDA, * ), T( LDT, * ), WORK( * ) +* .. +* +* ===================================================================== +* +* .. Parameters .. + REAL ONE, ZERO + PARAMETER ( ONE = 1.0E+0, ZERO = 0.0E+0 ) +* .. +* .. Local Scalars .. + LOGICAL LQUERY + INTEGER IINFO, LDC, LWORKOPT, LC, LW, NBLOCAL, J +* .. +* .. External Subroutines .. + EXTERNAL SCOPY, SLAMTSQR, SLASET, XERBLA +* .. +* .. Intrinsic Functions .. + INTRINSIC REAL, MAX, MIN +* .. +* .. Executable Statements .. +* +* Test the input parameters +* + LQUERY = LWORK.EQ.-1 + INFO = 0 + IF( M.LT.0 ) THEN + INFO = -1 + ELSE IF( N.LT.0 .OR. M.LT.N ) THEN + INFO = -2 + ELSE IF( MB.LE.N ) THEN + INFO = -3 + ELSE IF( NB.LT.1 ) THEN + INFO = -4 + ELSE IF( LDA.LT.MAX( 1, M ) ) THEN + INFO = -6 + ELSE IF( LDT.LT.MAX( 1, MIN( NB, N ) ) ) THEN + INFO = -8 + ELSE +* +* Test the input LWORK for the dimension of the array WORK. +* This workspace is used to store array C(LDC, N) and WORK(LWORK) +* in the call to DLAMTSQR. See the documentation for DLAMTSQR. +* + IF( LWORK.LT.2 .AND. (.NOT.LQUERY) ) THEN + INFO = -10 + ELSE +* +* Set block size for column blocks +* + NBLOCAL = MIN( NB, N ) +* +* LWORK = -1, then set the size for the array C(LDC,N) +* in DLAMTSQR call and set the optimal size of the work array +* WORK(LWORK) in DLAMTSQR call. +* + LDC = M + LC = LDC*N + LW = N * NBLOCAL +* + LWORKOPT = LC+LW +* + IF( ( LWORK.LT.MAX( 1, LWORKOPT ) ).AND.(.NOT.LQUERY) ) THEN + INFO = -10 + END IF + END IF +* + END IF +* +* Handle error in the input parameters and return workspace query. +* + IF( INFO.NE.0 ) THEN + CALL XERBLA( 'SORGTSQR', -INFO ) + RETURN + ELSE IF ( LQUERY ) THEN + WORK( 1 ) = REAL( LWORKOPT ) + RETURN + END IF +* +* Quick return if possible +* + IF( MIN( M, N ).EQ.0 ) THEN + WORK( 1 ) = REAL( LWORKOPT ) + RETURN + END IF +* +* (1) Form explicitly the tall-skinny M-by-N left submatrix Q1_in +* of M-by-M orthogonal matrix Q_in, which is implicitly stored in +* the subdiagonal part of input array A and in the input array T. +* Perform by the following operation using the routine DLAMTSQR. +* +* Q1_in = Q_in * ( I ), where I is a N-by-N identity matrix, +* ( 0 ) 0 is a (M-N)-by-N zero matrix. +* +* (1a) Form M-by-N matrix in the array WORK(1:LDC*N) with ones +* on the diagonal and zeros elsewhere. +* + CALL SLASET( 'F', M, N, ZERO, ONE, WORK, LDC ) +* +* (1b) On input, WORK(1:LDC*N) stores ( I ); +* ( 0 ) +* +* On output, WORK(1:LDC*N) stores Q1_in. +* + CALL SLAMTSQR( 'L', 'N', M, N, N, MB, NBLOCAL, A, LDA, T, LDT, + $ WORK, LDC, WORK( LC+1 ), LW, IINFO ) +* +* (2) Copy the result from the part of the work array (1:M,1:N) +* with the leading dimension LDC that starts at WORK(1) into +* the output array A(1:M,1:N) column-by-column. +* + DO J = 1, N + CALL SCOPY( M, WORK( (J-1)*LDC + 1 ), 1, A( 1, J ), 1 ) + END DO +* + WORK( 1 ) = REAL( LWORKOPT ) + RETURN +* +* End of SORGTSQR +* + END \ No newline at end of file diff --git a/lapack-netlib/SRC/sorhr_col.f b/lapack-netlib/SRC/sorhr_col.f new file mode 100644 index 000000000..38976245c --- /dev/null +++ b/lapack-netlib/SRC/sorhr_col.f @@ -0,0 +1,439 @@ +*> \brief \b SORHR_COL +* +* =========== DOCUMENTATION =========== +* +* Online html documentation available at +* http://www.netlib.org/lapack/explore-html/ +* +*> \htmlonly +*> Download SORHR_COL + dependencies +*> +*> [TGZ] +*> +*> [ZIP] +*> +*> [TXT] +*> +* Definition: +* =========== +* +* SUBROUTINE SORHR_COL( M, N, NB, A, LDA, T, LDT, D, INFO ) +* +* .. Scalar Arguments .. +* INTEGER INFO, LDA, LDT, M, N, NB +* .. +* .. Array Arguments .. +* REAL A( LDA, * ), D( * ), T( LDT, * ) +* .. +* +*> \par Purpose: +* ============= +*> +*> \verbatim +*> +*> SORHR_COL takes an M-by-N real matrix Q_in with orthonormal columns +*> as input, stored in A, and performs Householder Reconstruction (HR), +*> i.e. reconstructs Householder vectors V(i) implicitly representing +*> another M-by-N matrix Q_out, with the property that Q_in = Q_out*S, +*> where S is an N-by-N diagonal matrix with diagonal entries +*> equal to +1 or -1. The Householder vectors (columns V(i) of V) are +*> stored in A on output, and the diagonal entries of S are stored in D. +*> Block reflectors are also returned in T +*> (same output format as SGEQRT). +*> \endverbatim +* +* Arguments: +* ========== +* +*> \param[in] M +*> \verbatim +*> M is INTEGER +*> The number of rows of the matrix A. M >= 0. +*> \endverbatim +*> +*> \param[in] N +*> \verbatim +*> N is INTEGER +*> The number of columns of the matrix A. M >= N >= 0. +*> \endverbatim +*> +*> \param[in] NB +*> \verbatim +*> NB is INTEGER +*> The column block size to be used in the reconstruction +*> of Householder column vector blocks in the array A and +*> corresponding block reflectors in the array T. NB >= 1. +*> (Note that if NB > N, then N is used instead of NB +*> as the column block size.) +*> \endverbatim +*> +*> \param[in,out] A +*> \verbatim +*> A is REAL array, dimension (LDA,N) +*> +*> On entry: +*> +*> The array A contains an M-by-N orthonormal matrix Q_in, +*> i.e the columns of A are orthogonal unit vectors. +*> +*> On exit: +*> +*> The elements below the diagonal of A represent the unit +*> lower-trapezoidal matrix V of Householder column vectors +*> V(i). The unit diagonal entries of V are not stored +*> (same format as the output below the diagonal in A from +*> SGEQRT). The matrix T and the matrix V stored on output +*> in A implicitly define Q_out. +*> +*> The elements above the diagonal contain the factor U +*> of the "modified" LU-decomposition: +*> Q_in - ( S ) = V * U +*> ( 0 ) +*> where 0 is a (M-N)-by-(M-N) zero matrix. +*> \endverbatim +*> +*> \param[in] LDA +*> \verbatim +*> LDA is INTEGER +*> The leading dimension of the array A. LDA >= max(1,M). +*> \endverbatim +*> +*> \param[out] T +*> \verbatim +*> T is REAL array, +*> dimension (LDT, N) +*> +*> Let NOCB = Number_of_output_col_blocks +*> = CEIL(N/NB) +*> +*> On exit, T(1:NB, 1:N) contains NOCB upper-triangular +*> block reflectors used to define Q_out stored in compact +*> form as a sequence of upper-triangular NB-by-NB column +*> blocks (same format as the output T in SGEQRT). +*> The matrix T and the matrix V stored on output in A +*> implicitly define Q_out. NOTE: The lower triangles +*> below the upper-triangular blcoks will be filled with +*> zeros. See Further Details. +*> \endverbatim +*> +*> \param[in] LDT +*> \verbatim +*> LDT is INTEGER +*> The leading dimension of the array T. +*> LDT >= max(1,min(NB,N)). +*> \endverbatim +*> +*> \param[out] D +*> \verbatim +*> D is REAL array, dimension min(M,N). +*> The elements can be only plus or minus one. +*> +*> D(i) is constructed as D(i) = -SIGN(Q_in_i(i,i)), where +*> 1 <= i <= min(M,N), and Q_in_i is Q_in after performing +*> i-1 steps of “modified” Gaussian elimination. +*> See Further Details. +*> \endverbatim +*> +*> \param[out] INFO +*> \verbatim +*> INFO is INTEGER +*> = 0: successful exit +*> < 0: if INFO = -i, the i-th argument had an illegal value +*> \endverbatim +*> +*> \par Further Details: +* ===================== +*> +*> \verbatim +*> +*> The computed M-by-M orthogonal factor Q_out is defined implicitly as +*> a product of orthogonal matrices Q_out(i). Each Q_out(i) is stored in +*> the compact WY-representation format in the corresponding blocks of +*> matrices V (stored in A) and T. +*> +*> The M-by-N unit lower-trapezoidal matrix V stored in the M-by-N +*> matrix A contains the column vectors V(i) in NB-size column +*> blocks VB(j). For example, VB(1) contains the columns +*> V(1), V(2), ... V(NB). NOTE: The unit entries on +*> the diagonal of Y are not stored in A. +*> +*> The number of column blocks is +*> +*> NOCB = Number_of_output_col_blocks = CEIL(N/NB) +*> +*> where each block is of order NB except for the last block, which +*> is of order LAST_NB = N - (NOCB-1)*NB. +*> +*> For example, if M=6, N=5 and NB=2, the matrix V is +*> +*> +*> V = ( VB(1), VB(2), VB(3) ) = +*> +*> = ( 1 ) +*> ( v21 1 ) +*> ( v31 v32 1 ) +*> ( v41 v42 v43 1 ) +*> ( v51 v52 v53 v54 1 ) +*> ( v61 v62 v63 v54 v65 ) +*> +*> +*> For each of the column blocks VB(i), an upper-triangular block +*> reflector TB(i) is computed. These blocks are stored as +*> a sequence of upper-triangular column blocks in the NB-by-N +*> matrix T. The size of each TB(i) block is NB-by-NB, except +*> for the last block, whose size is LAST_NB-by-LAST_NB. +*> +*> For example, if M=6, N=5 and NB=2, the matrix T is +*> +*> T = ( TB(1), TB(2), TB(3) ) = +*> +*> = ( t11 t12 t13 t14 t15 ) +*> ( t22 t24 ) +*> +*> +*> The M-by-M factor Q_out is given as a product of NOCB +*> orthogonal M-by-M matrices Q_out(i). +*> +*> Q_out = Q_out(1) * Q_out(2) * ... * Q_out(NOCB), +*> +*> where each matrix Q_out(i) is given by the WY-representation +*> using corresponding blocks from the matrices V and T: +*> +*> Q_out(i) = I - VB(i) * TB(i) * (VB(i))**T, +*> +*> where I is the identity matrix. Here is the formula with matrix +*> dimensions: +*> +*> Q(i){M-by-M} = I{M-by-M} - +*> VB(i){M-by-INB} * TB(i){INB-by-INB} * (VB(i))**T {INB-by-M}, +*> +*> where INB = NB, except for the last block NOCB +*> for which INB=LAST_NB. +*> +*> ===== +*> NOTE: +*> ===== +*> +*> If Q_in is the result of doing a QR factorization +*> B = Q_in * R_in, then: +*> +*> B = (Q_out*S) * R_in = Q_out * (S * R_in) = O_out * R_out. +*> +*> So if one wants to interpret Q_out as the result +*> of the QR factorization of B, then corresponding R_out +*> should be obtained by R_out = S * R_in, i.e. some rows of R_in +*> should be multiplied by -1. +*> +*> For the details of the algorithm, see [1]. +*> +*> [1] "Reconstructing Householder vectors from tall-skinny QR", +*> G. Ballard, J. Demmel, L. Grigori, M. Jacquelin, H.D. Nguyen, +*> E. Solomonik, J. Parallel Distrib. Comput., +*> vol. 85, pp. 3-31, 2015. +*> \endverbatim +*> +* Authors: +* ======== +* +*> \author Univ. of Tennessee +*> \author Univ. of California Berkeley +*> \author Univ. of Colorado Denver +*> \author NAG Ltd. +* +*> \date November 2019 +* +*> \ingroup singleOTHERcomputational +* +*> \par Contributors: +* ================== +*> +*> \verbatim +*> +*> November 2019, Igor Kozachenko, +*> Computer Science Division, +*> University of California, Berkeley +*> +*> \endverbatim +* +* ===================================================================== + SUBROUTINE SORHR_COL( M, N, NB, A, LDA, T, LDT, D, INFO ) + IMPLICIT NONE +* +* -- LAPACK computational routine (version 3.9.0) -- +* -- LAPACK is a software package provided by Univ. of Tennessee, -- +* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- +* November 2019 +* +* .. Scalar Arguments .. + INTEGER INFO, LDA, LDT, M, N, NB +* .. +* .. Array Arguments .. + REAL A( LDA, * ), D( * ), T( LDT, * ) +* .. +* +* ===================================================================== +* +* .. Parameters .. + REAL ONE, ZERO + PARAMETER ( ONE = 1.0E+0, ZERO = 0.0E+0 ) +* .. +* .. Local Scalars .. + INTEGER I, IINFO, J, JB, JBTEMP1, JBTEMP2, JNB, + $ NPLUSONE +* .. +* .. External Subroutines .. + EXTERNAL SCOPY, SLAORHR_COL_GETRFNP, SSCAL, STRSM, XERBLA +* .. +* .. Intrinsic Functions .. + INTRINSIC MAX, MIN +* .. +* .. Executable Statements .. +* +* Test the input parameters +* + INFO = 0 + IF( M.LT.0 ) THEN + INFO = -1 + ELSE IF( N.LT.0 .OR. N.GT.M ) THEN + INFO = -2 + ELSE IF( NB.LT.1 ) THEN + INFO = -3 + ELSE IF( LDA.LT.MAX( 1, M ) ) THEN + INFO = -5 + ELSE IF( LDT.LT.MAX( 1, MIN( NB, N ) ) ) THEN + INFO = -7 + END IF +* +* Handle error in the input parameters. +* + IF( INFO.NE.0 ) THEN + CALL XERBLA( 'SORHR_COL', -INFO ) + RETURN + END IF +* +* Quick return if possible +* + IF( MIN( M, N ).EQ.0 ) THEN + RETURN + END IF +* +* On input, the M-by-N matrix A contains the orthogonal +* M-by-N matrix Q_in. +* +* (1) Compute the unit lower-trapezoidal V (ones on the diagonal +* are not stored) by performing the "modified" LU-decomposition. +* +* Q_in - ( S ) = V * U = ( V1 ) * U, +* ( 0 ) ( V2 ) +* +* where 0 is an (M-N)-by-N zero matrix. +* +* (1-1) Factor V1 and U. + + CALL SLAORHR_COL_GETRFNP( N, N, A, LDA, D, IINFO ) +* +* (1-2) Solve for V2. +* + IF( M.GT.N ) THEN + CALL STRSM( 'R', 'U', 'N', 'N', M-N, N, ONE, A, LDA, + $ A( N+1, 1 ), LDA ) + END IF +* +* (2) Reconstruct the block reflector T stored in T(1:NB, 1:N) +* as a sequence of upper-triangular blocks with NB-size column +* blocking. +* +* Loop over the column blocks of size NB of the array A(1:M,1:N) +* and the array T(1:NB,1:N), JB is the column index of a column +* block, JNB is the column block size at each step JB. +* + NPLUSONE = N + 1 + DO JB = 1, N, NB +* +* (2-0) Determine the column block size JNB. +* + JNB = MIN( NPLUSONE-JB, NB ) +* +* (2-1) Copy the upper-triangular part of the current JNB-by-JNB +* diagonal block U(JB) (of the N-by-N matrix U) stored +* in A(JB:JB+JNB-1,JB:JB+JNB-1) into the upper-triangular part +* of the current JNB-by-JNB block T(1:JNB,JB:JB+JNB-1) +* column-by-column, total JNB*(JNB+1)/2 elements. +* + JBTEMP1 = JB - 1 + DO J = JB, JB+JNB-1 + CALL SCOPY( J-JBTEMP1, A( JB, J ), 1, T( 1, J ), 1 ) + END DO +* +* (2-2) Perform on the upper-triangular part of the current +* JNB-by-JNB diagonal block U(JB) (of the N-by-N matrix U) stored +* in T(1:JNB,JB:JB+JNB-1) the following operation in place: +* (-1)*U(JB)*S(JB), i.e the result will be stored in the upper- +* triangular part of T(1:JNB,JB:JB+JNB-1). This multiplication +* of the JNB-by-JNB diagonal block U(JB) by the JNB-by-JNB +* diagonal block S(JB) of the N-by-N sign matrix S from the +* right means changing the sign of each J-th column of the block +* U(JB) according to the sign of the diagonal element of the block +* S(JB), i.e. S(J,J) that is stored in the array element D(J). +* + DO J = JB, JB+JNB-1 + IF( D( J ).EQ.ONE ) THEN + CALL SSCAL( J-JBTEMP1, -ONE, T( 1, J ), 1 ) + END IF + END DO +* +* (2-3) Perform the triangular solve for the current block +* matrix X(JB): +* +* X(JB) * (A(JB)**T) = B(JB), where: +* +* A(JB)**T is a JNB-by-JNB unit upper-triangular +* coefficient block, and A(JB)=V1(JB), which +* is a JNB-by-JNB unit lower-triangular block +* stored in A(JB:JB+JNB-1,JB:JB+JNB-1). +* The N-by-N matrix V1 is the upper part +* of the M-by-N lower-trapezoidal matrix V +* stored in A(1:M,1:N); +* +* B(JB) is a JNB-by-JNB upper-triangular right-hand +* side block, B(JB) = (-1)*U(JB)*S(JB), and +* B(JB) is stored in T(1:JNB,JB:JB+JNB-1); +* +* X(JB) is a JNB-by-JNB upper-triangular solution +* block, X(JB) is the upper-triangular block +* reflector T(JB), and X(JB) is stored +* in T(1:JNB,JB:JB+JNB-1). +* +* In other words, we perform the triangular solve for the +* upper-triangular block T(JB): +* +* T(JB) * (V1(JB)**T) = (-1)*U(JB)*S(JB). +* +* Even though the blocks X(JB) and B(JB) are upper- +* triangular, the routine STRSM will access all JNB**2 +* elements of the square T(1:JNB,JB:JB+JNB-1). Therefore, +* we need to set to zero the elements of the block +* T(1:JNB,JB:JB+JNB-1) below the diagonal before the call +* to STRSM. +* +* (2-3a) Set the elements to zero. +* + JBTEMP2 = JB - 2 + DO J = JB, JB+JNB-2 + DO I = J-JBTEMP2, NB + T( I, J ) = ZERO + END DO + END DO +* +* (2-3b) Perform the triangular solve. +* + CALL STRSM( 'R', 'L', 'T', 'U', JNB, JNB, ONE, + $ A( JB, JB ), LDA, T( 1, JB ), LDT ) +* + END DO +* + RETURN +* +* End of SORHR_COL +* + END \ No newline at end of file diff --git a/lapack-netlib/SRC/sporfsx.f b/lapack-netlib/SRC/sporfsx.f index 52fab6976..ce8c26569 100644 --- a/lapack-netlib/SRC/sporfsx.f +++ b/lapack-netlib/SRC/sporfsx.f @@ -135,7 +135,7 @@ *> \param[in,out] S *> \verbatim *> S is REAL array, dimension (N) -*> The row scale factors for A. If EQUED = 'Y', A is multiplied on +*> The scale factors for A. If EQUED = 'Y', A is multiplied on *> the left and right by diag(S). S is an input argument if FACT = *> 'F'; otherwise, S is an output argument. If FACT = 'F' and EQUED *> = 'Y', each element of S must be positive. If S is output, each @@ -263,7 +263,7 @@ *> information as described below. There currently are up to three *> pieces of information returned for each right-hand side. If *> componentwise accuracy is not requested (PARAMS(3) = 0.0), then -*> ERR_BNDS_COMP is not accessed. If N_ERR_BNDS .LT. 3, then at most +*> ERR_BNDS_COMP is not accessed. If N_ERR_BNDS < 3, then at most *> the first (:,N_ERR_BNDS) entries are returned. *> *> The first index in ERR_BNDS_COMP(i,:) corresponds to the ith @@ -299,14 +299,14 @@ *> \param[in] NPARAMS *> \verbatim *> NPARAMS is INTEGER -*> Specifies the number of parameters set in PARAMS. If .LE. 0, the +*> Specifies the number of parameters set in PARAMS. If <= 0, the *> PARAMS array is never referenced and default values are used. *> \endverbatim *> *> \param[in,out] PARAMS *> \verbatim *> PARAMS is REAL array, dimension NPARAMS -*> Specifies algorithm parameters. If an entry is .LT. 0.0, then +*> Specifies algorithm parameters. If an entry is < 0.0, then *> that entry will be filled with default value used for that *> parameter. Only positions up to NPARAMS are accessed; defaults *> are used for higher-numbered parameters. @@ -314,9 +314,9 @@ *> PARAMS(LA_LINRX_ITREF_I = 1) : Whether to perform iterative *> refinement or not. *> Default: 1.0 -*> = 0.0 : No refinement is performed, and no error bounds are +*> = 0.0: No refinement is performed, and no error bounds are *> computed. -*> = 1.0 : Use the double-precision refinement algorithm, +*> = 1.0: Use the double-precision refinement algorithm, *> possibly with doubled-single computations if the *> compilation environment does not support DOUBLE *> PRECISION. diff --git a/lapack-netlib/SRC/sposvxx.f b/lapack-netlib/SRC/sposvxx.f index 3cdfa749c..fa2c0d3f3 100644 --- a/lapack-netlib/SRC/sposvxx.f +++ b/lapack-netlib/SRC/sposvxx.f @@ -366,7 +366,7 @@ *> information as described below. There currently are up to three *> pieces of information returned for each right-hand side. If *> componentwise accuracy is not requested (PARAMS(3) = 0.0), then -*> ERR_BNDS_COMP is not accessed. If N_ERR_BNDS .LT. 3, then at most +*> ERR_BNDS_COMP is not accessed. If N_ERR_BNDS < 3, then at most *> the first (:,N_ERR_BNDS) entries are returned. *> *> The first index in ERR_BNDS_COMP(i,:) corresponds to the ith @@ -402,14 +402,14 @@ *> \param[in] NPARAMS *> \verbatim *> NPARAMS is INTEGER -*> Specifies the number of parameters set in PARAMS. If .LE. 0, the +*> Specifies the number of parameters set in PARAMS. If <= 0, the *> PARAMS array is never referenced and default values are used. *> \endverbatim *> *> \param[in,out] PARAMS *> \verbatim *> PARAMS is REAL array, dimension NPARAMS -*> Specifies algorithm parameters. If an entry is .LT. 0.0, then +*> Specifies algorithm parameters. If an entry is < 0.0, then *> that entry will be filled with default value used for that *> parameter. Only positions up to NPARAMS are accessed; defaults *> are used for higher-numbered parameters. @@ -417,9 +417,9 @@ *> PARAMS(LA_LINRX_ITREF_I = 1) : Whether to perform iterative *> refinement or not. *> Default: 1.0 -*> = 0.0 : No refinement is performed, and no error bounds are +*> = 0.0: No refinement is performed, and no error bounds are *> computed. -*> = 1.0 : Use the double-precision refinement algorithm, +*> = 1.0: Use the double-precision refinement algorithm, *> possibly with doubled-single computations if the *> compilation environment does not support DOUBLE *> PRECISION. diff --git a/lapack-netlib/SRC/ssb2st_kernels.f b/lapack-netlib/SRC/ssb2st_kernels.f index 54479f89e..08859169b 100644 --- a/lapack-netlib/SRC/ssb2st_kernels.f +++ b/lapack-netlib/SRC/ssb2st_kernels.f @@ -1,26 +1,26 @@ *> \brief \b SSB2ST_KERNELS * * @generated from zhb2st_kernels.f, fortran z -> s, Wed Dec 7 08:22:40 2016 -* +* * =========== DOCUMENTATION =========== * -* Online html documentation available at -* http://www.netlib.org/lapack/explore-html/ +* Online html documentation available at +* http://www.netlib.org/lapack/explore-html/ * *> \htmlonly -*> Download SSB2ST_KERNELS + dependencies -*> -*> [TGZ] -*> -*> [ZIP] -*> +*> Download SSB2ST_KERNELS + dependencies +*> +*> [TGZ] +*> +*> [ZIP] +*> *> [TXT] -*> \endhtmlonly +*> \endhtmlonly * * Definition: * =========== * -* SUBROUTINE SSB2ST_KERNELS( UPLO, WANTZ, TTYPE, +* SUBROUTINE SSB2ST_KERNELS( UPLO, WANTZ, TTYPE, * ST, ED, SWEEP, N, NB, IB, * A, LDA, V, TAU, LDVT, WORK) * @@ -32,9 +32,9 @@ * INTEGER TTYPE, ST, ED, SWEEP, N, NB, IB, LDA, LDVT * .. * .. Array Arguments .. -* REAL A( LDA, * ), V( * ), +* REAL A( LDA, * ), V( * ), * TAU( * ), WORK( * ) -* +* *> \par Purpose: * ============= *> @@ -124,7 +124,7 @@ *> LDVT is INTEGER. *> \endverbatim *> -*> \param[in] WORK +*> \param[out] WORK *> \verbatim *> WORK is REAL array. Workspace of size nb. *> \endverbatim @@ -150,7 +150,7 @@ *> http://doi.acm.org/10.1145/2063384.2063394 *> *> A. Haidar, J. Kurzak, P. Luszczek, 2013. -*> An improved parallel singular value algorithm and its implementation +*> An improved parallel singular value algorithm and its implementation *> for multicore hardware, In Proceedings of 2013 International Conference *> for High Performance Computing, Networking, Storage and Analysis (SC '13). *> Denver, Colorado, USA, 2013. @@ -158,16 +158,16 @@ *> http://doi.acm.org/10.1145/2503210.2503292 *> *> A. Haidar, R. Solca, S. Tomov, T. Schulthess and J. Dongarra. -*> A novel hybrid CPU-GPU generalized eigensolver for electronic structure +*> A novel hybrid CPU-GPU generalized eigensolver for electronic structure *> calculations based on fine-grained memory aware tasks. *> International Journal of High Performance Computing Applications. *> Volume 28 Issue 2, Pages 196-209, May 2014. -*> http://hpc.sagepub.com/content/28/2/196 +*> http://hpc.sagepub.com/content/28/2/196 *> *> \endverbatim *> * ===================================================================== - SUBROUTINE SSB2ST_KERNELS( UPLO, WANTZ, TTYPE, + SUBROUTINE SSB2ST_KERNELS( UPLO, WANTZ, TTYPE, $ ST, ED, SWEEP, N, NB, IB, $ A, LDA, V, TAU, LDVT, WORK) * @@ -184,7 +184,7 @@ INTEGER TTYPE, ST, ED, SWEEP, N, NB, IB, LDA, LDVT * .. * .. Array Arguments .. - REAL A( LDA, * ), V( * ), + REAL A( LDA, * ), V( * ), $ TAU( * ), WORK( * ) * .. * @@ -198,8 +198,8 @@ * .. Local Scalars .. LOGICAL UPPER INTEGER I, J1, J2, LM, LN, VPOS, TAUPOS, - $ DPOS, OFDPOS, AJETER - REAL CTMP + $ DPOS, OFDPOS, AJETER + REAL CTMP * .. * .. External Subroutines .. EXTERNAL SLARFG, SLARFX, SLARFY @@ -212,7 +212,7 @@ * .. * .. * .. Executable Statements .. -* +* AJETER = IB + LDVT UPPER = LSAME( UPLO, 'U' ) @@ -243,10 +243,10 @@ V( VPOS ) = ONE DO 10 I = 1, LM-1 V( VPOS+I ) = ( A( OFDPOS-I, ST+I ) ) - A( OFDPOS-I, ST+I ) = ZERO + A( OFDPOS-I, ST+I ) = ZERO 10 CONTINUE CTMP = ( A( OFDPOS, ST ) ) - CALL SLARFG( LM, CTMP, V( VPOS+1 ), 1, + CALL SLARFG( LM, CTMP, V( VPOS+1 ), 1, $ TAU( TAUPOS ) ) A( OFDPOS, ST ) = CTMP * @@ -284,14 +284,14 @@ * V( VPOS ) = ONE DO 30 I = 1, LM-1 - V( VPOS+I ) = + V( VPOS+I ) = $ ( A( DPOS-NB-I, J1+I ) ) A( DPOS-NB-I, J1+I ) = ZERO 30 CONTINUE CTMP = ( A( DPOS-NB, J1 ) ) CALL SLARFG( LM, CTMP, V( VPOS+1 ), 1, TAU( TAUPOS ) ) A( DPOS-NB, J1 ) = CTMP -* +* CALL SLARFX( 'Right', LN-1, LM, V( VPOS ), $ TAU( TAUPOS ), $ A( DPOS-NB+1, J1 ), LDA-1, WORK) @@ -299,9 +299,9 @@ ENDIF * * Lower case -* +* ELSE -* +* IF( WANTZ ) THEN VPOS = MOD( SWEEP-1, 2 ) * N + ST TAUPOS = MOD( SWEEP-1, 2 ) * N + ST @@ -316,9 +316,9 @@ V( VPOS ) = ONE DO 20 I = 1, LM-1 V( VPOS+I ) = A( OFDPOS+I, ST-1 ) - A( OFDPOS+I, ST-1 ) = ZERO + A( OFDPOS+I, ST-1 ) = ZERO 20 CONTINUE - CALL SLARFG( LM, A( OFDPOS, ST-1 ), V( VPOS+1 ), 1, + CALL SLARFG( LM, A( OFDPOS, ST-1 ), V( VPOS+1 ), 1, $ TAU( TAUPOS ) ) * LM = ED - ST + 1 @@ -345,7 +345,7 @@ LM = J2-J1+1 * IF( LM.GT.0) THEN - CALL SLARFX( 'Right', LM, LN, V( VPOS ), + CALL SLARFX( 'Right', LM, LN, V( VPOS ), $ TAU( TAUPOS ), A( DPOS+NB, ST ), $ LDA-1, WORK) * @@ -362,13 +362,13 @@ V( VPOS+I ) = A( DPOS+NB+I, ST ) A( DPOS+NB+I, ST ) = ZERO 40 CONTINUE - CALL SLARFG( LM, A( DPOS+NB, ST ), V( VPOS+1 ), 1, + CALL SLARFG( LM, A( DPOS+NB, ST ), V( VPOS+1 ), 1, $ TAU( TAUPOS ) ) * - CALL SLARFX( 'Left', LM, LN-1, V( VPOS ), + CALL SLARFX( 'Left', LM, LN-1, V( VPOS ), $ ( TAU( TAUPOS ) ), $ A( DPOS+NB-1, ST+1 ), LDA-1, WORK) - + ENDIF ENDIF ENDIF @@ -377,4 +377,4 @@ * * END OF SSB2ST_KERNELS * - END + END diff --git a/lapack-netlib/SRC/ssbgvx.f b/lapack-netlib/SRC/ssbgvx.f index 3408810bd..22f96729e 100644 --- a/lapack-netlib/SRC/ssbgvx.f +++ b/lapack-netlib/SRC/ssbgvx.f @@ -261,11 +261,11 @@ *> \param[out] INFO *> \verbatim *> INFO is INTEGER -*> = 0 : successful exit -*> < 0 : if INFO = -i, the i-th argument had an illegal value +*> = 0: successful exit +*> < 0: if INFO = -i, the i-th argument had an illegal value *> <= N: if INFO = i, then i eigenvectors failed to converge. *> Their indices are stored in IFAIL. -*> > N : SPBSTF returned an error code; i.e., +*> > N: SPBSTF returned an error code; i.e., *> if INFO = N + i, for 1 <= i <= N, then the leading *> minor of order i of B is not positive definite. *> The factorization of B could not be completed and diff --git a/lapack-netlib/SRC/sstemr.f b/lapack-netlib/SRC/sstemr.f index 228538161..d550f87e0 100644 --- a/lapack-netlib/SRC/sstemr.f +++ b/lapack-netlib/SRC/sstemr.f @@ -233,13 +233,13 @@ *> \param[in,out] TRYRAC *> \verbatim *> TRYRAC is LOGICAL -*> If TRYRAC.EQ..TRUE., indicates that the code should check whether +*> If TRYRAC = .TRUE., indicates that the code should check whether *> the tridiagonal matrix defines its eigenvalues to high relative *> accuracy. If so, the code uses relative-accuracy preserving *> algorithms that might be (a bit) slower depending on the matrix. *> If the matrix does not define its eigenvalues to high relative *> accuracy, the code can uses possibly faster algorithms. -*> If TRYRAC.EQ..FALSE., the code is not required to guarantee +*> If TRYRAC = .FALSE., the code is not required to guarantee *> relatively accurate eigenvalues and can use the fastest possible *> techniques. *> On exit, a .TRUE. TRYRAC will be set to .FALSE. if the matrix diff --git a/lapack-netlib/SRC/ssyconvf.f b/lapack-netlib/SRC/ssyconvf.f index d43b9473f..c6f08428f 100644 --- a/lapack-netlib/SRC/ssyconvf.f +++ b/lapack-netlib/SRC/ssyconvf.f @@ -291,7 +291,7 @@ * * Convert PERMUTATIONS and IPIV * -* Apply permutaions to submatrices of upper part of A +* Apply permutations to submatrices of upper part of A * in factorization order where i decreases from N to 1 * I = N @@ -344,7 +344,7 @@ * * Revert PERMUTATIONS and IPIV * -* Apply permutaions to submatrices of upper part of A +* Apply permutations to submatrices of upper part of A * in reverse factorization order where i increases from 1 to N * I = 1 @@ -435,7 +435,7 @@ * * Convert PERMUTATIONS and IPIV * -* Apply permutaions to submatrices of lower part of A +* Apply permutations to submatrices of lower part of A * in factorization order where k increases from 1 to N * I = 1 @@ -488,7 +488,7 @@ * * Revert PERMUTATIONS and IPIV * -* Apply permutaions to submatrices of lower part of A +* Apply permutations to submatrices of lower part of A * in reverse factorization order where i decreases from N to 1 * I = N diff --git a/lapack-netlib/SRC/ssyconvf_rook.f b/lapack-netlib/SRC/ssyconvf_rook.f index 833b9c632..a7e0d5258 100644 --- a/lapack-netlib/SRC/ssyconvf_rook.f +++ b/lapack-netlib/SRC/ssyconvf_rook.f @@ -282,7 +282,7 @@ * * Convert PERMUTATIONS * -* Apply permutaions to submatrices of upper part of A +* Apply permutations to submatrices of upper part of A * in factorization order where i decreases from N to 1 * I = N @@ -333,7 +333,7 @@ * * Revert PERMUTATIONS * -* Apply permutaions to submatrices of upper part of A +* Apply permutations to submatrices of upper part of A * in reverse factorization order where i increases from 1 to N * I = 1 @@ -423,7 +423,7 @@ * * Convert PERMUTATIONS * -* Apply permutaions to submatrices of lower part of A +* Apply permutations to submatrices of lower part of A * in factorization order where i increases from 1 to N * I = 1 @@ -474,7 +474,7 @@ * * Revert PERMUTATIONS * -* Apply permutaions to submatrices of lower part of A +* Apply permutations to submatrices of lower part of A * in reverse factorization order where i decreases from N to 1 * I = N diff --git a/lapack-netlib/SRC/ssyev_2stage.f b/lapack-netlib/SRC/ssyev_2stage.f index 166766919..5d354c1b3 100644 --- a/lapack-netlib/SRC/ssyev_2stage.f +++ b/lapack-netlib/SRC/ssyev_2stage.f @@ -317,7 +317,7 @@ IF( .NOT.WANTZ ) THEN CALL SSTERF( N, W, WORK( INDE ), INFO ) ELSE -* Not available in this release, and agrument checking should not +* Not available in this release, and argument checking should not * let it getting here RETURN CALL SORGTR( UPLO, N, A, LDA, WORK( INDTAU ), WORK( INDWRK ), diff --git a/lapack-netlib/SRC/ssyevd_2stage.f b/lapack-netlib/SRC/ssyevd_2stage.f index 8ab90b641..625713b85 100644 --- a/lapack-netlib/SRC/ssyevd_2stage.f +++ b/lapack-netlib/SRC/ssyevd_2stage.f @@ -385,7 +385,7 @@ IF( .NOT.WANTZ ) THEN CALL SSTERF( N, W, WORK( INDE ), INFO ) ELSE -* Not available in this release, and agrument checking should not +* Not available in this release, and argument checking should not * let it getting here RETURN CALL SSTEDC( 'I', N, W, WORK( INDE ), WORK( INDWRK ), N, diff --git a/lapack-netlib/SRC/ssyrfsx.f b/lapack-netlib/SRC/ssyrfsx.f index b5dd0b2df..bfb7b6005 100644 --- a/lapack-netlib/SRC/ssyrfsx.f +++ b/lapack-netlib/SRC/ssyrfsx.f @@ -271,7 +271,7 @@ *> information as described below. There currently are up to three *> pieces of information returned for each right-hand side. If *> componentwise accuracy is not requested (PARAMS(3) = 0.0), then -*> ERR_BNDS_COMP is not accessed. If N_ERR_BNDS .LT. 3, then at most +*> ERR_BNDS_COMP is not accessed. If N_ERR_BNDS < 3, then at most *> the first (:,N_ERR_BNDS) entries are returned. *> *> The first index in ERR_BNDS_COMP(i,:) corresponds to the ith @@ -307,14 +307,14 @@ *> \param[in] NPARAMS *> \verbatim *> NPARAMS is INTEGER -*> Specifies the number of parameters set in PARAMS. If .LE. 0, the +*> Specifies the number of parameters set in PARAMS. If <= 0, the *> PARAMS array is never referenced and default values are used. *> \endverbatim *> *> \param[in,out] PARAMS *> \verbatim *> PARAMS is REAL array, dimension NPARAMS -*> Specifies algorithm parameters. If an entry is .LT. 0.0, then +*> Specifies algorithm parameters. If an entry is < 0.0, then *> that entry will be filled with default value used for that *> parameter. Only positions up to NPARAMS are accessed; defaults *> are used for higher-numbered parameters. @@ -322,9 +322,9 @@ *> PARAMS(LA_LINRX_ITREF_I = 1) : Whether to perform iterative *> refinement or not. *> Default: 1.0 -*> = 0.0 : No refinement is performed, and no error bounds are +*> = 0.0: No refinement is performed, and no error bounds are *> computed. -*> = 1.0 : Use the double-precision refinement algorithm, +*> = 1.0: Use the double-precision refinement algorithm, *> possibly with doubled-single computations if the *> compilation environment does not support DOUBLE *> PRECISION. diff --git a/lapack-netlib/SRC/ssysv_aa.f b/lapack-netlib/SRC/ssysv_aa.f index e470f5883..7e58d1e75 100644 --- a/lapack-netlib/SRC/ssysv_aa.f +++ b/lapack-netlib/SRC/ssysv_aa.f @@ -42,7 +42,7 @@ *> matrices. *> *> Aasen's algorithm is used to factor A as -*> A = U * T * U**T, if UPLO = 'U', or +*> A = U**T * T * U, if UPLO = 'U', or *> A = L * T * L**T, if UPLO = 'L', *> where U (or L) is a product of permutation and unit upper (lower) *> triangular matrices, and T is symmetric tridiagonal. The factored @@ -86,7 +86,7 @@ *> *> On exit, if INFO = 0, the tridiagonal matrix T and the *> multipliers used to obtain the factor U or L from the -*> factorization A = U*T*U**T or A = L*T*L**T as computed by +*> factorization A = U**T*T*U or A = L*T*L**T as computed by *> SSYTRF. *> \endverbatim *> @@ -229,7 +229,7 @@ RETURN END IF * -* Compute the factorization A = U*T*U**T or A = L*T*L**T. +* Compute the factorization A = U**T*T*U or A = L*T*L**T. * CALL SSYTRF_AA( UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO ) IF( INFO.EQ.0 ) THEN diff --git a/lapack-netlib/SRC/ssysv_aa_2stage.f b/lapack-netlib/SRC/ssysv_aa_2stage.f index 43d937141..5e2e0e340 100644 --- a/lapack-netlib/SRC/ssysv_aa_2stage.f +++ b/lapack-netlib/SRC/ssysv_aa_2stage.f @@ -44,7 +44,7 @@ *> matrices. *> *> Aasen's 2-stage algorithm is used to factor A as -*> A = U * T * U**T, if UPLO = 'U', or +*> A = U**T * T * U, if UPLO = 'U', or *> A = L * T * L**T, if UPLO = 'L', *> where U (or L) is a product of permutation and unit upper (lower) *> triangular matrices, and T is symmetric and band. The matrix T is @@ -258,7 +258,7 @@ END IF * * -* Compute the factorization A = U*T*U**T or A = L*T*L**T. +* Compute the factorization A = U**T*T*U or A = L*T*L**T. * CALL SSYTRF_AA_2STAGE( UPLO, N, A, LDA, TB, LTB, IPIV, IPIV2, $ WORK, LWORK, INFO ) diff --git a/lapack-netlib/SRC/ssysvxx.f b/lapack-netlib/SRC/ssysvxx.f index 4762748c0..e2be0128b 100644 --- a/lapack-netlib/SRC/ssysvxx.f +++ b/lapack-netlib/SRC/ssysvxx.f @@ -377,7 +377,7 @@ *> information as described below. There currently are up to three *> pieces of information returned for each right-hand side. If *> componentwise accuracy is not requested (PARAMS(3) = 0.0), then -*> ERR_BNDS_COMP is not accessed. If N_ERR_BNDS .LT. 3, then at most +*> ERR_BNDS_COMP is not accessed. If N_ERR_BNDS < 3, then at most *> the first (:,N_ERR_BNDS) entries are returned. *> *> The first index in ERR_BNDS_COMP(i,:) corresponds to the ith @@ -413,14 +413,14 @@ *> \param[in] NPARAMS *> \verbatim *> NPARAMS is INTEGER -*> Specifies the number of parameters set in PARAMS. If .LE. 0, the +*> Specifies the number of parameters set in PARAMS. If <= 0, the *> PARAMS array is never referenced and default values are used. *> \endverbatim *> *> \param[in,out] PARAMS *> \verbatim *> PARAMS is REAL array, dimension NPARAMS -*> Specifies algorithm parameters. If an entry is .LT. 0.0, then +*> Specifies algorithm parameters. If an entry is < 0.0, then *> that entry will be filled with default value used for that *> parameter. Only positions up to NPARAMS are accessed; defaults *> are used for higher-numbered parameters. @@ -428,9 +428,9 @@ *> PARAMS(LA_LINRX_ITREF_I = 1) : Whether to perform iterative *> refinement or not. *> Default: 1.0 -*> = 0.0 : No refinement is performed, and no error bounds are +*> = 0.0: No refinement is performed, and no error bounds are *> computed. -*> = 1.0 : Use the double-precision refinement algorithm, +*> = 1.0: Use the double-precision refinement algorithm, *> possibly with doubled-single computations if the *> compilation environment does not support DOUBLE *> PRECISION. diff --git a/lapack-netlib/SRC/ssytf2_rk.f b/lapack-netlib/SRC/ssytf2_rk.f index bf113d1bd..400e48353 100644 --- a/lapack-netlib/SRC/ssytf2_rk.f +++ b/lapack-netlib/SRC/ssytf2_rk.f @@ -312,7 +312,7 @@ * * Factorize A as U*D*U**T using the upper triangle of A * -* Initilize the first entry of array E, where superdiagonal +* Initialize the first entry of array E, where superdiagonal * elements of D are stored * E( 1 ) = ZERO @@ -623,7 +623,7 @@ * * Factorize A as L*D*L**T using the lower triangle of A * -* Initilize the unused last entry of the subdiagonal array E. +* Initialize the unused last entry of the subdiagonal array E. * E( N ) = ZERO * diff --git a/lapack-netlib/SRC/ssytrd_2stage.f b/lapack-netlib/SRC/ssytrd_2stage.f index 7ddc0224e..d2502f483 100644 --- a/lapack-netlib/SRC/ssytrd_2stage.f +++ b/lapack-netlib/SRC/ssytrd_2stage.f @@ -123,23 +123,22 @@ *> *> \param[out] HOUS2 *> \verbatim -*> HOUS2 is REAL array, dimension LHOUS2, that -*> store the Householder representation of the stage2 +*> HOUS2 is REAL array, dimension (LHOUS2) +*> Stores the Householder representation of the stage2 *> band to tridiagonal. *> \endverbatim *> *> \param[in] LHOUS2 *> \verbatim *> LHOUS2 is INTEGER -*> The dimension of the array HOUS2. LHOUS2 = MAX(1, dimension) -*> If LWORK = -1, or LHOUS2=-1, +*> The dimension of the array HOUS2. +*> If LWORK = -1, or LHOUS2 = -1, *> then a query is assumed; the routine *> only calculates the optimal size of the HOUS2 array, returns *> this value as the first entry of the HOUS2 array, and no error *> message related to LHOUS2 is issued by XERBLA. -*> LHOUS2 = MAX(1, dimension) where -*> dimension = 4*N if VECT='N' -*> not available now if VECT='H' +*> If VECT='N', LHOUS2 = max(1, 4*n); +*> if VECT='V', option not yet available. *> \endverbatim *> *> \param[out] WORK diff --git a/lapack-netlib/SRC/ssytrd_sb2st.F b/lapack-netlib/SRC/ssytrd_sb2st.F index 0df1173e4..9893181cf 100644 --- a/lapack-netlib/SRC/ssytrd_sb2st.F +++ b/lapack-netlib/SRC/ssytrd_sb2st.F @@ -50,9 +50,9 @@ * Arguments: * ========== * -*> \param[in] STAGE +*> \param[in] STAGE1 *> \verbatim -*> STAGE is CHARACTER*1 +*> STAGE1 is CHARACTER*1 *> = 'N': "No": to mention that the stage 1 of the reduction *> from dense to band using the ssytrd_sy2sb routine *> was not called before this routine to reproduce AB. @@ -481,7 +481,7 @@ * * Call the kernel * -#if defined(_OPENMP) && _OPENMP >= 201307 +#if defined(_OPENMP) IF( TTYPE.NE.1 ) THEN !$OMP TASK DEPEND(in:WORK(MYID+SHIFT-1)) !$OMP$ DEPEND(in:WORK(MYID-1)) diff --git a/lapack-netlib/SRC/ssytrd_sy2sb.f b/lapack-netlib/SRC/ssytrd_sy2sb.f index 272876700..98169dc00 100644 --- a/lapack-netlib/SRC/ssytrd_sy2sb.f +++ b/lapack-netlib/SRC/ssytrd_sy2sb.f @@ -363,7 +363,7 @@ * * * Set the workspace of the triangular matrix T to zero once such a -* way everytime T is generated the upper/lower portion will be always zero +* way every time T is generated the upper/lower portion will be always zero * CALL SLASET( "A", LDT, KD, ZERO, ZERO, WORK( TPOS ), LDT ) * diff --git a/lapack-netlib/SRC/ssytrf.f b/lapack-netlib/SRC/ssytrf.f index 2c29475df..ae4550f28 100644 --- a/lapack-netlib/SRC/ssytrf.f +++ b/lapack-netlib/SRC/ssytrf.f @@ -39,7 +39,7 @@ *> the Bunch-Kaufman diagonal pivoting method. The form of the *> factorization is *> -*> A = U*D*U**T or A = L*D*L**T +*> A = U**T*D*U or A = L*D*L**T *> *> where U (or L) is a product of permutation and unit upper (lower) *> triangular matrices, and D is symmetric and block diagonal with @@ -144,7 +144,7 @@ *> *> \verbatim *> -*> If UPLO = 'U', then A = U*D*U**T, where +*> If UPLO = 'U', then A = U**T*D*U, where *> U = P(n)*U(n)* ... *P(k)U(k)* ..., *> i.e., U is a product of terms P(k)*U(k), where k decreases from n to *> 1 in steps of 1 or 2, and D is a block diagonal matrix with 1-by-1 @@ -262,7 +262,7 @@ * IF( UPPER ) THEN * -* Factorize A as U*D*U**T using the upper triangle of A +* Factorize A as U**T*D*U using the upper triangle of A * * K is the main loop index, decreasing from N to 1 in steps of * KB, where KB is the number of columns factorized by SLASYF; diff --git a/lapack-netlib/SRC/ssytrf_aa.f b/lapack-netlib/SRC/ssytrf_aa.f index 4aaa978ad..7f428561c 100644 --- a/lapack-netlib/SRC/ssytrf_aa.f +++ b/lapack-netlib/SRC/ssytrf_aa.f @@ -37,7 +37,7 @@ *> SSYTRF_AA computes the factorization of a real symmetric matrix A *> using the Aasen's algorithm. The form of the factorization is *> -*> A = U*T*U**T or A = L*T*L**T +*> A = U**T*T*U or A = L*T*L**T *> *> where U (or L) is a product of permutation and unit upper (lower) *> triangular matrices, and T is a symmetric tridiagonal matrix. @@ -223,7 +223,7 @@ IF( UPPER ) THEN * * ..................................................... -* Factorize A as L*D*L**T using the upper triangle of A +* Factorize A as U**T*D*U using the upper triangle of A * ..................................................... * * Copy first row A(1, 1:N) into H(1:n) (stored in WORK(1:N)) @@ -256,7 +256,7 @@ $ A( MAX(1, J), J+1 ), LDA, $ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ) ) * -* Ajust IPIV and apply it back (J-th step picks (J+1)-th pivot) +* Adjust IPIV and apply it back (J-th step picks (J+1)-th pivot) * DO J2 = J+2, MIN(N, J+JB+1) IPIV( J2 ) = IPIV( J2 ) + J @@ -375,7 +375,7 @@ $ A( J+1, MAX(1, J) ), LDA, $ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ) ) * -* Ajust IPIV and apply it back (J-th step picks (J+1)-th pivot) +* Adjust IPIV and apply it back (J-th step picks (J+1)-th pivot) * DO J2 = J+2, MIN(N, J+JB+1) IPIV( J2 ) = IPIV( J2 ) + J diff --git a/lapack-netlib/SRC/ssytrf_aa_2stage.f b/lapack-netlib/SRC/ssytrf_aa_2stage.f index 0e0f6edb7..03690815b 100644 --- a/lapack-netlib/SRC/ssytrf_aa_2stage.f +++ b/lapack-netlib/SRC/ssytrf_aa_2stage.f @@ -38,7 +38,7 @@ *> SSYTRF_AA_2STAGE computes the factorization of a real symmetric matrix A *> using the Aasen's algorithm. The form of the factorization is *> -*> A = U*T*U**T or A = L*T*L**T +*> A = U**T*T*U or A = L*T*L**T *> *> where U (or L) is a product of permutation and unit upper (lower) *> triangular matrices, and T is a symmetric band matrix with the @@ -275,7 +275,7 @@ IF( UPPER ) THEN * * ..................................................... -* Factorize A as L*D*L**T using the upper triangle of A +* Factorize A as U**T*D*U using the upper triangle of A * ..................................................... * DO J = 0, NT-1 @@ -442,12 +442,14 @@ c END IF * > Apply pivots to previous columns of L CALL SSWAP( K-1, A( (J+1)*NB+1, I1 ), 1, $ A( (J+1)*NB+1, I2 ), 1 ) -* > Swap A(I1+1:M, I1) with A(I2, I1+1:M) - CALL SSWAP( I2-I1-1, A( I1, I1+1 ), LDA, - $ A( I1+1, I2 ), 1 ) +* > Swap A(I1+1:M, I1) with A(I2, I1+1:M) + IF( I2.GT.(I1+1) ) + $ CALL SSWAP( I2-I1-1, A( I1, I1+1 ), LDA, + $ A( I1+1, I2 ), 1 ) * > Swap A(I2+1:M, I1) with A(I2+1:M, I2) - CALL SSWAP( N-I2, A( I1, I2+1 ), LDA, - $ A( I2, I2+1 ), LDA ) + IF( I2.LT.N ) + $ CALL SSWAP( N-I2, A( I1, I2+1 ), LDA, + $ A( I2, I2+1 ), LDA ) * > Swap A(I1, I1) with A(I2, I2) PIV = A( I1, I1 ) A( I1, I1 ) = A( I2, I2 ) @@ -616,11 +618,13 @@ c END IF CALL SSWAP( K-1, A( I1, (J+1)*NB+1 ), LDA, $ A( I2, (J+1)*NB+1 ), LDA ) * > Swap A(I1+1:M, I1) with A(I2, I1+1:M) - CALL SSWAP( I2-I1-1, A( I1+1, I1 ), 1, - $ A( I2, I1+1 ), LDA ) + IF( I2.GT.(I1+1) ) + $ CALL SSWAP( I2-I1-1, A( I1+1, I1 ), 1, + $ A( I2, I1+1 ), LDA ) * > Swap A(I2+1:M, I1) with A(I2+1:M, I2) - CALL SSWAP( N-I2, A( I2+1, I1 ), 1, - $ A( I2+1, I2 ), 1 ) + IF( I2.LT.N ) + $ CALL SSWAP( N-I2, A( I2+1, I1 ), 1, + $ A( I2+1, I2 ), 1 ) * > Swap A(I1, I1) with A(I2, I2) PIV = A( I1, I1 ) A( I1, I1 ) = A( I2, I2 ) diff --git a/lapack-netlib/SRC/ssytri2.f b/lapack-netlib/SRC/ssytri2.f index 4b9ea4e7b..897116c23 100644 --- a/lapack-netlib/SRC/ssytri2.f +++ b/lapack-netlib/SRC/ssytri2.f @@ -62,7 +62,7 @@ *> \param[in,out] A *> \verbatim *> A is REAL array, dimension (LDA,N) -*> On entry, the NB diagonal matrix D and the multipliers +*> On entry, the block diagonal matrix D and the multipliers *> used to obtain the factor U or L as computed by SSYTRF. *> *> On exit, if INFO = 0, the (symmetric) inverse of the original @@ -82,7 +82,7 @@ *> \param[in] IPIV *> \verbatim *> IPIV is INTEGER array, dimension (N) -*> Details of the interchanges and the NB structure of D +*> Details of the interchanges and the block structure of D *> as determined by SSYTRF. *> \endverbatim *> diff --git a/lapack-netlib/SRC/ssytrs_aa.f b/lapack-netlib/SRC/ssytrs_aa.f index b05c9f7e6..ed4377ae7 100644 --- a/lapack-netlib/SRC/ssytrs_aa.f +++ b/lapack-netlib/SRC/ssytrs_aa.f @@ -37,7 +37,7 @@ *> \verbatim *> *> SSYTRS_AA solves a system of linear equations A*X = B with a real -*> symmetric matrix A using the factorization A = U*T*U**T or +*> symmetric matrix A using the factorization A = U**T*T*U or *> A = L*T*L**T computed by SSYTRF_AA. *> \endverbatim * @@ -49,7 +49,7 @@ *> UPLO is CHARACTER*1 *> Specifies whether the details of the factorization are stored *> as an upper or lower triangular matrix. -*> = 'U': Upper triangular, form is A = U*T*U**T; +*> = 'U': Upper triangular, form is A = U**T*T*U; *> = 'L': Lower triangular, form is A = L*T*L**T. *> \endverbatim *> @@ -97,14 +97,16 @@ *> The leading dimension of the array B. LDB >= max(1,N). *> \endverbatim *> -*> \param[in] WORK +*> \param[out] WORK *> \verbatim -*> WORK is DOUBLE array, dimension (MAX(1,LWORK)) +*> WORK is REAL array, dimension (MAX(1,LWORK)) *> \endverbatim *> *> \param[in] LWORK *> \verbatim -*> LWORK is INTEGER, LWORK >= MAX(1,3*N-2). +*> LWORK is INTEGER +*> The dimension of the array WORK. LWORK >= max(1,3*N-2). +*> \endverbatim *> *> \param[out] INFO *> \verbatim @@ -198,24 +200,31 @@ * IF( UPPER ) THEN * -* Solve A*X = B, where A = U*T*U**T. +* Solve A*X = B, where A = U**T*T*U. * -* Pivot, P**T * B +* 1) Forward substitution with U**T * - K = 1 - DO WHILE ( K.LE.N ) - KP = IPIV( K ) - IF( KP.NE.K ) - $ CALL SSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) - K = K + 1 - END DO + IF( N.GT.1 ) THEN * -* Compute (U \P**T * B) -> B [ (U \P**T * B) ] +* Pivot, P**T * B -> B * - CALL STRSM('L', 'U', 'T', 'U', N-1, NRHS, ONE, A( 1, 2 ), LDA, - $ B( 2, 1 ), LDB) + K = 1 + DO WHILE ( K.LE.N ) + KP = IPIV( K ) + IF( KP.NE.K ) + $ CALL SSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) + K = K + 1 + END DO * -* Compute T \ B -> B [ T \ (U \P**T * B) ] +* Compute U**T \ B -> B [ (U**T \P**T * B) ] +* + CALL STRSM( 'L', 'U', 'T', 'U', N-1, NRHS, ONE, A( 1, 2 ), + $ LDA, B( 2, 1 ), LDB) + END IF +* +* 2) Solve with triangular matrix T +* +* Compute T \ B -> B [ T \ (U**T \P**T * B) ] * CALL SLACPY( 'F', 1, N, A(1, 1), LDA+1, WORK(N), 1) IF( N.GT.1 ) THEN @@ -224,41 +233,53 @@ END IF CALL SGTSV(N, NRHS, WORK(1), WORK(N), WORK(2*N), B, LDB, $ INFO) +* +* 3) Backward substitution with U +* + IF( N.GT.1 ) THEN * * -* Compute (U**T \ B) -> B [ U**T \ (T \ (U \P**T * B) ) ] +* Compute U \ B -> B [ U \ (T \ (U**T \P**T * B) ) ] * - CALL STRSM( 'L', 'U', 'N', 'U', N-1, NRHS, ONE, A( 1, 2 ), LDA, - $ B(2, 1), LDB) + CALL STRSM( 'L', 'U', 'N', 'U', N-1, NRHS, ONE, A( 1, 2 ), + $ LDA, B(2, 1), LDB) * -* Pivot, P * B [ P * (U**T \ (T \ (U \P**T * B) )) ] +* Pivot, P * B -> B [ P * (U \ (T \ (U**T \P**T * B) )) ] * - K = N - DO WHILE ( K.GE.1 ) - KP = IPIV( K ) - IF( KP.NE.K ) - $ CALL SSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) - K = K - 1 - END DO + K = N + DO WHILE ( K.GE.1 ) + KP = IPIV( K ) + IF( KP.NE.K ) + $ CALL SSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) + K = K - 1 + END DO + END IF * ELSE * * Solve A*X = B, where A = L*T*L**T. * -* Pivot, P**T * B +* 1) Forward substitution with L * - K = 1 - DO WHILE ( K.LE.N ) - KP = IPIV( K ) - IF( KP.NE.K ) - $ CALL SSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) - K = K + 1 - END DO + IF( N.GT.1 ) THEN * -* Compute (L \P**T * B) -> B [ (L \P**T * B) ] +* Pivot, P**T * B -> B * - CALL STRSM( 'L', 'L', 'N', 'U', N-1, NRHS, ONE, A( 2, 1), LDA, - $ B(2, 1), LDB) + K = 1 + DO WHILE ( K.LE.N ) + KP = IPIV( K ) + IF( KP.NE.K ) + $ CALL SSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) + K = K + 1 + END DO +* +* Compute L \ B -> B [ (L \P**T * B) ] +* + CALL STRSM( 'L', 'L', 'N', 'U', N-1, NRHS, ONE, A( 2, 1), + $ LDA, B(2, 1), LDB) + END IF +* +* 2) Solve with triangular matrix T * * Compute T \ B -> B [ T \ (L \P**T * B) ] * @@ -270,20 +291,25 @@ CALL SGTSV(N, NRHS, WORK(1), WORK(N), WORK(2*N), B, LDB, $ INFO) * -* Compute (L**T \ B) -> B [ L**T \ (T \ (L \P**T * B) ) ] +* 3) Backward substitution with L**T * - CALL STRSM( 'L', 'L', 'T', 'U', N-1, NRHS, ONE, A( 2, 1 ), LDA, - $ B( 2, 1 ), LDB) + IF( N.GT.1 ) THEN * -* Pivot, P * B [ P * (L**T \ (T \ (L \P**T * B) )) ] +* Compute L**T \ B -> B [ L**T \ (T \ (L \P**T * B) ) ] * - K = N - DO WHILE ( K.GE.1 ) - KP = IPIV( K ) - IF( KP.NE.K ) - $ CALL SSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) - K = K - 1 - END DO + CALL STRSM( 'L', 'L', 'T', 'U', N-1, NRHS, ONE, A( 2, 1 ), + $ LDA, B( 2, 1 ), LDB) +* +* Pivot, P * B -> B [ P * (L**T \ (T \ (L \P**T * B) )) ] +* + K = N + DO WHILE ( K.GE.1 ) + KP = IPIV( K ) + IF( KP.NE.K ) + $ CALL SSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) + K = K - 1 + END DO + END IF * END IF * diff --git a/lapack-netlib/SRC/ssytrs_aa_2stage.f b/lapack-netlib/SRC/ssytrs_aa_2stage.f index d271b9481..cf2da529d 100644 --- a/lapack-netlib/SRC/ssytrs_aa_2stage.f +++ b/lapack-netlib/SRC/ssytrs_aa_2stage.f @@ -36,7 +36,7 @@ *> \verbatim *> *> SSYTRS_AA_2STAGE solves a system of linear equations A*X = B with a real -*> symmetric matrix A using the factorization A = U*T*U**T or +*> symmetric matrix A using the factorization A = U**T*T*U or *> A = L*T*L**T computed by SSYTRF_AA_2STAGE. *> \endverbatim * @@ -48,7 +48,7 @@ *> UPLO is CHARACTER*1 *> Specifies whether the details of the factorization are stored *> as an upper or lower triangular matrix. -*> = 'U': Upper triangular, form is A = U*T*U**T; +*> = 'U': Upper triangular, form is A = U**T*T*U; *> = 'L': Lower triangular, form is A = L*T*L**T. *> \endverbatim *> @@ -208,15 +208,15 @@ * IF( UPPER ) THEN * -* Solve A*X = B, where A = U*T*U**T. +* Solve A*X = B, where A = U**T*T*U. * IF( N.GT.NB ) THEN * -* Pivot, P**T * B +* Pivot, P**T * B -> B * CALL SLASWP( NRHS, B, LDB, NB+1, N, IPIV, 1 ) * -* Compute (U**T \P**T * B) -> B [ (U**T \P**T * B) ] +* Compute (U**T \ B) -> B [ (U**T \P**T * B) ] * CALL STRSM( 'L', 'U', 'T', 'U', N-NB, NRHS, ONE, A(1, NB+1), $ LDA, B(NB+1, 1), LDB) @@ -234,7 +234,7 @@ CALL STRSM( 'L', 'U', 'N', 'U', N-NB, NRHS, ONE, A(1, NB+1), $ LDA, B(NB+1, 1), LDB) * -* Pivot, P * B [ P * (U \ (T \ (U**T \P**T * B) )) ] +* Pivot, P * B -> B [ P * (U \ (T \ (U**T \P**T * B) )) ] * CALL SLASWP( NRHS, B, LDB, NB+1, N, IPIV, -1 ) * @@ -246,11 +246,11 @@ * IF( N.GT.NB ) THEN * -* Pivot, P**T * B +* Pivot, P**T * B -> B * CALL SLASWP( NRHS, B, LDB, NB+1, N, IPIV, 1 ) * -* Compute (L \P**T * B) -> B [ (L \P**T * B) ] +* Compute (L \ B) -> B [ (L \P**T * B) ] * CALL STRSM( 'L', 'L', 'N', 'U', N-NB, NRHS, ONE, A(NB+1, 1), $ LDA, B(NB+1, 1), LDB) @@ -268,7 +268,7 @@ CALL STRSM( 'L', 'L', 'T', 'U', N-NB, NRHS, ONE, A(NB+1, 1), $ LDA, B(NB+1, 1), LDB) * -* Pivot, P * B [ P * (L**T \ (T \ (L \P**T * B) )) ] +* Pivot, P * B -> B [ P * (L**T \ (T \ (L \P**T * B) )) ] * CALL SLASWP( NRHS, B, LDB, NB+1, N, IPIV, -1 ) * diff --git a/lapack-netlib/SRC/stgsy2.f b/lapack-netlib/SRC/stgsy2.f index ca9946a7e..2814889fc 100644 --- a/lapack-netlib/SRC/stgsy2.f +++ b/lapack-netlib/SRC/stgsy2.f @@ -71,7 +71,7 @@ *> R * B**T + L * E**T = scale * -F *> *> This case is used to compute an estimate of Dif[(A, D), (B, E)] = -*> sigma_min(Z) using reverse communicaton with SLACON. +*> sigma_min(Z) using reverse communication with SLACON. *> *> STGSY2 also (IJOB >= 1) contributes to the computation in STGSYL *> of an upper bound on the separation between to matrix pairs. Then @@ -85,7 +85,7 @@ *> \param[in] TRANS *> \verbatim *> TRANS is CHARACTER*1 -*> = 'N', solve the generalized Sylvester equation (1). +*> = 'N': solve the generalized Sylvester equation (1). *> = 'T': solve the 'transposed' system (3). *> \endverbatim *> diff --git a/lapack-netlib/SRC/stgsyl.f b/lapack-netlib/SRC/stgsyl.f index cd597f37d..ff634b1de 100644 --- a/lapack-netlib/SRC/stgsyl.f +++ b/lapack-netlib/SRC/stgsyl.f @@ -88,20 +88,20 @@ *> \param[in] TRANS *> \verbatim *> TRANS is CHARACTER*1 -*> = 'N', solve the generalized Sylvester equation (1). -*> = 'T', solve the 'transposed' system (3). +*> = 'N': solve the generalized Sylvester equation (1). +*> = 'T': solve the 'transposed' system (3). *> \endverbatim *> *> \param[in] IJOB *> \verbatim *> IJOB is INTEGER *> Specifies what kind of functionality to be performed. -*> =0: solve (1) only. -*> =1: The functionality of 0 and 3. -*> =2: The functionality of 0 and 4. -*> =3: Only an estimate of Dif[(A,D), (B,E)] is computed. +*> = 0: solve (1) only. +*> = 1: The functionality of 0 and 3. +*> = 2: The functionality of 0 and 4. +*> = 3: Only an estimate of Dif[(A,D), (B,E)] is computed. *> (look ahead strategy IJOB = 1 is used). -*> =4: Only an estimate of Dif[(A,D), (B,E)] is computed. +*> = 4: Only an estimate of Dif[(A,D), (B,E)] is computed. *> ( SGECON on sub-systems is used ). *> Not referenced if TRANS = 'T'. *> \endverbatim diff --git a/lapack-netlib/SRC/stpmlqt.f b/lapack-netlib/SRC/stpmlqt.f index 565dadd0c..8fc7823c2 100644 --- a/lapack-netlib/SRC/stpmlqt.f +++ b/lapack-netlib/SRC/stpmlqt.f @@ -94,7 +94,7 @@ *> *> \param[in] V *> \verbatim -*> V is REAL array, dimension (LDA,K) +*> V is REAL array, dimension (LDV,K) *> The i-th row must contain the vector which defines the *> elementary reflector H(i), for i = 1,2,...,k, as returned by *> DTPLQT in B. See Further Details. diff --git a/lapack-netlib/SRC/stpmqrt.f b/lapack-netlib/SRC/stpmqrt.f index b1813b7dd..6a5cbb981 100644 --- a/lapack-netlib/SRC/stpmqrt.f +++ b/lapack-netlib/SRC/stpmqrt.f @@ -94,7 +94,7 @@ *> *> \param[in] V *> \verbatim -*> V is REAL array, dimension (LDA,K) +*> V is REAL array, dimension (LDV,K) *> The i-th column must contain the vector which defines the *> elementary reflector H(i), for i = 1,2,...,k, as returned by *> CTPQRT in B. See Further Details. diff --git a/lapack-netlib/SRC/stprfb.f b/lapack-netlib/SRC/stprfb.f index 66e67252f..fcd164183 100644 --- a/lapack-netlib/SRC/stprfb.f +++ b/lapack-netlib/SRC/stprfb.f @@ -152,8 +152,8 @@ *> \verbatim *> LDA is INTEGER *> The leading dimension of the array A. -*> If SIDE = 'L', LDC >= max(1,K); -*> If SIDE = 'R', LDC >= max(1,M). +*> If SIDE = 'L', LDA >= max(1,K); +*> If SIDE = 'R', LDA >= max(1,M). *> \endverbatim *> *> \param[in,out] B