Skip to content

Commit 8f8a45e

Browse files
authored
Avoid intermediate integer overflow in ?GESVJ/?GSVJ0 (Reference-LAPACK PR 1454) (#6205)
* Avoid integer overflow in calculation of EMPTSW (Reference-LAPACK PR 1454) * Add tests for ?GESVJ integer overflow protection (Reference-LAPACK PR 1454)
1 parent 0440606 commit 8f8a45e

25 files changed

Lines changed: 901 additions & 739 deletions

‎Makefile‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -406,7 +406,7 @@ endif
406406

407407
lapack-test :
408408
(cd $(NETLIB_LAPACK_DIR)/TESTING && rm -f x* *.out)
409-
$(MAKE) -j 1 -C $(NETLIB_LAPACK_DIR)/TESTING/EIG xeigtstc xeigtstd xeigtsts xeigtstz test_cgges3_scaling test_sgges3_scaling test_dgges3_scaling test_zgges3_scaling
409+
$(MAKE) -j 1 -C $(NETLIB_LAPACK_DIR)/TESTING/EIG xeigtstc xeigtstd xeigtsts xeigtstz test_cgges3_scaling test_sgges3_scaling test_dgges3_scaling test_zgges3_scaling test_cgesvj_emptsw test_dgesvj_emptsw test_sgesvj_emptsw test_zgesvj_emptsw
410410

411411
$(MAKE) -j 1 -C $(NETLIB_LAPACK_DIR)/TESTING/LIN xlintstc xlintstd xlintstds xlintstrfd xlintstrfz xlintsts xlintstz xlintstzc xlintstrfs xlintstrfc
412412

‎lapack-netlib/SRC/cgesvj.c‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1087,7 +1087,7 @@ f"> */
10871087

10881088
/* Row-cyclic Jacobi SVD algorithm with column pivoting */
10891089

1090-
emptsw = *n * (*n - 1) / 2;
1090+
emptsw = (*n / 2) * (*n - 1 + *n % 2);
10911091
notrot = 0;
10921092
i__1 = *n;
10931093
for (q = 1; q <= i__1; ++q) {

‎lapack-netlib/SRC/cgesvj.f‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -714,7 +714,7 @@ SUBROUTINE CGESVJ( JOBA, JOBU, JOBV, M, N, A, LDA, SVA, MV, V,
714714
*
715715
* Row-cyclic Jacobi SVD algorithm with column pivoting
716716
*
717-
EMPTSW = ( N*( N-1 ) ) / 2
717+
EMPTSW = ( N / 2 )*( N-1 + MOD( N, 2 ) )
718718
NOTROT = 0
719719

720720
DO 1868 q = 1, N

‎lapack-netlib/SRC/cgsvj0.c‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -638,7 +638,7 @@ f"> */
638638
roottol = sqrt(*tol);
639639

640640

641-
emptsw = *n * (*n - 1) / 2;
641+
emptsw = (*n / 2) * (*n - 1 + *n % 2);
642642
notrot = 0;
643643

644644

‎lapack-netlib/SRC/cgsvj0.f‎

Lines changed: 32 additions & 19 deletions
Original file line numberDiff line numberDiff line change
@@ -5,15 +5,13 @@
55
* Online html documentation available at
66
* http://www.netlib.org/lapack/explore-html/
77
*
8-
*> \htmlonly
98
*> Download CGSVJ0 + dependencies
109
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cgsvj0.f">
1110
*> [TGZ]</a>
1211
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cgsvj0.f">
1312
*> [ZIP]</a>
1413
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cgsvj0.f">
1514
*> [TXT]</a>
16-
*> \endhtmlonly
1715
*
1816
* Definition:
1917
* ===========
@@ -193,7 +191,7 @@
193191
*> \author Univ. of Colorado Denver
194192
*> \author NAG Ltd.
195193
*
196-
*> \ingroup complexOTHERcomputational
194+
*> \ingroup gsvj0
197195
*
198196
*> \par Further Details:
199197
* =====================
@@ -321,7 +319,7 @@ SUBROUTINE CGSVJ0( JOBV, M, N, A, LDA, D, SVA, MV, V, LDV, EPS,
321319
*
322320
* .. Row-cyclic Jacobi SVD algorithm with column pivoting ..
323321
*
324-
EMPTSW = ( N*( N-1 ) ) / 2
322+
EMPTSW = ( N / 2 )*( N-1 + MOD( N, 2 ) )
325323
NOTROT = 0
326324
*
327325
* .. Row-cyclic pivot strategy with de Rijk's pivoting ..
@@ -499,7 +497,8 @@ SUBROUTINE CGSVJ0( JOBV, M, N, A, LDA, D, SVA, MV, V, LDV, EPS,
499497
T = HALF / THETA
500498
CS = ONE
501499

502-
CALL CROT( M, A(1,p), 1, A(1,q), 1,
500+
CALL CROT( M, A(1,p), 1, A(1,q),
501+
$ 1,
503502
$ CS, CONJG(OMPQ)*T )
504503
IF ( RSVEC ) THEN
505504
CALL CROT( MVL, V(1,p), 1,
@@ -528,7 +527,8 @@ SUBROUTINE CGSVJ0( JOBV, M, N, A, LDA, D, SVA, MV, V, LDV, EPS,
528527
AAPP = AAPP*SQRT( MAX( ZERO,
529528
$ ONE-T*AQOAP*AAPQ1 ) )
530529
*
531-
CALL CROT( M, A(1,p), 1, A(1,q), 1,
530+
CALL CROT( M, A(1,p), 1, A(1,q),
531+
$ 1,
532532
$ CS, CONJG(OMPQ)*SN )
533533
IF ( RSVEC ) THEN
534534
CALL CROT( MVL, V(1,p), 1,
@@ -541,14 +541,17 @@ SUBROUTINE CGSVJ0( JOBV, M, N, A, LDA, D, SVA, MV, V, LDV, EPS,
541541
* .. have to use modified Gram-Schmidt like transformation
542542
CALL CCOPY( M, A( 1, p ), 1,
543543
$ WORK, 1 )
544-
CALL CLASCL( 'G', 0, 0, AAPP, ONE, M,
544+
CALL CLASCL( 'G', 0, 0, AAPP, ONE,
545+
$ M,
545546
$ 1, WORK, LDA,
546547
$ IERR )
547-
CALL CLASCL( 'G', 0, 0, AAQQ, ONE, M,
548+
CALL CLASCL( 'G', 0, 0, AAQQ, ONE,
549+
$ M,
548550
$ 1, A( 1, q ), LDA, IERR )
549551
CALL CAXPY( M, -AAPQ, WORK, 1,
550552
$ A( 1, q ), 1 )
551-
CALL CLASCL( 'G', 0, 0, ONE, AAQQ, M,
553+
CALL CLASCL( 'G', 0, 0, ONE, AAQQ,
554+
$ M,
552555
$ 1, A( 1, q ), LDA, IERR )
553556
SVA( q ) = AAQQ*SQRT( MAX( ZERO,
554557
$ ONE-AAPQ1*AAPQ1 ) )
@@ -563,7 +566,8 @@ SUBROUTINE CGSVJ0( JOBV, M, N, A, LDA, D, SVA, MV, V, LDV, EPS,
563566
$ THEN
564567
IF( ( AAQQ.LT.ROOTBIG ) .AND.
565568
$ ( AAQQ.GT.ROOTSFMIN ) ) THEN
566-
SVA( q ) = SCNRM2( M, A( 1, q ), 1 )
569+
SVA( q ) = SCNRM2( M, A( 1, q ),
570+
$ 1 )
567571
ELSE
568572
T = ZERO
569573
AAQQ = ONE
@@ -715,7 +719,8 @@ SUBROUTINE CGSVJ0( JOBV, M, N, A, LDA, D, SVA, MV, V, LDV, EPS,
715719
IF( ABS( THETA ).GT.BIGTHETA ) THEN
716720
T = HALF / THETA
717721
CS = ONE
718-
CALL CROT( M, A(1,p), 1, A(1,q), 1,
722+
CALL CROT( M, A(1,p), 1, A(1,q),
723+
$ 1,
719724
$ CS, CONJG(OMPQ)*T )
720725
IF( RSVEC ) THEN
721726
CALL CROT( MVL, V(1,p), 1,
@@ -742,7 +747,8 @@ SUBROUTINE CGSVJ0( JOBV, M, N, A, LDA, D, SVA, MV, V, LDV, EPS,
742747
AAPP = AAPP*SQRT( MAX( ZERO,
743748
$ ONE-T*AQOAP*AAPQ1 ) )
744749
*
745-
CALL CROT( M, A(1,p), 1, A(1,q), 1,
750+
CALL CROT( M, A(1,p), 1, A(1,q),
751+
$ 1,
746752
$ CS, CONJG(OMPQ)*SN )
747753
IF( RSVEC ) THEN
748754
CALL CROT( MVL, V(1,p), 1,
@@ -756,15 +762,18 @@ SUBROUTINE CGSVJ0( JOBV, M, N, A, LDA, D, SVA, MV, V, LDV, EPS,
756762
IF( AAPP.GT.AAQQ ) THEN
757763
CALL CCOPY( M, A( 1, p ), 1,
758764
$ WORK, 1 )
759-
CALL CLASCL( 'G', 0, 0, AAPP, ONE,
765+
CALL CLASCL( 'G', 0, 0, AAPP,
766+
$ ONE,
760767
$ M, 1, WORK,LDA,
761768
$ IERR )
762-
CALL CLASCL( 'G', 0, 0, AAQQ, ONE,
769+
CALL CLASCL( 'G', 0, 0, AAQQ,
770+
$ ONE,
763771
$ M, 1, A( 1, q ), LDA,
764772
$ IERR )
765773
CALL CAXPY( M, -AAPQ, WORK,
766774
$ 1, A( 1, q ), 1 )
767-
CALL CLASCL( 'G', 0, 0, ONE, AAQQ,
775+
CALL CLASCL( 'G', 0, 0, ONE,
776+
$ AAQQ,
768777
$ M, 1, A( 1, q ), LDA,
769778
$ IERR )
770779
SVA( q ) = AAQQ*SQRT( MAX( ZERO,
@@ -773,15 +782,18 @@ SUBROUTINE CGSVJ0( JOBV, M, N, A, LDA, D, SVA, MV, V, LDV, EPS,
773782
ELSE
774783
CALL CCOPY( M, A( 1, q ), 1,
775784
$ WORK, 1 )
776-
CALL CLASCL( 'G', 0, 0, AAQQ, ONE,
785+
CALL CLASCL( 'G', 0, 0, AAQQ,
786+
$ ONE,
777787
$ M, 1, WORK,LDA,
778788
$ IERR )
779-
CALL CLASCL( 'G', 0, 0, AAPP, ONE,
789+
CALL CLASCL( 'G', 0, 0, AAPP,
790+
$ ONE,
780791
$ M, 1, A( 1, p ), LDA,
781792
$ IERR )
782793
CALL CAXPY( M, -CONJG(AAPQ),
783794
$ WORK, 1, A( 1, p ), 1 )
784-
CALL CLASCL( 'G', 0, 0, ONE, AAPP,
795+
CALL CLASCL( 'G', 0, 0, ONE,
796+
$ AAPP,
785797
$ M, 1, A( 1, p ), LDA,
786798
$ IERR )
787799
SVA( p ) = AAPP*SQRT( MAX( ZERO,
@@ -797,7 +809,8 @@ SUBROUTINE CGSVJ0( JOBV, M, N, A, LDA, D, SVA, MV, V, LDV, EPS,
797809
$ THEN
798810
IF( ( AAQQ.LT.ROOTBIG ) .AND.
799811
$ ( AAQQ.GT.ROOTSFMIN ) ) THEN
800-
SVA( q ) = SCNRM2( M, A( 1, q ), 1)
812+
SVA( q ) = SCNRM2( M, A( 1, q ),
813+
$ 1)
801814
ELSE
802815
T = ZERO
803816
AAQQ = ONE

‎lapack-netlib/SRC/dgesvj.c‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -992,7 +992,7 @@ f"> */
992992

993993
/* Row-cyclic Jacobi SVD algorithm with column pivoting */
994994

995-
emptsw = *n * (*n - 1) / 2;
995+
emptsw = (*n / 2) * (*n - 1 + *n % 2);
996996
notrot = 0;
997997
fastr[0] = 0.;
998998

‎lapack-netlib/SRC/dgesvj.f‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -695,7 +695,7 @@ SUBROUTINE DGESVJ( JOBA, JOBU, JOBV, M, N, A, LDA, SVA, MV, V,
695695
*
696696
* Row-cyclic Jacobi SVD algorithm with column pivoting
697697
*
698-
EMPTSW = ( N*( N-1 ) ) / 2
698+
EMPTSW = ( N / 2 )*( N-1 + MOD( N, 2 ) )
699699
NOTROT = 0
700700
FASTR( 1 ) = ZERO
701701
*

‎lapack-netlib/SRC/dgsvj0.c‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -574,7 +574,7 @@ f"> */
574574

575575
/* -#- Row-cyclic Jacobi SVD algorithm with column pivoting -#- */
576576

577-
emptsw = *n * (*n - 1) / 2;
577+
emptsw = (*n / 2) * (*n - 1 + *n % 2);
578578
notrot = 0;
579579
fastr[0] = 0.;
580580

0 commit comments

Comments
 (0)