Skip to content

Commit 9d9a6b7

Browse files
xGGQRCS: check for sufficiently large workspace
1 parent 2ad6db2 commit 9d9a6b7

File tree

4 files changed

+50
-4
lines changed

4 files changed

+50
-4
lines changed

SRC/cggqrcs.f

Lines changed: 14 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -374,7 +374,7 @@ RECURSIVE SUBROUTINE CGGQRCS( JOBU1, JOBU2, JOBX, M, N, P, L,
374374
* .. Local Scalars ..
375375
LOGICAL WANTU1, WANTU2, WANTX, LQUERY
376376
INTEGER I, J, K, K1, LMAX, IG, IG11, IG21, IG22,
377-
$ IVT, IVT12, LDG, LDX, LDVT, LWKOPT
377+
$ IVT, IVT12, LDG, LDX, LDVT, LWKMIN, LWKOPT
378378
REAL BASE, NAN, NORMA, NORMB, NORMG, TOL, ULP, UNFL,
379379
$ THETA, IOTA, W
380380
COMPLEX CNAN
@@ -477,26 +477,32 @@ RECURSIVE SUBROUTINE CGGQRCS( JOBU1, JOBU2, JOBX, M, N, P, L,
477477
* Compute workspace
478478
*
479479
IF( INFO.EQ.0 ) THEN
480+
LWKMIN = 0
480481
LWKOPT = 0
481482
*
482483
CALL CGEQP3( M + P, N, WORK( IG ), LDG, IWORK, WORK, WORK, -1,
483484
$ RWORK, INFO )
485+
LWKMIN = MAX( LWKMIN, N + 1 )
484486
LWKOPT = MAX( LWKOPT, INT( WORK( 1 ) ) )
485487
*
486488
CALL CUNGQR( M + P, LMAX, LMAX, WORK( IG ), LDG, WORK, WORK,
487489
$ -1, INFO )
490+
LWKMIN = MAX( LWKMIN, LMAX )
488491
LWKOPT = MAX( LWKOPT, INT( WORK( 1 ) ) )
489492
*
490493
CALL CUNCSD2BY1( JOBU1, JOBU2, JOBX, M + P, M, LMAX,
491494
$ WORK( IG ), LDG, WORK( IG ), LDG,
492495
$ ALPHA,
493496
$ U1, LDU1, U2, LDU2, WORK( IVT ), LDVT,
494497
$ WORK, -1, RWORK, LRWORK, IWORK, INFO )
498+
LWKMIN = MAX( LWKMIN, INT( WORK( 1 ) ) )
495499
LWKOPT = MAX( LWKOPT, INT( WORK( 1 ) ) )
496500
* The matrix (A, B) must be stored sequentially for CUNGQR
501+
LWKMIN = LWKMIN + IVT
497502
LWKOPT = LWKOPT + IVT
498503
* 2-by-1 CSD matrix V1 must be stored
499504
IF( WANTX ) THEN
505+
LWKMIN = LWKMIN + LDVT*N
500506
LWKOPT = LWKOPT + LDVT*N
501507
END IF
502508
* Adjust CUNCSD2BY1 LRWORK for case with maximum memory
@@ -509,6 +515,13 @@ RECURSIVE SUBROUTINE CGGQRCS( JOBU1, JOBU2, JOBX, M, N, P, L,
509515
$ - 8 * MAX( 0, MIN( M, P, N, M+P-N ) )
510516
$ + 8 * MIN( M, P, N )
511517
LRWKOPT = MAX( 2*N, LRWORK2BY1 )
518+
* Check workspace size
519+
IF( LWORK.LT.LWKMIN .AND. .NOT.LQUERY ) THEN
520+
INFO = -20
521+
END IF
522+
IF( LRWORK.LT.LRWKOPT .AND. .NOT.LQUERY ) THEN
523+
INFO = -22
524+
END IF
512525
*
513526
WORK( 1 ) = CMPLX( REAL( LWKOPT ), 0.0E0 )
514527
RWORK( 1 ) = REAL( LRWKOPT )

SRC/dggqrcs.f

Lines changed: 11 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -351,7 +351,7 @@ RECURSIVE SUBROUTINE DGGQRCS( JOBU1, JOBU2, JOBX, M, N, P, L,
351351
* .. Local Scalars ..
352352
LOGICAL WANTU1, WANTU2, WANTX, LQUERY
353353
INTEGER I, J, K, K1, LMAX, IG, IG11, IG21, IG22,
354-
$ IVT, IVT12, LDG, LDX, LDVT, LWKOPT
354+
$ IVT, IVT12, LDG, LDX, LDVT, LWKMIN, LWKOPT
355355
DOUBLE PRECISION BASE, NAN, NORMA, NORMB, NORMG, TOL, ULP, UNFL,
356356
$ THETA, IOTA, W
357357
* ..
@@ -452,28 +452,38 @@ RECURSIVE SUBROUTINE DGGQRCS( JOBU1, JOBU2, JOBX, M, N, P, L,
452452
* Compute workspace
453453
*
454454
IF( INFO.EQ.0 ) THEN
455+
LWKMIN = 0
455456
LWKOPT = 0
456457
*
457458
CALL DGEQP3( M + P, N, WORK( IG ), LDG, IWORK, ALPHA, WORK, -1,
458459
$ INFO )
460+
LWKMIN = MAX( LWKMIN, 3 * N + 1 )
459461
LWKOPT = MAX( LWKOPT, INT( WORK( 1 ) ) )
460462
*
461463
CALL DORGQR( M + P, LMAX, LMAX, WORK( IG ), LDG, ALPHA, WORK,
462464
$ -1, INFO )
465+
LWKMIN = MAX( LWKMIN, LMAX )
463466
LWKOPT = MAX( LWKOPT, INT( WORK( 1 ) ) )
464467
*
465468
CALL DORCSD2BY1( JOBU1, JOBU2, JOBX, M + P, M, LMAX,
466469
$ WORK( IG ), LDG, WORK( IG ), LDG,
467470
$ ALPHA,
468471
$ U1, LDU1, U2, LDU2, WORK( IVT ), LDVT,
469472
$ WORK, -1, IWORK, INFO )
473+
LWKMIN = MAX( LWKMIN, INT( WORK( 1 ) ) )
470474
LWKOPT = MAX( LWKOPT, INT( WORK( 1 ) ) )
471475
* The matrix (A, B) must be stored sequentially for DORGQR
476+
LWKMIN = LWKMIN + IVT
472477
LWKOPT = LWKOPT + IVT
473478
* 2-by-1 CSD matrix V1 must be stored
474479
IF( WANTX ) THEN
480+
LWKMIN = LWKMIN + LDVT*N
475481
LWKOPT = LWKOPT + LDVT*N
476482
END IF
483+
* Check for minimum workspace size
484+
IF( LWORK.LT.LWKMIN .AND. .NOT.LQUERY ) THEN
485+
INFO = -20
486+
END IF
477487
*
478488
WORK( 1 ) = DBLE( LWKOPT )
479489
END IF

SRC/sggqrcs.f

Lines changed: 11 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -351,7 +351,7 @@ RECURSIVE SUBROUTINE SGGQRCS( JOBU1, JOBU2, JOBX, M, N, P, L,
351351
* .. Local Scalars ..
352352
LOGICAL WANTU1, WANTU2, WANTX, LQUERY
353353
INTEGER I, J, K, K1, LMAX, IG, IG11, IG21, IG22,
354-
$ IVT, IVT12, LDG, LDX, LDVT, LWKOPT
354+
$ IVT, IVT12, LDG, LDX, LDVT, LWKMIN, LWKOPT
355355
REAL BASE, NAN, NORMA, NORMB, NORMG, TOL, ULP, UNFL,
356356
$ THETA, IOTA, W
357357
* ..
@@ -452,28 +452,38 @@ RECURSIVE SUBROUTINE SGGQRCS( JOBU1, JOBU2, JOBX, M, N, P, L,
452452
* Compute workspace
453453
*
454454
IF( INFO.EQ.0 ) THEN
455+
LWKMIN = 0
455456
LWKOPT = 0
456457
*
457458
CALL SGEQP3( M + P, N, WORK( IG ), LDG, IWORK, ALPHA, WORK, -1,
458459
$ INFO )
460+
LWKMIN = MAX( LWKMIN, 3 * N + 1 )
459461
LWKOPT = MAX( LWKOPT, INT( WORK( 1 ) ) )
460462
*
461463
CALL SORGQR( M + P, LMAX, LMAX, WORK( IG ), LDG, ALPHA, WORK,
462464
$ -1, INFO )
465+
LWKMIN = MAX( LWKMIN, LMAX )
463466
LWKOPT = MAX( LWKOPT, INT( WORK( 1 ) ) )
464467
*
465468
CALL SORCSD2BY1( JOBU1, JOBU2, JOBX, M + P, M, LMAX,
466469
$ WORK( IG ), LDG, WORK( IG ), LDG,
467470
$ ALPHA,
468471
$ U1, LDU1, U2, LDU2, WORK( IVT ), LDVT,
469472
$ WORK, -1, IWORK, INFO )
473+
LWKMIN = MAX( LWKMIN, INT( WORK( 1 ) ) )
470474
LWKOPT = MAX( LWKOPT, INT( WORK( 1 ) ) )
471475
* The matrix (A, B) must be stored sequentially for SORGQR
476+
LWKMIN = LWKMIN + IVT
472477
LWKOPT = LWKOPT + IVT
473478
* 2-by-1 CSD matrix V1 must be stored
474479
IF( WANTX ) THEN
480+
LWKMIN = LWKMIN + LDVT*N
475481
LWKOPT = LWKOPT + LDVT*N
476482
END IF
483+
* Check workspace size
484+
IF( LWORK.LT.LWKMIN .AND. .NOT.LQUERY ) THEN
485+
INFO = -20
486+
END IF
477487
*
478488
WORK( 1 ) = REAL( LWKOPT )
479489
END IF

SRC/zggqrcs.f

Lines changed: 14 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -374,7 +374,7 @@ RECURSIVE SUBROUTINE ZGGQRCS( JOBU1, JOBU2, JOBX, M, N, P, L,
374374
* .. Local Scalars ..
375375
LOGICAL WANTU1, WANTU2, WANTX, LQUERY
376376
INTEGER I, J, K, K1, LMAX, IG, IG11, IG21, IG22,
377-
$ IVT, IVT12, LDG, LDX, LDVT, LWKOPT
377+
$ IVT, IVT12, LDG, LDX, LDVT, LWKMIN, LWKOPT
378378
DOUBLE PRECISION BASE, NAN, NORMA, NORMB, NORMG, TOL, ULP, UNFL,
379379
$ THETA, IOTA, W
380380
COMPLEX*16 ZNAN
@@ -477,26 +477,32 @@ RECURSIVE SUBROUTINE ZGGQRCS( JOBU1, JOBU2, JOBX, M, N, P, L,
477477
* Compute workspace
478478
*
479479
IF( INFO.EQ.0 ) THEN
480+
LWKMIN = 0
480481
LWKOPT = 0
481482
*
482483
CALL ZGEQP3( M + P, N, WORK( IG ), LDG, IWORK, WORK, WORK, -1,
483484
$ RWORK, INFO )
485+
LWKMIN = MAX( LWKMIN, N + 1 )
484486
LWKOPT = MAX( LWKOPT, INT( WORK( 1 ) ) )
485487
*
486488
CALL ZUNGQR( M + P, LMAX, LMAX, WORK( IG ), LDG, WORK, WORK,
487489
$ -1, INFO )
490+
LWKMIN = MAX( LWKMIN, LMAX )
488491
LWKOPT = MAX( LWKOPT, INT( WORK( 1 ) ) )
489492
*
490493
CALL ZUNCSD2BY1( JOBU1, JOBU2, JOBX, M + P, M, LMAX,
491494
$ WORK( IG ), LDG, WORK( IG ), LDG,
492495
$ ALPHA,
493496
$ U1, LDU1, U2, LDU2, WORK( IVT ), LDVT,
494497
$ WORK, -1, RWORK, LRWORK, IWORK, INFO )
498+
LWKMIN = MAX( LWKMIN, INT( WORK( 1 ) ) )
495499
LWKOPT = MAX( LWKOPT, INT( WORK( 1 ) ) )
496500
* The matrix (A, B) must be stored sequentially for ZUNGQR
501+
LWKMIN = LWKMIN + IVT
497502
LWKOPT = LWKOPT + IVT
498503
* 2-by-1 CSD matrix V1 must be stored
499504
IF( WANTX ) THEN
505+
LWKMIN = LWKMIN + LDVT*N
500506
LWKOPT = LWKOPT + LDVT*N
501507
END IF
502508
* Adjust ZUNCSD2BY1 LRWORK for case with maximum memory
@@ -509,6 +515,13 @@ RECURSIVE SUBROUTINE ZGGQRCS( JOBU1, JOBU2, JOBX, M, N, P, L,
509515
$ - 8 * MAX( 0, MIN( M, P, N, M+P-N ) )
510516
$ + 8 * MIN( M, P, N )
511517
LRWKOPT = MAX( 2*N, LRWORK2BY1 )
518+
* Check workspace size
519+
IF( LWORK.LT.LWKMIN .AND. .NOT.LQUERY ) THEN
520+
INFO = -20
521+
END IF
522+
IF( LRWORK.LT.LRWKOPT .AND. .NOT.LQUERY ) THEN
523+
INFO = -22
524+
END IF
512525
*
513526
WORK( 1 ) = DCMPLX( DBLE( LWKOPT ), 0.0D0 )
514527
RWORK( 1 ) = DBLE( LRWKOPT )

0 commit comments

Comments
 (0)