Skip to content

Commit

Permalink
Handle corner cases of LWORK (Reference-LAPACK PR 942)
Browse files Browse the repository at this point in the history
  • Loading branch information
martin-frbg authored Dec 23, 2023
1 parent 5c11b2f commit 0814491
Show file tree
Hide file tree
Showing 47 changed files with 783 additions and 533 deletions.
32 changes: 20 additions & 12 deletions lapack-netlib/SRC/cgebrd.f
Original file line number Diff line number Diff line change
Expand Up @@ -123,7 +123,8 @@
*> \param[in] LWORK
*> \verbatim
*> LWORK is INTEGER
*> The length of the array WORK. LWORK >= max(1,M,N).
*> The length of the array WORK.
*> LWORK >= 1, if MIN(M,N) = 0, and LWORK >= MAX(M,N), otherwise.
*> For optimum performance LWORK >= (M+N)*NB, where NB
*> is the optimal blocksize.
*>
Expand All @@ -148,7 +149,7 @@
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*
*> \ingroup complexGEcomputational
*> \ingroup gebrd
*
*> \par Further Details:
* =====================
Expand Down Expand Up @@ -225,8 +226,8 @@ SUBROUTINE CGEBRD( M, N, A, LDA, D, E, TAUQ, TAUP, WORK, LWORK,
* ..
* .. Local Scalars ..
LOGICAL LQUERY
INTEGER I, IINFO, J, LDWRKX, LDWRKY, LWKOPT, MINMN, NB,
$ NBMIN, NX, WS
INTEGER I, IINFO, J, LDWRKX, LDWRKY, LWKMIN, LWKOPT,
$ MINMN, NB, NBMIN, NX, WS
* ..
* .. External Subroutines ..
EXTERNAL CGEBD2, CGEMM, CLABRD, XERBLA
Expand All @@ -236,24 +237,32 @@ SUBROUTINE CGEBRD( M, N, A, LDA, D, E, TAUQ, TAUP, WORK, LWORK,
* ..
* .. External Functions ..
INTEGER ILAENV
EXTERNAL ILAENV
REAL SROUNDUP_LWORK
EXTERNAL ILAENV, SROUNDUP_LWORK
* ..
* .. Executable Statements ..
*
* Test the input parameters
*
INFO = 0
NB = MAX( 1, ILAENV( 1, 'CGEBRD', ' ', M, N, -1, -1 ) )
LWKOPT = ( M+N )*NB
WORK( 1 ) = REAL( LWKOPT )
MINMN = MIN( M, N )
IF( MINMN.EQ.0 ) THEN
LWKMIN = 1
LWKOPT = 1
ELSE
LWKMIN = MAX( M, N )
NB = MAX( 1, ILAENV( 1, 'CGEBRD', ' ', M, N, -1, -1 ) )
LWKOPT = ( M+N )*NB
END IF
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
LQUERY = ( LWORK.EQ.-1 )
IF( M.LT.0 ) THEN
INFO = -1
ELSE IF( N.LT.0 ) THEN
INFO = -2
ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
INFO = -4
ELSE IF( LWORK.LT.MAX( 1, M, N ) .AND. .NOT.LQUERY ) THEN
ELSE IF( LWORK.LT.LWKMIN .AND. .NOT.LQUERY ) THEN
INFO = -10
END IF
IF( INFO.LT.0 ) THEN
Expand All @@ -265,7 +274,6 @@ SUBROUTINE CGEBRD( M, N, A, LDA, D, E, TAUQ, TAUP, WORK, LWORK,
*
* Quick return if possible
*
MINMN = MIN( M, N )
IF( MINMN.EQ.0 ) THEN
WORK( 1 ) = 1
RETURN
Expand All @@ -284,7 +292,7 @@ SUBROUTINE CGEBRD( M, N, A, LDA, D, E, TAUQ, TAUP, WORK, LWORK,
* Determine when to switch from blocked to unblocked code.
*
IF( NX.LT.MINMN ) THEN
WS = ( M+N )*NB
WS = LWKOPT
IF( LWORK.LT.WS ) THEN
*
* Not enough work space for the optimal NB, consider using
Expand Down Expand Up @@ -343,7 +351,7 @@ SUBROUTINE CGEBRD( M, N, A, LDA, D, E, TAUQ, TAUP, WORK, LWORK,
*
CALL CGEBD2( M-I+1, N-I+1, A( I, I ), LDA, D( I ), E( I ),
$ TAUQ( I ), TAUP( I ), WORK, IINFO )
WORK( 1 ) = WS
WORK( 1 ) = SROUNDUP_LWORK( WS )
RETURN
*
* End of CGEBRD
Expand Down
20 changes: 13 additions & 7 deletions lapack-netlib/SRC/cgehrd.f
Original file line number Diff line number Diff line change
Expand Up @@ -89,7 +89,7 @@
*>
*> \param[out] WORK
*> \verbatim
*> WORK is COMPLEX array, dimension (LWORK)
*> WORK is COMPLEX array, dimension (MAX(1,LWORK))
*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
*> \endverbatim
*>
Expand Down Expand Up @@ -222,13 +222,19 @@ SUBROUTINE CGEHRD( N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO )
INFO = -8
END IF
*
NH = IHI - ILO + 1
IF( INFO.EQ.0 ) THEN
*
* Compute the workspace requirements
*
NB = MIN( NBMAX, ILAENV( 1, 'CGEHRD', ' ', N, ILO, IHI, -1 ) )
LWKOPT = N*NB + TSIZE
WORK( 1 ) = SROUNDUP_LWORK(LWKOPT)
IF( NH.LE.1 ) THEN
LWKOPT = 1
ELSE
NB = MIN( NBMAX, ILAENV( 1, 'DGEHRD', ' ', N, ILO, IHI,
$ -1 ) )
LWKOPT = N*NB + TSIZE
END IF
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
END IF
*
IF( INFO.NE.0 ) THEN
Expand All @@ -249,7 +255,6 @@ SUBROUTINE CGEHRD( N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO )
*
* Quick return if possible
*
NH = IHI - ILO + 1
IF( NH.LE.1 ) THEN
WORK( 1 ) = 1
RETURN
Expand All @@ -269,7 +274,7 @@ SUBROUTINE CGEHRD( N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO )
*
* Determine if workspace is large enough for blocked code
*
IF( LWORK.LT.N*NB+TSIZE ) THEN
IF( LWORK.LT.LWKOPT ) THEN
*
* Not enough workspace to use optimal NB: determine the
* minimum value of NB, and reduce NB or force use of
Expand Down Expand Up @@ -345,7 +350,8 @@ SUBROUTINE CGEHRD( N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO )
* Use unblocked code to reduce the rest of the matrix
*
CALL CGEHD2( N, I, IHI, A, LDA, TAU, WORK, IINFO )
WORK( 1 ) = SROUNDUP_LWORK(LWKOPT)
*
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
*
RETURN
*
Expand Down
8 changes: 4 additions & 4 deletions lapack-netlib/SRC/cgelq.f
Original file line number Diff line number Diff line change
Expand Up @@ -98,7 +98,7 @@
*> \param[in] LWORK
*> \verbatim
*> LWORK is INTEGER
*> The dimension of the array WORK.
*> The dimension of the array WORK. LWORK >= 1.
*> If LWORK = -1 or -2, then a workspace query is assumed. The routine
*> only calculates the sizes of the T and WORK arrays, returns these
*> values as the first entries of the T and WORK arrays, and no error
Expand Down Expand Up @@ -295,9 +295,9 @@ SUBROUTINE CGELQ( M, N, A, LDA, T, TSIZE, WORK, LWORK,
T( 2 ) = MB
T( 3 ) = NB
IF( MINW ) THEN
WORK( 1 ) = SROUNDUP_LWORK(LWMIN)
WORK( 1 ) = SROUNDUP_LWORK( LWMIN )
ELSE
WORK( 1 ) = SROUNDUP_LWORK(LWREQ)
WORK( 1 ) = SROUNDUP_LWORK( LWREQ )
END IF
END IF
IF( INFO.NE.0 ) THEN
Expand All @@ -322,7 +322,7 @@ SUBROUTINE CGELQ( M, N, A, LDA, T, TSIZE, WORK, LWORK,
$ LWORK, INFO )
END IF
*
WORK( 1 ) = SROUNDUP_LWORK(LWREQ)
WORK( 1 ) = SROUNDUP_LWORK( LWREQ )
*
RETURN
*
Expand Down
20 changes: 13 additions & 7 deletions lapack-netlib/SRC/cgelqf.f
Original file line number Diff line number Diff line change
Expand Up @@ -93,7 +93,8 @@
*> \param[in] LWORK
*> \verbatim
*> LWORK is INTEGER
*> The dimension of the array WORK. LWORK >= max(1,M).
*> The dimension of the array WORK.
*> LWORK >= 1, if MIN(M,N) = 0, and LWORK >= M, otherwise.
*> For optimum performance LWORK >= M*NB, where NB is the
*> optimal blocksize.
*>
Expand Down Expand Up @@ -175,29 +176,34 @@ SUBROUTINE CGELQF( M, N, A, LDA, TAU, WORK, LWORK, INFO )
* Test the input arguments
*
INFO = 0
K = MIN( M, N )
NB = ILAENV( 1, 'CGELQF', ' ', M, N, -1, -1 )
LWKOPT = M*NB
WORK( 1 ) = SROUNDUP_LWORK(LWKOPT)
LQUERY = ( LWORK.EQ.-1 )
IF( M.LT.0 ) THEN
INFO = -1
ELSE IF( N.LT.0 ) THEN
INFO = -2
ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
INFO = -4
ELSE IF( LWORK.LT.MAX( 1, M ) .AND. .NOT.LQUERY ) THEN
INFO = -7
ELSE IF( .NOT.LQUERY ) THEN
IF( LWORK.LE.0 .OR. ( N.GT.0 .AND. LWORK.LT.MAX( 1, M ) ) )
$ INFO = -7
END IF
IF( INFO.NE.0 ) THEN
CALL XERBLA( 'CGELQF', -INFO )
RETURN
ELSE IF( LQUERY ) THEN
IF( K.EQ.0 ) THEN
LWKOPT = 1
ELSE
LWKOPT = M*NB
END IF
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
RETURN
END IF
*
* Quick return if possible
*
K = MIN( M, N )
IF( K.EQ.0 ) THEN
WORK( 1 ) = 1
RETURN
Expand Down Expand Up @@ -267,7 +273,7 @@ SUBROUTINE CGELQF( M, N, A, LDA, TAU, WORK, LWORK, INFO )
$ CALL CGELQ2( M-I+1, N-I+1, A( I, I ), LDA, TAU( I ), WORK,
$ IINFO )
*
WORK( 1 ) = SROUNDUP_LWORK(IWS)
WORK( 1 ) = SROUNDUP_LWORK( IWS )
RETURN
*
* End of CGELQF
Expand Down
35 changes: 23 additions & 12 deletions lapack-netlib/SRC/cgemlq.f
Original file line number Diff line number Diff line change
Expand Up @@ -110,16 +110,17 @@
*>
*> \param[out] WORK
*> \verbatim
*> (workspace) COMPLEX array, dimension (MAX(1,LWORK))
*> (workspace) COMPLEX array, dimension (MAX(1,LWORK))
*> On exit, if INFO = 0, WORK(1) returns the minimal LWORK.
*> \endverbatim
*>
*> \param[in] LWORK
*> \verbatim
*> LWORK is INTEGER
*> The dimension of the array WORK.
*> The dimension of the array WORK. LWORK >= 1.
*> If LWORK = -1, then a workspace query is assumed. The routine
*> only calculates the size of the WORK array, returns this
*> value as WORK(1), and no error message related to WORK
*> value as WORK(1), and no error message related to WORK
*> is issued by XERBLA.
*> \endverbatim
*>
Expand All @@ -143,7 +144,7 @@
*>
*> \verbatim
*>
*> These details are particular for this LAPACK implementation. Users should not
*> These details are particular for this LAPACK implementation. Users should not
*> take them for granted. These details may change in the future, and are not likely
*> true for another LAPACK implementation. These details are relevant if one wants
*> to try to understand the code. They are not part of the interface.
Expand All @@ -159,11 +160,13 @@
*> block sizes MB and NB returned by ILAENV, CGELQ will use either
*> CLASWLQ (if the matrix is wide-and-short) or CGELQT to compute
*> the LQ factorization.
*> This version of CGEMLQ will use either CLAMSWLQ or CGEMLQT to
*> This version of CGEMLQ will use either CLAMSWLQ or CGEMLQT to
*> multiply matrix Q by another matrix.
*> Further Details in CLAMSWLQ or CGEMLQT.
*> \endverbatim
*>
*> \ingroup gemlq
*>
* =====================================================================
SUBROUTINE CGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
$ C, LDC, WORK, LWORK, INFO )
Expand All @@ -185,11 +188,12 @@ SUBROUTINE CGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
* ..
* .. Local Scalars ..
LOGICAL LEFT, RIGHT, TRAN, NOTRAN, LQUERY
INTEGER MB, NB, LW, NBLCKS, MN
INTEGER MB, NB, LW, NBLCKS, MN, MINMNK, LWMIN
* ..
* .. External Functions ..
LOGICAL LSAME
EXTERNAL LSAME
REAL SROUNDUP_LWORK
EXTERNAL LSAME, SROUNDUP_LWORK
* ..
* .. External Subroutines ..
EXTERNAL CLAMSWLQ, CGEMLQT, XERBLA
Expand All @@ -201,7 +205,7 @@ SUBROUTINE CGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
*
* Test the input arguments
*
LQUERY = LWORK.EQ.-1
LQUERY = ( LWORK.EQ.-1 )
NOTRAN = LSAME( TRANS, 'N' )
TRAN = LSAME( TRANS, 'C' )
LEFT = LSAME( SIDE, 'L' )
Expand All @@ -216,6 +220,13 @@ SUBROUTINE CGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
LW = M * MB
MN = N
END IF
*
MINMNK = MIN( M, N, K )
IF( MINMNK.EQ.0 ) THEN
LWMIN = 1
ELSE
LWMIN = MAX( 1, LW )
END IF
*
IF( ( NB.GT.K ) .AND. ( MN.GT.K ) ) THEN
IF( MOD( MN - K, NB - K ) .EQ. 0 ) THEN
Expand Down Expand Up @@ -244,12 +255,12 @@ SUBROUTINE CGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
INFO = -9
ELSE IF( LDC.LT.MAX( 1, M ) ) THEN
INFO = -11
ELSE IF( ( LWORK.LT.MAX( 1, LW ) ) .AND. ( .NOT.LQUERY ) ) THEN
ELSE IF( ( LWORK.LT.LWMIN ) .AND. ( .NOT.LQUERY ) ) THEN
INFO = -13
END IF
*
IF( INFO.EQ.0 ) THEN
WORK( 1 ) = REAL( LW )
WORK( 1 ) = SROUNDUP_LWORK( LWMIN )
END IF
*
IF( INFO.NE.0 ) THEN
Expand All @@ -261,7 +272,7 @@ SUBROUTINE CGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
*
* Quick return if possible
*
IF( MIN( M, N, K ).EQ.0 ) THEN
IF( MINMNK.EQ.0 ) THEN
RETURN
END IF
*
Expand All @@ -274,7 +285,7 @@ SUBROUTINE CGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
$ MB, C, LDC, WORK, LWORK, INFO )
END IF
*
WORK( 1 ) = REAL( LW )
WORK( 1 ) = SROUNDUP_LWORK( LWMIN )
*
RETURN
*
Expand Down
Loading

0 comments on commit 0814491

Please sign in to comment.