Skip to content

Commit 45ef0d7

Browse files
martin-frbgmartin-frbg
authored
Handle corner cases of LWORK (Reference-LAPACK PR 942)
1 parent c082669 commit 45ef0d7

42 files changed

Lines changed: 733 additions & 504 deletions

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

‎lapack-netlib/SRC/zgebrd.f‎

Lines changed: 17 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -122,7 +122,8 @@
122122
*> \param[in] LWORK
123123
*> \verbatim
124124
*> LWORK is INTEGER
125-
*> The length of the array WORK. LWORK >= max(1,M,N).
125+
*> The length of the array WORK.
126+
*> LWORK >= 1, if MIN(M,N) = 0, and LWORK >= MAX(M,N), otherwise.
126127
*> For optimum performance LWORK >= (M+N)*NB, where NB
127128
*> is the optimal blocksize.
128129
*>
@@ -147,7 +148,7 @@
147148
*> \author Univ. of Colorado Denver
148149
*> \author NAG Ltd.
149150
*
150-
*> \ingroup complex16GEcomputational
151+
*> \ingroup gebrd
151152
*
152153
*> \par Further Details:
153154
* =====================
@@ -223,8 +224,8 @@ SUBROUTINE ZGEBRD( M, N, A, LDA, D, E, TAUQ, TAUP, WORK, LWORK,
223224
* ..
224225
* .. Local Scalars ..
225226
LOGICAL LQUERY
226-
INTEGER I, IINFO, J, LDWRKX, LDWRKY, LWKOPT, MINMN, NB,
227-
$ NBMIN, NX, WS
227+
INTEGER I, IINFO, J, LDWRKX, LDWRKY, LWKMIN, LWKOPT,
228+
$ MINMN, NB, NBMIN, NX, WS
228229
* ..
229230
* .. External Subroutines ..
230231
EXTERNAL XERBLA, ZGEBD2, ZGEMM, ZLABRD
@@ -241,17 +242,25 @@ SUBROUTINE ZGEBRD( M, N, A, LDA, D, E, TAUQ, TAUP, WORK, LWORK,
241242
* Test the input parameters
242243
*
243244
INFO = 0
244-
NB = MAX( 1, ILAENV( 1, 'ZGEBRD', ' ', M, N, -1, -1 ) )
245-
LWKOPT = ( M+N )*NB
245+
MINMN = MIN( M, N )
246+
IF( MINMN.EQ.0 ) THEN
247+
LWKMIN = 1
248+
LWKOPT = 1
249+
ELSE
250+
LWKMIN = MAX( M, N )
251+
NB = MAX( 1, ILAENV( 1, 'ZGEBRD', ' ', M, N, -1, -1 ) )
252+
LWKOPT = ( M+N )*NB
253+
END IF
246254
WORK( 1 ) = DBLE( LWKOPT )
255+
*
247256
LQUERY = ( LWORK.EQ.-1 )
248257
IF( M.LT.0 ) THEN
249258
INFO = -1
250259
ELSE IF( N.LT.0 ) THEN
251260
INFO = -2
252261
ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
253262
INFO = -4
254-
ELSE IF( LWORK.LT.MAX( 1, M, N ) .AND. .NOT.LQUERY ) THEN
263+
ELSE IF( LWORK.LT.LWKMIN .AND. .NOT.LQUERY ) THEN
255264
INFO = -10
256265
END IF
257266
IF( INFO.LT.0 ) THEN
@@ -263,7 +272,6 @@ SUBROUTINE ZGEBRD( M, N, A, LDA, D, E, TAUQ, TAUP, WORK, LWORK,
263272
*
264273
* Quick return if possible
265274
*
266-
MINMN = MIN( M, N )
267275
IF( MINMN.EQ.0 ) THEN
268276
WORK( 1 ) = 1
269277
RETURN
@@ -282,7 +290,7 @@ SUBROUTINE ZGEBRD( M, N, A, LDA, D, E, TAUQ, TAUP, WORK, LWORK,
282290
* Determine when to switch from blocked to unblocked code.
283291
*
284292
IF( NX.LT.MINMN ) THEN
285-
WS = ( M+N )*NB
293+
WS = LWKOPT
286294
IF( LWORK.LT.WS ) THEN
287295
*
288296
* Not enough work space for the optimal NB, consider using

‎lapack-netlib/SRC/zgehrd.f‎

Lines changed: 14 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -89,7 +89,7 @@
8989
*>
9090
*> \param[out] WORK
9191
*> \verbatim
92-
*> WORK is COMPLEX*16 array, dimension (LWORK)
92+
*> WORK is COMPLEX*16 array, dimension (MAX(1,LWORK))
9393
*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
9494
*> \endverbatim
9595
*>
@@ -120,7 +120,7 @@
120120
*> \author Univ. of Colorado Denver
121121
*> \author NAG Ltd.
122122
*
123-
*> \ingroup complex16GEcomputational
123+
*> \ingroup gehrd
124124
*
125125
*> \par Further Details:
126126
* =====================
@@ -173,7 +173,7 @@ SUBROUTINE ZGEHRD( N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO )
173173
INTEGER IHI, ILO, INFO, LDA, LWORK, N
174174
* ..
175175
* .. Array Arguments ..
176-
COMPLEX*16 A( LDA, * ), TAU( * ), WORK( * )
176+
COMPLEX*16 A( LDA, * ), TAU( * ), WORK( * )
177177
* ..
178178
*
179179
* =====================================================================
@@ -182,15 +182,15 @@ SUBROUTINE ZGEHRD( N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO )
182182
INTEGER NBMAX, LDT, TSIZE
183183
PARAMETER ( NBMAX = 64, LDT = NBMAX+1,
184184
$ TSIZE = LDT*NBMAX )
185-
COMPLEX*16 ZERO, ONE
185+
COMPLEX*16 ZERO, ONE
186186
PARAMETER ( ZERO = ( 0.0D+0, 0.0D+0 ),
187187
$ ONE = ( 1.0D+0, 0.0D+0 ) )
188188
* ..
189189
* .. Local Scalars ..
190190
LOGICAL LQUERY
191191
INTEGER I, IB, IINFO, IWT, J, LDWORK, LWKOPT, NB,
192192
$ NBMIN, NH, NX
193-
COMPLEX*16 EI
193+
COMPLEX*16 EI
194194
* ..
195195
* .. External Subroutines ..
196196
EXTERNAL ZAXPY, ZGEHD2, ZGEMM, ZLAHR2, ZLARFB, ZTRMM,
@@ -221,12 +221,18 @@ SUBROUTINE ZGEHRD( N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO )
221221
INFO = -8
222222
END IF
223223
*
224+
NH = IHI - ILO + 1
224225
IF( INFO.EQ.0 ) THEN
225226
*
226227
* Compute the workspace requirements
227228
*
228-
NB = MIN( NBMAX, ILAENV( 1, 'ZGEHRD', ' ', N, ILO, IHI, -1 ) )
229-
LWKOPT = N*NB + TSIZE
229+
IF( NH.LE.1 ) THEN
230+
LWKOPT = 1
231+
ELSE
232+
NB = MIN( NBMAX, ILAENV( 1, 'ZGEHRD', ' ', N, ILO, IHI,
233+
$ -1 ) )
234+
LWKOPT = N*NB + TSIZE
235+
END IF
230236
WORK( 1 ) = LWKOPT
231237
ENDIF
232238
*
@@ -248,7 +254,6 @@ SUBROUTINE ZGEHRD( N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO )
248254
*
249255
* Quick return if possible
250256
*
251-
NH = IHI - ILO + 1
252257
IF( NH.LE.1 ) THEN
253258
WORK( 1 ) = 1
254259
RETURN
@@ -268,7 +273,7 @@ SUBROUTINE ZGEHRD( N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO )
268273
*
269274
* Determine if workspace is large enough for blocked code
270275
*
271-
IF( LWORK.LT.N*NB+TSIZE ) THEN
276+
IF( LWORK.LT.LWKOPT ) THEN
272277
*
273278
* Not enough workspace to use optimal NB: determine the
274279
* minimum value of NB, and reduce NB or force use of

‎lapack-netlib/SRC/zgelq.f‎

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -98,7 +98,7 @@
9898
*> \param[in] LWORK
9999
*> \verbatim
100100
*> LWORK is INTEGER
101-
*> The dimension of the array WORK.
101+
*> The dimension of the array WORK. LWORK >= 1.
102102
*> If LWORK = -1 or -2, then a workspace query is assumed. The routine
103103
*> only calculates the sizes of the T and WORK arrays, returns these
104104
*> values as the first entries of the T and WORK arrays, and no error
@@ -166,6 +166,8 @@
166166
*> the LQ factorization.
167167
*> \endverbatim
168168
*>
169+
*> \ingroup gelq
170+
*>
169171
* =====================================================================
170172
SUBROUTINE ZGELQ( M, N, A, LDA, T, TSIZE, WORK, LWORK,
171173
$ INFO )

‎lapack-netlib/SRC/zgelqf.f‎

Lines changed: 13 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -93,7 +93,8 @@
9393
*> \param[in] LWORK
9494
*> \verbatim
9595
*> LWORK is INTEGER
96-
*> The dimension of the array WORK. LWORK >= max(1,M).
96+
*> The dimension of the array WORK.
97+
*> LWORK >= 1, if MIN(M,N) = 0, and LWORK >= M, otherwise.
9798
*> For optimum performance LWORK >= M*NB, where NB is the
9899
*> optimal blocksize.
99100
*>
@@ -118,7 +119,7 @@
118119
*> \author Univ. of Colorado Denver
119120
*> \author NAG Ltd.
120121
*
121-
*> \ingroup complex16GEcomputational
122+
*> \ingroup gelqf
122123
*
123124
*> \par Further Details:
124125
* =====================
@@ -174,29 +175,34 @@ SUBROUTINE ZGELQF( M, N, A, LDA, TAU, WORK, LWORK, INFO )
174175
* Test the input arguments
175176
*
176177
INFO = 0
178+
K = MIN( M, N )
177179
NB = ILAENV( 1, 'ZGELQF', ' ', M, N, -1, -1 )
178-
LWKOPT = M*NB
179-
WORK( 1 ) = LWKOPT
180180
LQUERY = ( LWORK.EQ.-1 )
181181
IF( M.LT.0 ) THEN
182182
INFO = -1
183183
ELSE IF( N.LT.0 ) THEN
184184
INFO = -2
185185
ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
186186
INFO = -4
187-
ELSE IF( LWORK.LT.MAX( 1, M ) .AND. .NOT.LQUERY ) THEN
188-
INFO = -7
187+
ELSE IF( .NOT.LQUERY ) THEN
188+
IF( LWORK.LE.0 .OR. ( N.GT.0 .AND. LWORK.LT.MAX( 1, M ) ) )
189+
$ INFO = -7
189190
END IF
190191
IF( INFO.NE.0 ) THEN
191192
CALL XERBLA( 'ZGELQF', -INFO )
192193
RETURN
193194
ELSE IF( LQUERY ) THEN
195+
IF( K.EQ.0 ) THEN
196+
LWKOPT = 1
197+
ELSE
198+
LWKOPT = M*NB
199+
END IF
200+
WORK( 1 ) = LWKOPT
194201
RETURN
195202
END IF
196203
*
197204
* Quick return if possible
198205
*
199-
K = MIN( M, N )
200206
IF( K.EQ.0 ) THEN
201207
WORK( 1 ) = 1
202208
RETURN

‎lapack-netlib/SRC/zgemlq.f‎

Lines changed: 19 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -109,16 +109,17 @@
109109
*>
110110
*> \param[out] WORK
111111
*> \verbatim
112-
*> (workspace) COMPLEX*16 array, dimension (MAX(1,LWORK))
112+
*> (workspace) COMPLEX*16 array, dimension (MAX(1,LWORK))
113+
*> On exit, if INFO = 0, WORK(1) returns the minimal LWORK.
113114
*> \endverbatim
114115
*>
115116
*> \param[in] LWORK
116117
*> \verbatim
117118
*> LWORK is INTEGER
118-
*> The dimension of the array WORK.
119+
*> The dimension of the array WORK. LWORK >= 1.
119120
*> If LWORK = -1, then a workspace query is assumed. The routine
120121
*> only calculates the size of the WORK array, returns this
121-
*> value as WORK(1), and no error message related to WORK
122+
*> value as WORK(1), and no error message related to WORK
122123
*> is issued by XERBLA.
123124
*> \endverbatim
124125
*>
@@ -142,7 +143,7 @@
142143
*>
143144
*> \verbatim
144145
*>
145-
*> These details are particular for this LAPACK implementation. Users should not
146+
*> These details are particular for this LAPACK implementation. Users should not
146147
*> take them for granted. These details may change in the future, and are not likely
147148
*> true for another LAPACK implementation. These details are relevant if one wants
148149
*> to try to understand the code. They are not part of the interface.
@@ -158,11 +159,13 @@
158159
*> block sizes MB and NB returned by ILAENV, ZGELQ will use either
159160
*> ZLASWLQ (if the matrix is wide-and-short) or ZGELQT to compute
160161
*> the LQ factorization.
161-
*> This version of ZGEMLQ will use either ZLAMSWLQ or ZGEMLQT to
162+
*> This version of ZGEMLQ will use either ZLAMSWLQ or ZGEMLQT to
162163
*> multiply matrix Q by another matrix.
163164
*> Further Details in ZLAMSWLQ or ZGEMLQT.
164165
*> \endverbatim
165166
*>
167+
*> \ingroup gemlq
168+
*>
166169
* =====================================================================
167170
SUBROUTINE ZGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
168171
$ C, LDC, WORK, LWORK, INFO )
@@ -184,7 +187,7 @@ SUBROUTINE ZGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
184187
* ..
185188
* .. Local Scalars ..
186189
LOGICAL LEFT, RIGHT, TRAN, NOTRAN, LQUERY
187-
INTEGER MB, NB, LW, NBLCKS, MN
190+
INTEGER MB, NB, LW, NBLCKS, MN, MINMNK, LWMIN
188191
* ..
189192
* .. External Functions ..
190193
LOGICAL LSAME
@@ -200,7 +203,7 @@ SUBROUTINE ZGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
200203
*
201204
* Test the input arguments
202205
*
203-
LQUERY = LWORK.EQ.-1
206+
LQUERY = ( LWORK.EQ.-1 )
204207
NOTRAN = LSAME( TRANS, 'N' )
205208
TRAN = LSAME( TRANS, 'C' )
206209
LEFT = LSAME( SIDE, 'L' )
@@ -215,6 +218,13 @@ SUBROUTINE ZGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
215218
LW = M * MB
216219
MN = N
217220
END IF
221+
*
222+
MINMNK = MIN( M, N, K )
223+
IF( MINMNK.EQ.0 ) THEN
224+
LWMIN = 1
225+
ELSE
226+
LWMIN = MAX( 1, LW )
227+
END IF
218228
*
219229
IF( ( NB.GT.K ) .AND. ( MN.GT.K ) ) THEN
220230
IF( MOD( MN - K, NB - K ) .EQ. 0 ) THEN
@@ -243,7 +253,7 @@ SUBROUTINE ZGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
243253
INFO = -9
244254
ELSE IF( LDC.LT.MAX( 1, M ) ) THEN
245255
INFO = -11
246-
ELSE IF( ( LWORK.LT.MAX( 1, LW ) ) .AND. ( .NOT.LQUERY ) ) THEN
256+
ELSE IF( ( LWORK.LT.LWMIN ) .AND. ( .NOT.LQUERY ) ) THEN
247257
INFO = -13
248258
END IF
249259
*
@@ -260,7 +270,7 @@ SUBROUTINE ZGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
260270
*
261271
* Quick return if possible
262272
*
263-
IF( MIN( M, N, K ).EQ.0 ) THEN
273+
IF( MINMNK.EQ.0 ) THEN
264274
RETURN
265275
END IF
266276
*

0 commit comments

Comments
 (0)