Skip to content

Commit 9dc2480

Browse files
dklyuchinskiydklyuchinskiy
committed
handle and document corner cases of lwork in lapack, align all precisions
1 parent 7d15f83 commit 9dc2480

Some content is hidden

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

96 files changed

+883
-688
lines changed

‎SRC/cgebrd.f

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -252,7 +252,7 @@ SUBROUTINE CGEBRD( M, N, A, LDA, D, E, TAUQ, TAUP, WORK, LWORK,
252252
ELSE
253253
LWKMIN = MAX( M, N )
254254
NB = MAX( 1, ILAENV( 1, 'CGEBRD', ' ', M, N, -1, -1 ) )
255-
LWKOPT = MAX( 1, ( M+N )*NB )
255+
LWKOPT = ( M+N )*NB
256256
END IF
257257
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
258258
LQUERY = ( LWORK.EQ.-1 )
@@ -292,7 +292,7 @@ SUBROUTINE CGEBRD( M, N, A, LDA, D, E, TAUQ, TAUP, WORK, LWORK,
292292
* Determine when to switch from blocked to unblocked code.
293293
*
294294
IF( NX.LT.MINMN ) THEN
295-
WS = ( M+N )*NB
295+
WS = LWKOPT
296296
IF( LWORK.LT.WS ) THEN
297297
*
298298
* Not enough work space for the optimal NB, consider using

‎SRC/cgehrd.f

Lines changed: 2 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -222,12 +222,12 @@ SUBROUTINE CGEHRD( N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO )
222222
INFO = -8
223223
END IF
224224
*
225+
NH = IHI - ILO + 1
225226
IF( INFO.EQ.0 ) THEN
226227
*
227228
* Compute the workspace requirements
228229
*
229-
230-
IF( N.EQ.0 ) THEN
230+
IF( NH.LE.1 ) THEN
231231
LWKOPT = 1
232232
ELSE
233233
NB = MIN( NBMAX, ILAENV( 1, 'DGEHRD', ' ', N, ILO, IHI,
@@ -255,7 +255,6 @@ SUBROUTINE CGEHRD( N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO )
255255
*
256256
* Quick return if possible
257257
*
258-
NH = IHI - ILO + 1
259258
IF( NH.LE.1 ) THEN
260259
WORK( 1 ) = 1
261260
RETURN

‎SRC/cgelqf.f

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -185,7 +185,7 @@ SUBROUTINE CGELQF( M, N, A, LDA, TAU, WORK, LWORK, INFO )
185185
INFO = -2
186186
ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
187187
INFO = -4
188-
ELSE IF ( .NOT.LQUERY ) THEN
188+
ELSE IF( .NOT.LQUERY ) THEN
189189
IF( LWORK.LE.0 .OR. ( N.GT.0 .AND. LWORK.LT.MAX( 1, M ) ) )
190190
$ INFO = -7
191191
END IF

‎SRC/cgemlq.f

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -110,8 +110,8 @@
110110
*>
111111
*> \param[out] WORK
112112
*> \verbatim
113-
*> (workspace) COMPLEX array, dimension (MAX(1,LWORK))
114-
*> On exit, if INFO = 0, WORK(1) returns the minimal LWORK.
113+
*> (workspace) COMPLEX array, dimension (MAX(1,LWORK))
114+
*> On exit, if INFO = 0, WORK(1) returns the minimal LWORK.
115115
*> \endverbatim
116116
*>
117117
*> \param[in] LWORK
@@ -227,7 +227,7 @@ SUBROUTINE CGEMLQ( SIDE, TRANS, M, N, K, A, LDA, T, TSIZE,
227227
ELSE
228228
LWMIN = MAX( 1, LW )
229229
END IF
230-
230+
*
231231
IF( ( NB.GT.K ) .AND. ( MN.GT.K ) ) THEN
232232
IF( MOD( MN - K, NB - K ) .EQ. 0 ) THEN
233233
NBLCKS = ( MN - K ) / ( NB - K )

‎SRC/cgemqr.f

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -111,8 +111,8 @@
111111
*>
112112
*> \param[out] WORK
113113
*> \verbatim
114-
*> (workspace) COMPLEX array, dimension (MAX(1,LWORK))
115-
*> On exit, if INFO = 0, WORK(1) returns the minimal LWORK.
114+
*> (workspace) COMPLEX array, dimension (MAX(1,LWORK))
115+
*> On exit, if INFO = 0, WORK(1) returns the minimal LWORK.
116116
*> \endverbatim
117117
*>
118118
*> \param[in] LWORK

‎SRC/cgeqlf.f

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -192,7 +192,7 @@ SUBROUTINE CGEQLF( M, N, A, LDA, TAU, WORK, LWORK, INFO )
192192
*
193193
IF( .NOT.LQUERY ) THEN
194194
IF( LWORK.LE.0 .OR. ( M.GT.0 .AND. LWORK.LT.MAX( 1, N ) ) )
195-
$ INFO = -7
195+
$ INFO = -7
196196
END IF
197197
END IF
198198
*

‎SRC/cgeqp3rk.f

Lines changed: 12 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -428,7 +428,7 @@
428428
*> \verbatim
429429
*> LWORK is INTEGER
430430
*> The dimension of the array WORK.
431-
*> LWORK >= 1, if MIN(M,N) = 0,
431+
*> LWORK >= 1, if MIN(M,N) = 0, and
432432
*> LWORK >= N+NRHS-1, otherwise.
433433
*> For optimal performance LWORK >= NB*( N+NRHS+1 ),
434434
*> where NB is the optimal block size for CGEQP3RK returned
@@ -628,8 +628,9 @@ SUBROUTINE CGEQP3RK( M, N, NRHS, KMAX, ABSTOL, RELTOL, A, LDA,
628628
* .. External Functions ..
629629
LOGICAL SISNAN
630630
INTEGER ISAMAX, ILAENV
631-
REAL SLAMCH, SCNRM2
632-
EXTERNAL SISNAN, SLAMCH, SCNRM2, ISAMAX, ILAENV
631+
REAL SLAMCH, SCNRM2, SROUNDUP_LWORK
632+
EXTERNAL SISNAN, SLAMCH, SCNRM2, ISAMAX, ILAENV,
633+
$ SROUNDUP_LWORK
633634
* ..
634635
* .. Intrinsic Functions ..
635636
INTRINSIC CMPLX, MAX, MIN
@@ -704,7 +705,7 @@ SUBROUTINE CGEQP3RK( M, N, NRHS, KMAX, ABSTOL, RELTOL, A, LDA,
704705
*
705706
LWKOPT = 2*N + NB*( N+NRHS+1 )
706707
END IF
707-
WORK( 1 ) = CMPLX( LWKOPT )
708+
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
708709
*
709710
IF( ( LWORK.LT.IWS ) .AND. .NOT.LQUERY ) THEN
710711
INFO = -15
@@ -727,7 +728,7 @@ SUBROUTINE CGEQP3RK( M, N, NRHS, KMAX, ABSTOL, RELTOL, A, LDA,
727728
K = 0
728729
MAXC2NRMK = ZERO
729730
RELMAXC2NRMK = ZERO
730-
WORK( 1 ) = CMPLX( LWKOPT )
731+
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
731732
RETURN
732733
END IF
733734
*
@@ -779,7 +780,7 @@ SUBROUTINE CGEQP3RK( M, N, NRHS, KMAX, ABSTOL, RELTOL, A, LDA,
779780
*
780781
* Array TAU is not set and contains undefined elements.
781782
*
782-
WORK( 1 ) = CMPLX( LWKOPT )
783+
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
783784
RETURN
784785
END IF
785786
*
@@ -798,7 +799,7 @@ SUBROUTINE CGEQP3RK( M, N, NRHS, KMAX, ABSTOL, RELTOL, A, LDA,
798799
TAU( J ) = CZERO
799800
END DO
800801
*
801-
WORK( 1 ) = CMPLX( LWKOPT )
802+
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
802803
RETURN
803804
*
804805
END IF
@@ -829,7 +830,7 @@ SUBROUTINE CGEQP3RK( M, N, NRHS, KMAX, ABSTOL, RELTOL, A, LDA,
829830
DO J = 1, MINMN
830831
TAU( J ) = CZERO
831832
END DO
832-
WORK( 1 ) = CMPLX( LWKOPT )
833+
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
833834
RETURN
834835
END IF
835836
*
@@ -874,7 +875,7 @@ SUBROUTINE CGEQP3RK( M, N, NRHS, KMAX, ABSTOL, RELTOL, A, LDA,
874875
TAU( J ) = CZERO
875876
END DO
876877
*
877-
WORK( 1 ) = CMPLX( LWKOPT )
878+
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
878879
RETURN
879880
END IF
880881
*
@@ -992,7 +993,7 @@ SUBROUTINE CGEQP3RK( M, N, NRHS, KMAX, ABSTOL, RELTOL, A, LDA,
992993
*
993994
* Return from the routine.
994995
*
995-
WORK( 1 ) = CMPLX( LWKOPT )
996+
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
996997
*
997998
RETURN
998999
*
@@ -1083,7 +1084,7 @@ SUBROUTINE CGEQP3RK( M, N, NRHS, KMAX, ABSTOL, RELTOL, A, LDA,
10831084
*
10841085
END IF
10851086
*
1086-
WORK( 1 ) = CMPLX( LWKOPT )
1087+
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
10871088
*
10881089
RETURN
10891090
*

‎SRC/cgeqrfp.f

Lines changed: 5 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -184,12 +184,12 @@ SUBROUTINE CGEQRFP( M, N, A, LDA, TAU, WORK, LWORK, INFO )
184184
INFO = 0
185185
NB = ILAENV( 1, 'CGEQRF', ' ', M, N, -1, -1 )
186186
K = MIN( M, N )
187-
IF ( K.EQ.0 ) THEN
188-
LWKMIN = 1
189-
LWKOPT = 1
187+
IF( K.EQ.0 ) THEN
188+
LWKMIN = 1
189+
LWKOPT = 1
190190
ELSE
191-
LWKMIN = N
192-
LWKOPT = N*NB
191+
LWKMIN = N
192+
LWKOPT = N*NB
193193
END IF
194194
WORK( 1 ) = SROUNDUP_LWORK( LWKOPT )
195195
*

‎SRC/cgesvj.f

Lines changed: 21 additions & 18 deletions
Original file line numberDiff line numberDiff line change
@@ -208,16 +208,17 @@
208208
*> \verbatim
209209
*> CWORK is COMPLEX array, dimension (max(1,LWORK))
210210
*> Used as workspace.
211-
*> If on entry LWORK = -1, then a workspace query is assumed and
212-
*> no computation is done; CWORK(1) is set to the minial (and optimal)
213-
*> length of CWORK.
214211
*> \endverbatim
215212
*>
216213
*> \param[in] LWORK
217214
*> \verbatim
218215
*> LWORK is INTEGER.
219216
*> Length of CWORK.
220-
*> LWORK >= 1, if MIN(M,N) = 0, and LWORK >= MAX(1,M+N), otherwise.
217+
*> LWORK >= 1, if MIN(M,N) = 0, and LWORK >= M+N, otherwise.
218+
*>
219+
*> If on entry LWORK = -1, then a workspace query is assumed and
220+
*> no computation is done; CWORK(1) is set to the minial (and optimal)
221+
*> length of CWORK.
221222
*> \endverbatim
222223
*>
223224
*> \param[in,out] RWORK
@@ -248,15 +249,17 @@
248249
*> RWORK(6) = the largest absolute value over all sines of the
249250
*> Jacobi rotation angles in the last sweep. It can be
250251
*> useful for a post festum analysis.
251-
*> If on entry LRWORK = -1, then a workspace query is assumed and
252-
*> no computation is done; RWORK(1) is set to the minial (and optimal)
253-
*> length of RWORK.
254252
*> \endverbatim
255253
*>
256254
*> \param[in] LRWORK
257255
*> \verbatim
258256
*> LRWORK is INTEGER
259-
*> Length of RWORK, LRWORK >= MAX(6,N).
257+
*> Length of RWORK.
258+
*> LRWORK >= 1, if MIN(M,N) = 0, and LRWORK >= MAX(6,N), otherwise
259+
*>
260+
*> If on entry LRWORK = -1, then a workspace query is assumed and
261+
*> no computation is done; RWORK(1) is set to the minial (and optimal)
262+
*> length of RWORK.
260263
*> \endverbatim
261264
*>
262265
*> \param[out] INFO
@@ -400,8 +403,8 @@ SUBROUTINE CGESVJ( JOBA, JOBU, JOBV, M, N, A, LDA, SVA, MV, V,
400403
INTEGER ISAMAX
401404
EXTERNAL ISAMAX
402405
* from LAPACK
403-
REAL SLAMCH
404-
EXTERNAL SLAMCH
406+
REAL SLAMCH, SROUNDUP_LWORK
407+
EXTERNAL SLAMCH, SROUNDUP_LWORK
405408
LOGICAL LSAME
406409
EXTERNAL LSAME
407410
* ..
@@ -423,19 +426,17 @@ SUBROUTINE CGESVJ( JOBA, JOBU, JOBV, M, N, A, LDA, SVA, MV, V,
423426
APPLV = LSAME( JOBV, 'A' )
424427
UPPER = LSAME( JOBA, 'U' )
425428
LOWER = LSAME( JOBA, 'L' )
426-
429+
*
427430
MINMN = MIN( M, N )
428431
IF( MINMN.EQ.0 ) THEN
429-
LWMIN = 1
430-
LRWMIN = 6
432+
LWMIN = 1
433+
LRWMIN = 1
431434
ELSE
432-
LWMIN = M + N
435+
LWMIN = M + N
433436
LRWMIN = MAX( 6, N )
434437
END IF
435-
CWORK(1) = LWMIN
436-
RWORK(1) = LRWMIN
437438
*
438-
LQUERY = ( LWORK .EQ. -1 ) .OR. ( LRWORK .EQ. -1 )
439+
LQUERY = ( LWORK.EQ.-1 ) .OR. ( LRWORK.EQ.-1 )
439440
IF( .NOT.( UPPER .OR. LOWER .OR. LSAME( JOBA, 'G' ) ) ) THEN
440441
INFO = -1
441442
ELSE IF( .NOT.( LSVEC .OR. UCTOL .OR. LSAME( JOBU, 'N' ) ) ) THEN
@@ -467,7 +468,9 @@ SUBROUTINE CGESVJ( JOBA, JOBU, JOBV, M, N, A, LDA, SVA, MV, V,
467468
IF( INFO.NE.0 ) THEN
468469
CALL XERBLA( 'CGESVJ', -INFO )
469470
RETURN
470-
ELSE IF ( LQUERY ) THEN
471+
ELSE IF( LQUERY ) THEN
472+
CWORK( 1 ) = SROUNDUP_LWORK( LWMIN )
473+
RWORK( 1 ) = SROUNDUP_LWORK( LRWMIN )
471474
RETURN
472475
END IF
473476
*

‎SRC/cgetsqrhrt.f

Lines changed: 11 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -131,13 +131,15 @@
131131
*> \param[in] LWORK
132132
*> \verbatim
133133
*> The dimension of the array WORK.
134+
*> If MIN(M,N) = 0, LWORK >= 1, else
134135
*> LWORK >= MAX( 1, LWT + LW1, MAX( LWT+N*N+LW2, LWT+N*N+N ) ),
135136
*> where
136137
*> NUM_ALL_ROW_BLOCKS = CEIL((M-N)/(MB1-N)),
137138
*> NB1LOCAL = MIN(NB1,N).
138139
*> LWT = NUM_ALL_ROW_BLOCKS * N * NB1LOCAL,
139140
*> LW1 = NB1LOCAL * N,
140141
*> LW2 = NB1LOCAL * MAX( NB1LOCAL, ( N - NB1LOCAL ) ).
142+
*>
141143
*> If LWORK = -1, then a workspace query is assumed.
142144
*> The routine only calculates the optimal size of the WORK
143145
*> array, returns this value as the first entry of the WORK
@@ -200,6 +202,10 @@ SUBROUTINE CGETSQRHRT( M, N, MB1, NB1, NB2, A, LDA, T, LDT, WORK,
200202
INTEGER I, IINFO, J, LW1, LW2, LWT, LDWT, LWORKOPT,
201203
$ NB1LOCAL, NB2LOCAL, NUM_ALL_ROW_BLOCKS
202204
* ..
205+
* .. External Functions ..
206+
REAL SROUNDUP_LWORK
207+
EXTERNAL SROUNDUP_LWORK
208+
* ..
203209
* .. External Subroutines ..
204210
EXTERNAL CCOPY, CLATSQR, CUNGTSQR_ROW, CUNHR_COL,
205211
$ XERBLA
@@ -212,7 +218,7 @@ SUBROUTINE CGETSQRHRT( M, N, MB1, NB1, NB2, A, LDA, T, LDT, WORK,
212218
* Test the input arguments
213219
*
214220
INFO = 0
215-
LQUERY = ( LWORK.EQ.-1 )
221+
LQUERY = ( LWORK.EQ.-1 )
216222
IF( M.LT.0 ) THEN
217223
INFO = -1
218224
ELSE IF( N.LT.0 .OR. M.LT.N ) THEN
@@ -225,7 +231,7 @@ SUBROUTINE CGETSQRHRT( M, N, MB1, NB1, NB2, A, LDA, T, LDT, WORK,
225231
INFO = -5
226232
ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
227233
INFO = -7
228-
ELSE IF( LDT.LT.MAX( 1, MIN( NB2, N ) ) ) THEN
234+
ELSE IF( LDT.LT.MAX( 1, MIN( NB2, N ) ) ) THEN
229235
INFO = -9
230236
ELSE
231237
*
@@ -278,14 +284,14 @@ SUBROUTINE CGETSQRHRT( M, N, MB1, NB1, NB2, A, LDA, T, LDT, WORK,
278284
CALL XERBLA( 'CGETSQRHRT', -INFO )
279285
RETURN
280286
ELSE IF ( LQUERY ) THEN
281-
WORK( 1 ) = CMPLX( LWORKOPT )
287+
WORK( 1 ) = SROUNDUP_LWORK( LWORKOPT )
282288
RETURN
283289
END IF
284290
*
285291
* Quick return if possible
286292
*
287293
IF( MIN( M, N ).EQ.0 ) THEN
288-
WORK( 1 ) = CMPLX( LWORKOPT )
294+
WORK( 1 ) = SROUNDUP_LWORK( LWORKOPT )
289295
RETURN
290296
END IF
291297
*
@@ -342,7 +348,7 @@ SUBROUTINE CGETSQRHRT( M, N, MB1, NB1, NB2, A, LDA, T, LDT, WORK,
342348
END IF
343349
END DO
344350
*
345-
WORK( 1 ) = CMPLX( LWORKOPT )
351+
WORK( 1 ) = SROUNDUP_LWORK( LWORKOPT )
346352
RETURN
347353
*
348354
* End of CGETSQRHRT

0 commit comments

Comments
 (0)