From 326bb27863ec19244bb3e86a664d74d7e94a58b9 Mon Sep 17 00:00:00 2001 From: Federico Perini Date: Sat, 15 Jun 2024 09:39:59 -0500 Subject: [PATCH 1/3] bound-check calls to `stdlib_*laset` from `*gesdd` --- src/stdlib_linalg_lapack_c.fypp | 8 ++++---- src/stdlib_linalg_lapack_d.fypp | 8 ++++---- src/stdlib_linalg_lapack_q.fypp | 8 ++++---- src/stdlib_linalg_lapack_s.fypp | 8 ++++---- src/stdlib_linalg_lapack_w.fypp | 8 ++++---- src/stdlib_linalg_lapack_z.fypp | 8 ++++---- 6 files changed, 24 insertions(+), 24 deletions(-) diff --git a/src/stdlib_linalg_lapack_c.fypp b/src/stdlib_linalg_lapack_c.fypp index 65dbcb36d..b058d1354 100644 --- a/src/stdlib_linalg_lapack_c.fypp +++ b/src/stdlib_linalg_lapack_c.fypp @@ -65732,7 +65732,7 @@ module stdlib_linalg_lapack_c call stdlib_cgeqrf( m, n, a, lda, work( itau ), work( nwork ),lwork-nwork+1, & ierr ) ! zero out below r - call stdlib_claset( 'L', n-1, n-1, czero, czero, a( 2, 1 ),lda ) + if (n>1) call stdlib_claset( 'L', n-1, n-1, czero, czero, a( 2, 1 ),lda ) ie = 1 itauq = 1 itaup = itauq + n @@ -65918,7 +65918,7 @@ module stdlib_linalg_lapack_c call stdlib_cungqr( m, m, n, u, ldu, work( itau ),work( nwork ), lwork-nwork+& 1, ierr ) ! produce r in a, zeroing out below it - call stdlib_claset( 'L', n-1, n-1, czero, czero, a( 2, 1 ),lda ) + if (n>1) call stdlib_claset( 'L', n-1, n-1, czero, czero, a( 2, 1 ),lda ) ie = 1 itauq = itau itaup = itauq + n @@ -66294,7 +66294,7 @@ module stdlib_linalg_lapack_c call stdlib_cgelqf( m, n, a, lda, work( itau ), work( nwork ),lwork-nwork+1, & ierr ) ! zero out above l - call stdlib_claset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) + if (m>1) call stdlib_claset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) ie = 1 itauq = 1 itaup = itauq + m @@ -66485,7 +66485,7 @@ module stdlib_linalg_lapack_c call stdlib_cunglq( n, n, m, vt, ldvt, work( itau ),work( nwork ), lwork-& nwork+1, ierr ) ! produce l in a, zeroing out above it - call stdlib_claset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) + if (m>1) call stdlib_claset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) ie = 1 itauq = itau itaup = itauq + m diff --git a/src/stdlib_linalg_lapack_d.fypp b/src/stdlib_linalg_lapack_d.fypp index 557311435..c9d2dc121 100644 --- a/src/stdlib_linalg_lapack_d.fypp +++ b/src/stdlib_linalg_lapack_d.fypp @@ -75410,7 +75410,7 @@ module stdlib_linalg_lapack_d call stdlib_dgeqrf( m, n, a, lda, work( itau ), work( nwork ),lwork - nwork + & 1, ierr ) ! zero out below r - call stdlib_dlaset( 'L', n-1, n-1, zero, zero, a( 2, 1 ), lda ) + if (n>1) call stdlib_dlaset( 'L', n-1, n-1, zero, zero, a( 2, 1 ), lda ) ie = 1 itauq = ie + n itaup = itauq + n @@ -75560,7 +75560,7 @@ module stdlib_linalg_lapack_d call stdlib_dorgqr( m, m, n, u, ldu, work( itau ),work( nwork ), lwork - & nwork + 1, ierr ) ! produce r in a, zeroing out other entries - call stdlib_dlaset( 'L', n-1, n-1, zero, zero, a( 2, 1 ), lda ) + if (n>1) call stdlib_dlaset( 'L', n-1, n-1, zero, zero, a( 2, 1 ), lda ) ie = itau itauq = ie + n itaup = itauq + n @@ -75726,7 +75726,7 @@ module stdlib_linalg_lapack_d call stdlib_dgelqf( m, n, a, lda, work( itau ), work( nwork ),lwork - nwork + & 1, ierr ) ! zero out above l - call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) + if (m>1) call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) ie = 1 itauq = ie + m itaup = itauq + m @@ -75880,7 +75880,7 @@ module stdlib_linalg_lapack_d call stdlib_dorglq( n, n, m, vt, ldvt, work( itau ),work( nwork ), lwork - & nwork + 1, ierr ) ! produce l in a, zeroing out other entries - call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) + if (m>1) call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) ie = itau itauq = ie + m itaup = itauq + m diff --git a/src/stdlib_linalg_lapack_q.fypp b/src/stdlib_linalg_lapack_q.fypp index 6dd252a47..e05e1a553 100644 --- a/src/stdlib_linalg_lapack_q.fypp +++ b/src/stdlib_linalg_lapack_q.fypp @@ -10442,7 +10442,7 @@ module stdlib_linalg_lapack_q call stdlib_qgeqrf( m, n, a, lda, work( itau ), work( nwork ),lwork - nwork + & 1, ierr ) ! zero out below r - call stdlib_qlaset( 'L', n-1, n-1, zero, zero, a( 2, 1 ), lda ) + if (n>1) call stdlib_qlaset( 'L', n-1, n-1, zero, zero, a( 2, 1 ), lda ) ie = 1 itauq = ie + n itaup = itauq + n @@ -10592,7 +10592,7 @@ module stdlib_linalg_lapack_q call stdlib_qorgqr( m, m, n, u, ldu, work( itau ),work( nwork ), lwork - & nwork + 1, ierr ) ! produce r in a, zeroing out other entries - call stdlib_qlaset( 'L', n-1, n-1, zero, zero, a( 2, 1 ), lda ) + if (n>1) call stdlib_qlaset( 'L', n-1, n-1, zero, zero, a( 2, 1 ), lda ) ie = itau itauq = ie + n itaup = itauq + n @@ -10758,7 +10758,7 @@ module stdlib_linalg_lapack_q call stdlib_qgelqf( m, n, a, lda, work( itau ), work( nwork ),lwork - nwork + & 1, ierr ) ! zero out above l - call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) + if (m>1) call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) ie = 1 itauq = ie + m itaup = itauq + m @@ -10912,7 +10912,7 @@ module stdlib_linalg_lapack_q call stdlib_qorglq( n, n, m, vt, ldvt, work( itau ),work( nwork ), lwork - & nwork + 1, ierr ) ! produce l in a, zeroing out other entries - call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) + if (m>1) call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) ie = itau itauq = ie + m itaup = itauq + m diff --git a/src/stdlib_linalg_lapack_s.fypp b/src/stdlib_linalg_lapack_s.fypp index 61c96b29c..73aadd4fe 100644 --- a/src/stdlib_linalg_lapack_s.fypp +++ b/src/stdlib_linalg_lapack_s.fypp @@ -72903,7 +72903,7 @@ module stdlib_linalg_lapack_s call stdlib_sgeqrf( m, n, a, lda, work( itau ), work( nwork ),lwork - nwork + & 1, ierr ) ! zero out below r - call stdlib_slaset( 'L', n-1, n-1, zero, zero, a( 2, 1 ), lda ) + if (n>1) call stdlib_slaset( 'L', n-1, n-1, zero, zero, a( 2, 1 ), lda ) ie = 1 itauq = ie + n itaup = itauq + n @@ -73053,7 +73053,7 @@ module stdlib_linalg_lapack_s call stdlib_sorgqr( m, m, n, u, ldu, work( itau ),work( nwork ), lwork - & nwork + 1, ierr ) ! produce r in a, zeroing out other entries - call stdlib_slaset( 'L', n-1, n-1, zero, zero, a( 2, 1 ), lda ) + if (n>1) call stdlib_slaset( 'L', n-1, n-1, zero, zero, a( 2, 1 ), lda ) ie = itau itauq = ie + n itaup = itauq + n @@ -73219,7 +73219,7 @@ module stdlib_linalg_lapack_s call stdlib_sgelqf( m, n, a, lda, work( itau ), work( nwork ),lwork - nwork + & 1, ierr ) ! zero out above l - call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) + if (m>1) call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) ie = 1 itauq = ie + m itaup = itauq + m @@ -73373,7 +73373,7 @@ module stdlib_linalg_lapack_s call stdlib_sorglq( n, n, m, vt, ldvt, work( itau ),work( nwork ), lwork - & nwork + 1, ierr ) ! produce l in a, zeroing out other entries - call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) + if (m>1) call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) ie = itau itauq = ie + m itaup = itauq + m diff --git a/src/stdlib_linalg_lapack_w.fypp b/src/stdlib_linalg_lapack_w.fypp index ea71d6446..7dc08c33e 100644 --- a/src/stdlib_linalg_lapack_w.fypp +++ b/src/stdlib_linalg_lapack_w.fypp @@ -10802,7 +10802,7 @@ module stdlib_linalg_lapack_w call stdlib_wgeqrf( m, n, a, lda, work( itau ), work( nwork ),lwork-nwork+1, & ierr ) ! zero out below r - call stdlib_wlaset( 'L', n-1, n-1, czero, czero, a( 2, 1 ),lda ) + if (n>1) call stdlib_wlaset( 'L', n-1, n-1, czero, czero, a( 2, 1 ),lda ) ie = 1 itauq = 1 itaup = itauq + n @@ -10988,7 +10988,7 @@ module stdlib_linalg_lapack_w call stdlib_wungqr( m, m, n, u, ldu, work( itau ),work( nwork ), lwork-nwork+& 1, ierr ) ! produce r in a, zeroing out below it - call stdlib_wlaset( 'L', n-1, n-1, czero, czero, a( 2, 1 ),lda ) + if (n>1) call stdlib_wlaset( 'L', n-1, n-1, czero, czero, a( 2, 1 ),lda ) ie = 1 itauq = itau itaup = itauq + n @@ -11364,7 +11364,7 @@ module stdlib_linalg_lapack_w call stdlib_wgelqf( m, n, a, lda, work( itau ), work( nwork ),lwork-nwork+1, & ierr ) ! zero out above l - call stdlib_wlaset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) + if (m>1) call stdlib_wlaset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) ie = 1 itauq = 1 itaup = itauq + m @@ -11555,7 +11555,7 @@ module stdlib_linalg_lapack_w call stdlib_wunglq( n, n, m, vt, ldvt, work( itau ),work( nwork ), lwork-& nwork+1, ierr ) ! produce l in a, zeroing out above it - call stdlib_wlaset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) + if (m>1) call stdlib_wlaset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) ie = 1 itauq = itau itaup = itauq + m diff --git a/src/stdlib_linalg_lapack_z.fypp b/src/stdlib_linalg_lapack_z.fypp index e28d61961..eb37665dd 100644 --- a/src/stdlib_linalg_lapack_z.fypp +++ b/src/stdlib_linalg_lapack_z.fypp @@ -66211,7 +66211,7 @@ module stdlib_linalg_lapack_z call stdlib_zgeqrf( m, n, a, lda, work( itau ), work( nwork ),lwork-nwork+1, & ierr ) ! zero out below r - call stdlib_zlaset( 'L', n-1, n-1, czero, czero, a( 2, 1 ),lda ) + if (n>1) call stdlib_zlaset( 'L', n-1, n-1, czero, czero, a( 2, 1 ),lda ) ie = 1 itauq = 1 itaup = itauq + n @@ -66397,7 +66397,7 @@ module stdlib_linalg_lapack_z call stdlib_zungqr( m, m, n, u, ldu, work( itau ),work( nwork ), lwork-nwork+& 1, ierr ) ! produce r in a, zeroing out below it - call stdlib_zlaset( 'L', n-1, n-1, czero, czero, a( 2, 1 ),lda ) + if (n>1) call stdlib_zlaset( 'L', n-1, n-1, czero, czero, a( 2, 1 ),lda ) ie = 1 itauq = itau itaup = itauq + n @@ -66773,7 +66773,7 @@ module stdlib_linalg_lapack_z call stdlib_zgelqf( m, n, a, lda, work( itau ), work( nwork ),lwork-nwork+1, & ierr ) ! zero out above l - call stdlib_zlaset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) + if (m>1) call stdlib_zlaset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) ie = 1 itauq = 1 itaup = itauq + m @@ -66964,7 +66964,7 @@ module stdlib_linalg_lapack_z call stdlib_zunglq( n, n, m, vt, ldvt, work( itau ),work( nwork ), lwork-& nwork+1, ierr ) ! produce l in a, zeroing out above it - call stdlib_zlaset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) + if (m>1) call stdlib_zlaset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) ie = 1 itauq = itau itaup = itauq + m From e832a103c90f9c8cea5575c618f2d4549b75d1f4 Mon Sep 17 00:00:00 2001 From: Federico Perini Date: Sat, 15 Jun 2024 09:47:20 -0500 Subject: [PATCH 2/3] fix test extend to all routines; add test add test extend checks to all routines where applicable Revert "fix test" This reverts commit de9a53d024aba03bbcb1741a674ccd578c8d9236. Reapply "fix test" This reverts commit ffe726d1bb9861eb15a7c46cc40d74de3cb9c0bb. Revert "Reapply "fix test"" This reverts commit 720ce1199e9bf965ad43834fd1df3af17e1ae57d. add checks everywhere applicable --- src/stdlib_linalg_lapack_c.fypp | 48 ++++++++++++++++----------------- src/stdlib_linalg_lapack_d.fypp | 20 +++++++------- src/stdlib_linalg_lapack_q.fypp | 22 +++++++-------- src/stdlib_linalg_lapack_s.fypp | 46 +++++++++++++++---------------- src/stdlib_linalg_lapack_w.fypp | 48 ++++++++++++++++----------------- src/stdlib_linalg_lapack_z.fypp | 48 ++++++++++++++++----------------- 6 files changed, 116 insertions(+), 116 deletions(-) diff --git a/src/stdlib_linalg_lapack_c.fypp b/src/stdlib_linalg_lapack_c.fypp index b058d1354..450701c90 100644 --- a/src/stdlib_linalg_lapack_c.fypp +++ b/src/stdlib_linalg_lapack_c.fypp @@ -68327,7 +68327,7 @@ module stdlib_linalg_lapack_c call stdlib_cgelqf( m, n, a, lda, work( itau ), work( iwork ),lwork-iwork+1, & ierr ) ! zero out above l - call stdlib_claset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) + if (m>1) call stdlib_claset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) ie = 1 itauq = 1 itaup = itauq + m @@ -68483,7 +68483,7 @@ module stdlib_linalg_lapack_c 1, ierr ) ! copy l to u, zeroing about above it call stdlib_clacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_claset( 'U', m-1, m-1, czero, czero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_claset( 'U', m-1, m-1, czero, czero, u( 1, 2 ),ldu ) ! generate q in a ! (cworkspace: need m*m+2*m, prefer m*m+m+m*nb) ! (rworkspace: 0) @@ -68540,7 +68540,7 @@ module stdlib_linalg_lapack_c 1, ierr ) ! copy l to u, zeroing out above it call stdlib_clacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_claset( 'U', m-1, m-1, czero, czero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_claset( 'U', m-1, m-1, czero, czero, u( 1, 2 ),ldu ) ! generate q in a ! (cworkspace: need 2*m, prefer m+m*nb) ! (rworkspace: 0) @@ -68654,7 +68654,7 @@ module stdlib_linalg_lapack_c itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_claset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) + if (m>1) call stdlib_claset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) ! bidiagonalize l in a ! (cworkspace: need 3*m, prefer 2*m+2*m*nb) ! (rworkspace: need m) @@ -68774,7 +68774,7 @@ module stdlib_linalg_lapack_c itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_claset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) + if (m>1) call stdlib_claset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) ! bidiagonalize l in a ! (cworkspace: need 3*m, prefer 2*m+2*m*nb) ! (rworkspace: need m) @@ -68882,7 +68882,7 @@ module stdlib_linalg_lapack_c lwork-iwork+1, ierr ) ! copy l to u, zeroing out above it call stdlib_clacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_claset( 'U', m-1, m-1, czero, czero,u( 1, 2 ), ldu ) + if (m>1) call stdlib_claset( 'U', m-1, m-1, czero, czero,u( 1, 2 ), ldu ) ie = 1 itauq = itau itaup = itauq + m @@ -68995,7 +68995,7 @@ module stdlib_linalg_lapack_c itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_claset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) + if (m>1) call stdlib_claset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) ! bidiagonalize l in a ! (cworkspace: need 3*m, prefer 2*m+2*m*nb) ! (rworkspace: need m) @@ -69117,7 +69117,7 @@ module stdlib_linalg_lapack_c itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_claset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) + if (m>1) call stdlib_claset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) ! bidiagonalize l in a ! (cworkspace: need 3*m, prefer 2*m+2*m*nb) ! (rworkspace: need m) @@ -69228,7 +69228,7 @@ module stdlib_linalg_lapack_c lwork-iwork+1, ierr ) ! copy l to u, zeroing out above it call stdlib_clacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_claset( 'U', m-1, m-1, czero, czero,u( 1, 2 ), ldu ) + if (m>1) call stdlib_claset( 'U', m-1, m-1, czero, czero,u( 1, 2 ), ldu ) ie = 1 itauq = itau itaup = itauq + m @@ -70098,7 +70098,7 @@ module stdlib_linalg_lapack_c v(q,p) = conjg(u(p,nr+q)) end do end do - call stdlib_claset('U',nr-1,nr-1,czero,czero,v(1,2),ldv) + if (nr>1) call stdlib_claset('U',nr-1,nr-1,czero,czero,v(1,2),ldv) call stdlib_cgesvd( 'S', 'O', nr, nr, v, ldv, s, u, ldu,v,ldv, cwork(n+nr+& 1),lcwork-n-nr,rwork, info ) call stdlib_claset('A',n-nr,nr,czero,czero,v(nr+1,1),ldv) @@ -75163,7 +75163,7 @@ module stdlib_linalg_lapack_c end do end do else - call stdlib_claset( 'U', nr-1,nr-1, czero,czero, a(1,2),lda ) + if (nr>1) call stdlib_claset( 'U', nr-1,nr-1, czero,czero, a(1,2),lda ) end if ! Second Preconditioning Using The Qr Factorization call stdlib_cgeqrf( n,nr, a,lda, cwork, cwork(n+1),lwork-n, ierr ) @@ -75188,7 +75188,7 @@ module stdlib_linalg_lapack_c end do end do else - call stdlib_claset( 'U', nr-1, nr-1, czero, czero, a(1,2), lda ) + if (nr>1) call stdlib_claset( 'U', nr-1, nr-1, czero, czero, a(1,2), lda ) end if ! .. and one-sided jacobi rotations are started on a lower ! triangular matrix (plus perturbation which is ignored in @@ -75206,7 +75206,7 @@ module stdlib_linalg_lapack_c call stdlib_ccopy( n-p+1, a(p,p), lda, v(p,p), 1 ) call stdlib_clacgv( n-p+1, v(p,p), 1 ) end do - call stdlib_claset( 'U', nr-1,nr-1, czero, czero, v(1,2), ldv ) + if (nr>1) call stdlib_claset( 'U', nr-1,nr-1, czero, czero, v(1,2), ldv ) call stdlib_cgesvj( 'L','U','N', n, nr, v, ldv, sva, nr, a, lda,cwork, lwork, & rwork, lrwork, info ) scalem = rwork(1) @@ -75214,17 +75214,17 @@ module stdlib_linalg_lapack_c else ! .. two more qr factorizations ( one qrf is not enough, two require ! accumulated product of jacobi rotations, three are perfect ) - call stdlib_claset( 'L', nr-1,nr-1, czero, czero, a(2,1), lda ) + if (nr>1) call stdlib_claset( 'L', nr-1,nr-1, czero, czero, a(2,1), lda ) call stdlib_cgelqf( nr,n, a, lda, cwork, cwork(n+1), lwork-n, ierr) call stdlib_clacpy( 'L', nr, nr, a, lda, v, ldv ) - call stdlib_claset( 'U', nr-1,nr-1, czero, czero, v(1,2), ldv ) + if (nr>1) call stdlib_claset( 'U', nr-1,nr-1, czero, czero, v(1,2), ldv ) call stdlib_cgeqrf( nr, nr, v, ldv, cwork(n+1), cwork(2*n+1),lwork-2*n, ierr ) do p = 1, nr call stdlib_ccopy( nr-p+1, v(p,p), ldv, v(p,p), 1 ) call stdlib_clacgv( nr-p+1, v(p,p), 1 ) end do - call stdlib_claset('U', nr-1, nr-1, czero, czero, v(1,2), ldv) + if (nr>1) call stdlib_claset('U', nr-1, nr-1, czero, czero, v(1,2), ldv) call stdlib_cgesvj( 'L', 'U','N', nr, nr, v,ldv, sva, nr, u,ldu, cwork(n+1), & lwork-n, rwork, lrwork, info ) scalem = rwork(1) @@ -75247,7 +75247,7 @@ module stdlib_linalg_lapack_c call stdlib_clacpy( 'A', n, n, v, ldv, u, ldu ) end if else if ( jracc .and. (.not. lsvec) .and. ( nr== n ) ) then - call stdlib_claset( 'L', n-1,n-1, czero, czero, a(2,1), lda ) + if (n>1) call stdlib_claset( 'L', n-1,n-1, czero, czero, a(2,1), lda ) call stdlib_cgesvj( 'U','N','V', n, n, a, lda, sva, n, v, ldv,cwork, lwork, rwork, & lrwork, info ) scalem = rwork(1) @@ -75261,14 +75261,14 @@ module stdlib_linalg_lapack_c call stdlib_ccopy( n-p+1, a(p,p), lda, u(p,p), 1 ) call stdlib_clacgv( n-p+1, u(p,p), 1 ) end do - call stdlib_claset( 'U', nr-1, nr-1, czero, czero, u(1,2), ldu ) + if (nr>1) call stdlib_claset( 'U', nr-1, nr-1, czero, czero, u(1,2), ldu ) call stdlib_cgeqrf( n, nr, u, ldu, cwork(n+1), cwork(2*n+1),lwork-2*n, ierr ) do p = 1, nr - 1 call stdlib_ccopy( nr-p, u(p,p+1), ldu, u(p+1,p), 1 ) call stdlib_clacgv( n-p+1, u(p,p), 1 ) end do - call stdlib_claset( 'U', nr-1, nr-1, czero, czero, u(1,2), ldu ) + if (nr>1) call stdlib_claset( 'U', nr-1, nr-1, czero, czero, u(1,2), ldu ) call stdlib_cgesvj( 'L', 'U', 'N', nr,nr, u, ldu, sva, nr, a,lda, cwork(n+1), lwork-& n, rwork, lrwork, info ) scalem = rwork(1) @@ -75327,7 +75327,7 @@ module stdlib_linalg_lapack_c end do end do else - call stdlib_claset( 'U', nr-1, nr-1, czero, czero, v(1,2), ldv ) + if (nr>1) call stdlib_claset( 'U', nr-1, nr-1, czero, czero, v(1,2), ldv ) end if ! estimate the row scaled condition number of r1 ! (if r1 is rectangular, n > nr, then the condition number @@ -75409,7 +75409,7 @@ module stdlib_linalg_lapack_c end do end do else - call stdlib_claset( 'L',nr-1,nr-1,czero,czero,v(2,1),ldv ) + if (nr>1) call stdlib_claset( 'L',nr-1,nr-1,czero,czero,v(2,1),ldv ) end if ! now, compute r2 = l3 * q3, the lq factorization. call stdlib_cgelqf( nr, nr, v, ldv, cwork(2*n+n*nr+1),cwork(2*n+n*nr+nr+1), & @@ -75443,7 +75443,7 @@ module stdlib_linalg_lapack_c end do end do else - call stdlib_claset( 'U', nr-1,nr-1, czero,czero, v(1,2), ldv ) + if (nr>1) call stdlib_claset( 'U', nr-1,nr-1, czero,czero, v(1,2), ldv ) end if ! second preconditioning finished; continue with jacobi svd ! the input matrix is lower trinagular. @@ -75662,7 +75662,7 @@ module stdlib_linalg_lapack_c end do end do else - call stdlib_claset( 'U', nr-1, nr-1, czero, czero, v(1,2), ldv ) + if (nr>1) call stdlib_claset( 'U', nr-1, nr-1, czero, czero, v(1,2), ldv ) end if call stdlib_cgeqrf( n, nr, v, ldv, cwork(n+1), cwork(2*n+1),lwork-2*n, ierr ) @@ -75681,7 +75681,7 @@ module stdlib_linalg_lapack_c end do end do else - call stdlib_claset('U', nr-1, nr-1, czero, czero, u(1,2), ldu ) + if (nr>1) call stdlib_claset('U', nr-1, nr-1, czero, czero, u(1,2), ldu ) end if call stdlib_cgesvj( 'L', 'U', 'V', nr, nr, u, ldu, sva,n, v, ldv, cwork(2*n+n*nr+1),& lwork-2*n-n*nr,rwork, lrwork, info ) diff --git a/src/stdlib_linalg_lapack_d.fypp b/src/stdlib_linalg_lapack_d.fypp index c9d2dc121..f5e6b1985 100644 --- a/src/stdlib_linalg_lapack_d.fypp +++ b/src/stdlib_linalg_lapack_d.fypp @@ -77363,7 +77363,7 @@ module stdlib_linalg_lapack_d call stdlib_dgelqf( m, n, a, lda, work( itau ), work( iwork ),lwork-iwork+1, & ierr ) ! zero out above l - call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) + if (m>1) call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) ie = 1 itauq = ie + m itaup = itauq + m @@ -77506,7 +77506,7 @@ module stdlib_linalg_lapack_d 1, ierr ) ! copy l to u, zeroing about above it call stdlib_dlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_dlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_dlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) ! generate q in a ! (workspace: need m*m + 2*m, prefer m*m + m + m*nb) call stdlib_dorglq( m, n, m, a, lda, work( itau ),work( iwork ), lwork-& @@ -77556,7 +77556,7 @@ module stdlib_linalg_lapack_d 1, ierr ) ! copy l to u, zeroing out above it call stdlib_dlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_dlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_dlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) ! generate q in a ! (workspace: need 2*m, prefer m + m*nb) call stdlib_dorglq( m, n, m, a, lda, work( itau ),work( iwork ), lwork-& @@ -77657,7 +77657,7 @@ module stdlib_linalg_lapack_d itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) + if (m>1) call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) ! bidiagonalize l in a ! (workspace: need 4*m, prefer 3*m + 2*m*nb) call stdlib_dgebrd( m, m, a, lda, s, work( ie ),work( itauq ), work( & @@ -77764,7 +77764,7 @@ module stdlib_linalg_lapack_d itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) + if (m>1) call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) ! bidiagonalize l in a ! (workspace: need 4*m, prefer 3*m + 2*m*nb) call stdlib_dgebrd( m, m, a, lda, s, work( ie ),work( itauq ), work( & @@ -77859,7 +77859,7 @@ module stdlib_linalg_lapack_d lwork-iwork+1, ierr ) ! copy l to u, zeroing out above it call stdlib_dlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_dlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_dlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) ie = itau itauq = ie + m itaup = itauq + m @@ -77960,7 +77960,7 @@ module stdlib_linalg_lapack_d itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) + if (m>1) call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) ! bidiagonalize l in a ! (workspace: need 4*m, prefer 3*m + 2*m*nb) call stdlib_dgebrd( m, m, a, lda, s, work( ie ),work( itauq ), work( & @@ -78070,7 +78070,7 @@ module stdlib_linalg_lapack_d itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) + if (m>1) call stdlib_dlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) ! bidiagonalize l in a ! (workspace: need 4*m, prefer 3*m + 2*m*nb) call stdlib_dgebrd( m, m, a, lda, s, work( ie ),work( itauq ), work( & @@ -78168,7 +78168,7 @@ module stdlib_linalg_lapack_d lwork-iwork+1, ierr ) ! copy l to u, zeroing out above it call stdlib_dlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_dlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_dlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) ie = itau itauq = ie + m itaup = itauq + m @@ -79029,7 +79029,7 @@ module stdlib_linalg_lapack_d v(q,p) = u(p,nr+q) end do end do - call stdlib_dlaset('U',nr-1,nr-1,zero,zero,v(1,2),ldv) + if (nr>1) call stdlib_dlaset('U',nr-1,nr-1,zero,zero,v(1,2),ldv) call stdlib_dgesvd( 'S', 'O', nr, nr, v, ldv, s, u, ldu,v,ldv, work(n+nr+1)& ,lwork-n-nr, info ) call stdlib_dlaset('A',n-nr,nr,zero,zero,v(nr+1,1),ldv) diff --git a/src/stdlib_linalg_lapack_q.fypp b/src/stdlib_linalg_lapack_q.fypp index e05e1a553..fb9d0c137 100644 --- a/src/stdlib_linalg_lapack_q.fypp +++ b/src/stdlib_linalg_lapack_q.fypp @@ -12443,7 +12443,7 @@ module stdlib_linalg_lapack_q call stdlib_qgelqf( m, n, a, lda, work( itau ), work( iwork ),lwork-iwork+1, & ierr ) ! zero out above l - call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) + if (m>1) call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) ie = 1 itauq = ie + m itaup = itauq + m @@ -12586,7 +12586,7 @@ module stdlib_linalg_lapack_q 1, ierr ) ! copy l to u, zeroing about above it call stdlib_qlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_qlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_qlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) ! generate q in a ! (workspace: need m*m + 2*m, prefer m*m + m + m*nb) call stdlib_qorglq( m, n, m, a, lda, work( itau ),work( iwork ), lwork-& @@ -12636,7 +12636,7 @@ module stdlib_linalg_lapack_q 1, ierr ) ! copy l to u, zeroing out above it call stdlib_qlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_qlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_qlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) ! generate q in a ! (workspace: need 2*m, prefer m + m*nb) call stdlib_qorglq( m, n, m, a, lda, work( itau ),work( iwork ), lwork-& @@ -12737,7 +12737,7 @@ module stdlib_linalg_lapack_q itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) + if (m>1) call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) ! bidiagonalize l in a ! (workspace: need 4*m, prefer 3*m + 2*m*nb) call stdlib_qgebrd( m, m, a, lda, s, work( ie ),work( itauq ), work( & @@ -12844,7 +12844,7 @@ module stdlib_linalg_lapack_q itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) + if (m>1) call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) ! bidiagonalize l in a ! (workspace: need 4*m, prefer 3*m + 2*m*nb) call stdlib_qgebrd( m, m, a, lda, s, work( ie ),work( itauq ), work( & @@ -12939,7 +12939,7 @@ module stdlib_linalg_lapack_q lwork-iwork+1, ierr ) ! copy l to u, zeroing out above it call stdlib_qlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_qlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_qlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) ie = itau itauq = ie + m itaup = itauq + m @@ -13040,7 +13040,7 @@ module stdlib_linalg_lapack_q itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) + if (m>1) call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) ! bidiagonalize l in a ! (workspace: need 4*m, prefer 3*m + 2*m*nb) call stdlib_qgebrd( m, m, a, lda, s, work( ie ),work( itauq ), work( & @@ -13150,7 +13150,7 @@ module stdlib_linalg_lapack_q itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) + if (m>1) call stdlib_qlaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) ! bidiagonalize l in a ! (workspace: need 4*m, prefer 3*m + 2*m*nb) call stdlib_qgebrd( m, m, a, lda, s, work( ie ),work( itauq ), work( & @@ -13248,7 +13248,7 @@ module stdlib_linalg_lapack_q lwork-iwork+1, ierr ) ! copy l to u, zeroing out above it call stdlib_qlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_qlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_qlaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) ie = itau itauq = ie + m itaup = itauq + m @@ -14109,7 +14109,7 @@ module stdlib_linalg_lapack_q v(q,p) = u(p,nr+q) end do end do - call stdlib_qlaset('U',nr-1,nr-1,zero,zero,v(1,2),ldv) + if (nr>1) call stdlib_qlaset('U',nr-1,nr-1,zero,zero,v(1,2),ldv) call stdlib_qgesvd( 'S', 'O', nr, nr, v, ldv, s, u, ldu,v,ldv, work(n+nr+1)& ,lwork-n-nr, info ) call stdlib_qlaset('A',n-nr,nr,zero,zero,v(nr+1,1),ldv) @@ -77581,7 +77581,7 @@ module stdlib_linalg_lapack_q go to 70 end if ! set lower triangle of b-part to zero - call stdlib_qlaset( 'LOWER', m-1, m-1, zero, zero, t(2,1), ldst ) + if (m>1) call stdlib_qlaset( 'LOWER', m-1, m-1, zero, zero, t(2,1), ldst ) if( wands ) then ! strong stability test: ! f-norm((a-ql**h*s*qr)) <= o(eps*f-norm((a))) diff --git a/src/stdlib_linalg_lapack_s.fypp b/src/stdlib_linalg_lapack_s.fypp index 73aadd4fe..bd8353c6c 100644 --- a/src/stdlib_linalg_lapack_s.fypp +++ b/src/stdlib_linalg_lapack_s.fypp @@ -56700,7 +56700,7 @@ module stdlib_linalg_lapack_s go to 70 end if ! set lower triangle of b-part to zero - call stdlib_slaset( 'LOWER', m-1, m-1, zero, zero, t(2,1), ldst ) + if (m>1) call stdlib_slaset( 'LOWER', m-1, m-1, zero, zero, t(2,1), ldst ) if( wands ) then ! strong stability test: ! f-norm((a-ql**h*s*qr)) <= o(eps*f-norm((a))) @@ -71232,7 +71232,7 @@ module stdlib_linalg_lapack_s end do end do else - call stdlib_slaset( 'U', nr-1,nr-1, zero,zero, a(1,2),lda ) + if (nr>1) call stdlib_slaset( 'U', nr-1,nr-1, zero,zero, a(1,2),lda ) end if ! Second Preconditioning Using The Qr Factorization call stdlib_sgeqrf( n,nr, a,lda, work, work(n+1),lwork-n, ierr ) @@ -71255,7 +71255,7 @@ module stdlib_linalg_lapack_s end do end do else - call stdlib_slaset( 'U', nr-1, nr-1, zero, zero, a(1,2), lda ) + if (nr>1) call stdlib_slaset( 'U', nr-1, nr-1, zero, zero, a(1,2), lda ) end if ! .. and one-sided jacobi rotations are started on a lower ! triangular matrix (plus perturbation which is ignored in @@ -71271,7 +71271,7 @@ module stdlib_linalg_lapack_s do p = 1, nr call stdlib_scopy( n-p+1, a(p,p), lda, v(p,p), 1 ) end do - call stdlib_slaset( 'UPPER', nr-1, nr-1, zero, zero, v(1,2), ldv ) + if (nr>1) call stdlib_slaset( 'UPPER', nr-1, nr-1, zero, zero, v(1,2), ldv ) call stdlib_sgesvj( 'L','U','N', n, nr, v,ldv, sva, nr, a,lda,work, lwork, info ) scalem = work(1) @@ -71279,16 +71279,16 @@ module stdlib_linalg_lapack_s else ! .. two more qr factorizations ( one qrf is not enough, two require ! accumulated product of jacobi rotations, three are perfect ) - call stdlib_slaset( 'LOWER', nr-1, nr-1, zero, zero, a(2,1), lda ) + if (nr>1) call stdlib_slaset( 'LOWER', nr-1, nr-1, zero, zero, a(2,1), lda ) call stdlib_sgelqf( nr, n, a, lda, work, work(n+1), lwork-n, ierr) call stdlib_slacpy( 'LOWER', nr, nr, a, lda, v, ldv ) - call stdlib_slaset( 'UPPER', nr-1, nr-1, zero, zero, v(1,2), ldv ) + if (nr>1) call stdlib_slaset( 'UPPER', nr-1, nr-1, zero, zero, v(1,2), ldv ) call stdlib_sgeqrf( nr, nr, v, ldv, work(n+1), work(2*n+1),lwork-2*n, ierr ) do p = 1, nr call stdlib_scopy( nr-p+1, v(p,p), ldv, v(p,p), 1 ) end do - call stdlib_slaset( 'UPPER', nr-1, nr-1, zero, zero, v(1,2), ldv ) + if (nr>1) call stdlib_slaset( 'UPPER', nr-1, nr-1, zero, zero, v(1,2), ldv ) call stdlib_sgesvj( 'LOWER', 'U','N', nr, nr, v,ldv, sva, nr, u,ldu, work(n+1), & lwork-n, info ) scalem = work(n+1) @@ -71315,12 +71315,12 @@ module stdlib_linalg_lapack_s do p = 1, nr call stdlib_scopy( n-p+1, a(p,p), lda, u(p,p), 1 ) end do - call stdlib_slaset( 'UPPER', nr-1, nr-1, zero, zero, u(1,2), ldu ) + if (nr>1) call stdlib_slaset( 'UPPER', nr-1, nr-1, zero, zero, u(1,2), ldu ) call stdlib_sgeqrf( n, nr, u, ldu, work(n+1), work(2*n+1),lwork-2*n, ierr ) do p = 1, nr - 1 call stdlib_scopy( nr-p, u(p,p+1), ldu, u(p+1,p), 1 ) end do - call stdlib_slaset( 'UPPER', nr-1, nr-1, zero, zero, u(1,2), ldu ) + if (nr>1) call stdlib_slaset( 'UPPER', nr-1, nr-1, zero, zero, u(1,2), ldu ) call stdlib_sgesvj( 'LOWER', 'U', 'N', nr,nr, u, ldu, sva, nr, a,lda, work(n+1), & lwork-n, info ) scalem = work(n+1) @@ -71377,7 +71377,7 @@ module stdlib_linalg_lapack_s end do end do else - call stdlib_slaset( 'U', nr-1, nr-1, zero, zero, v(1,2), ldv ) + if (nr>1) call stdlib_slaset( 'U', nr-1, nr-1, zero, zero, v(1,2), ldv ) end if ! estimate the row scaled condition number of r1 ! (if r1 is rectangular, n > nr, then the condition number @@ -71453,7 +71453,7 @@ module stdlib_linalg_lapack_s end do end do else - call stdlib_slaset( 'L',nr-1,nr-1,zero,zero,v(2,1),ldv ) + if (nr>1) call stdlib_slaset( 'L',nr-1,nr-1,zero,zero,v(2,1),ldv ) end if ! now, compute r2 = l3 * q3, the lq factorization. call stdlib_sgelqf( nr, nr, v, ldv, work(2*n+n*nr+1),work(2*n+n*nr+nr+1), & @@ -71487,7 +71487,7 @@ module stdlib_linalg_lapack_s end do end do else - call stdlib_slaset( 'U', nr-1,nr-1, zero,zero, v(1,2), ldv ) + if (nr>1) call stdlib_slaset( 'U', nr-1,nr-1, zero,zero, v(1,2), ldv ) end if ! second preconditioning finished; continue with jacobi svd ! the input matrix is lower trinagular. @@ -71704,7 +71704,7 @@ module stdlib_linalg_lapack_s end do end do else - call stdlib_slaset( 'U', nr-1, nr-1, zero, zero, v(1,2), ldv ) + if (nr>1) call stdlib_slaset( 'U', nr-1, nr-1, zero, zero, v(1,2), ldv ) end if call stdlib_sgeqrf( n, nr, v, ldv, work(n+1), work(2*n+1),lwork-2*n, ierr ) call stdlib_slacpy( 'L', n, nr, v, ldv, work(2*n+1), n ) @@ -71720,7 +71720,7 @@ module stdlib_linalg_lapack_s end do end do else - call stdlib_slaset('U', nr-1, nr-1, zero, zero, u(1,2), ldu ) + if (nr>1) call stdlib_slaset('U', nr-1, nr-1, zero, zero, u(1,2), ldu ) end if call stdlib_sgesvj( 'L', 'U', 'V', nr, nr, u, ldu, sva,n, v, ldv, work(2*n+n*nr+1), & lwork-2*n-n*nr, info ) @@ -74857,7 +74857,7 @@ module stdlib_linalg_lapack_s call stdlib_sgelqf( m, n, a, lda, work( itau ), work( iwork ),lwork-iwork+1, & ierr ) ! zero out above l - call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) + if (m>1) call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ), lda ) ie = 1 itauq = ie + m itaup = itauq + m @@ -75000,7 +75000,7 @@ module stdlib_linalg_lapack_s 1, ierr ) ! copy l to u, zeroing about above it call stdlib_slacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_slaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_slaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) ! generate q in a ! (workspace: need m*m+2*m, prefer m*m+m+m*nb) call stdlib_sorglq( m, n, m, a, lda, work( itau ),work( iwork ), lwork-& @@ -75050,7 +75050,7 @@ module stdlib_linalg_lapack_s 1, ierr ) ! copy l to u, zeroing out above it call stdlib_slacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_slaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_slaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) ! generate q in a ! (workspace: need 2*m, prefer m+m*nb) call stdlib_sorglq( m, n, m, a, lda, work( itau ),work( iwork ), lwork-& @@ -75151,7 +75151,7 @@ module stdlib_linalg_lapack_s itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) + if (m>1) call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) ! bidiagonalize l in a ! (workspace: need 4*m, prefer 3*m+2*m*nb) call stdlib_sgebrd( m, m, a, lda, s, work( ie ),work( itauq ), work( & @@ -75258,7 +75258,7 @@ module stdlib_linalg_lapack_s itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) + if (m>1) call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) ! bidiagonalize l in a ! (workspace: need 4*m, prefer 3*m+2*m*nb) call stdlib_sgebrd( m, m, a, lda, s, work( ie ),work( itauq ), work( & @@ -75353,7 +75353,7 @@ module stdlib_linalg_lapack_s lwork-iwork+1, ierr ) ! copy l to u, zeroing out above it call stdlib_slacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_slaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_slaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) ie = itau itauq = ie + m itaup = itauq + m @@ -75454,7 +75454,7 @@ module stdlib_linalg_lapack_s itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) + if (m>1) call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) ! bidiagonalize l in a ! (workspace: need 4*m, prefer 3*m+2*m*nb) call stdlib_sgebrd( m, m, a, lda, s, work( ie ),work( itauq ), work( & @@ -75564,7 +75564,7 @@ module stdlib_linalg_lapack_s itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) + if (m>1) call stdlib_slaset( 'U', m-1, m-1, zero, zero, a( 1, 2 ),lda ) ! bidiagonalize l in a ! (workspace: need 4*m, prefer 3*m+2*m*nb) call stdlib_sgebrd( m, m, a, lda, s, work( ie ),work( itauq ), work( & @@ -75662,7 +75662,7 @@ module stdlib_linalg_lapack_s lwork-iwork+1, ierr ) ! copy l to u, zeroing out above it call stdlib_slacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_slaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_slaset( 'U', m-1, m-1, zero, zero, u( 1, 2 ),ldu ) ie = itau itauq = ie + m itaup = itauq + m diff --git a/src/stdlib_linalg_lapack_w.fypp b/src/stdlib_linalg_lapack_w.fypp index 7dc08c33e..8c2323529 100644 --- a/src/stdlib_linalg_lapack_w.fypp +++ b/src/stdlib_linalg_lapack_w.fypp @@ -6429,7 +6429,7 @@ module stdlib_linalg_lapack_w end do end do else - call stdlib_wlaset( 'U', nr-1,nr-1, czero,czero, a(1,2),lda ) + if (nr>1) call stdlib_wlaset( 'U', nr-1,nr-1, czero,czero, a(1,2),lda ) end if ! Second Preconditioning Using The Qr Factorization call stdlib_wgeqrf( n,nr, a,lda, cwork, cwork(n+1),lwork-n, ierr ) @@ -6454,7 +6454,7 @@ module stdlib_linalg_lapack_w end do end do else - call stdlib_wlaset( 'U', nr-1, nr-1, czero, czero, a(1,2), lda ) + if (nr>1) call stdlib_wlaset( 'U', nr-1, nr-1, czero, czero, a(1,2), lda ) end if ! .. and one-sided jacobi rotations are started on a lower ! triangular matrix (plus perturbation which is ignored in @@ -6472,7 +6472,7 @@ module stdlib_linalg_lapack_w call stdlib_wcopy( n-p+1, a(p,p), lda, v(p,p), 1 ) call stdlib_wlacgv( n-p+1, v(p,p), 1 ) end do - call stdlib_wlaset( 'U', nr-1,nr-1, czero, czero, v(1,2), ldv ) + if (nr>1) call stdlib_wlaset( 'U', nr-1,nr-1, czero, czero, v(1,2), ldv ) call stdlib_wgesvj( 'L','U','N', n, nr, v, ldv, sva, nr, a, lda,cwork, lwork, & rwork, lrwork, info ) scalem = rwork(1) @@ -6480,17 +6480,17 @@ module stdlib_linalg_lapack_w else ! .. two more qr factorizations ( one qrf is not enough, two require ! accumulated product of jacobi rotations, three are perfect ) - call stdlib_wlaset( 'L', nr-1,nr-1, czero, czero, a(2,1), lda ) + if (nr>1) call stdlib_wlaset( 'L', nr-1,nr-1, czero, czero, a(2,1), lda ) call stdlib_wgelqf( nr,n, a, lda, cwork, cwork(n+1), lwork-n, ierr) call stdlib_wlacpy( 'L', nr, nr, a, lda, v, ldv ) - call stdlib_wlaset( 'U', nr-1,nr-1, czero, czero, v(1,2), ldv ) + if (nr>1) call stdlib_wlaset( 'U', nr-1,nr-1, czero, czero, v(1,2), ldv ) call stdlib_wgeqrf( nr, nr, v, ldv, cwork(n+1), cwork(2*n+1),lwork-2*n, ierr ) do p = 1, nr call stdlib_wcopy( nr-p+1, v(p,p), ldv, v(p,p), 1 ) call stdlib_wlacgv( nr-p+1, v(p,p), 1 ) end do - call stdlib_wlaset('U', nr-1, nr-1, czero, czero, v(1,2), ldv) + if (nr>1) call stdlib_wlaset('U', nr-1, nr-1, czero, czero, v(1,2), ldv) call stdlib_wgesvj( 'L', 'U','N', nr, nr, v,ldv, sva, nr, u,ldu, cwork(n+1), & lwork-n, rwork, lrwork, info ) scalem = rwork(1) @@ -6513,7 +6513,7 @@ module stdlib_linalg_lapack_w call stdlib_wlacpy( 'A', n, n, v, ldv, u, ldu ) end if else if ( jracc .and. (.not. lsvec) .and. ( nr== n ) ) then - call stdlib_wlaset( 'L', n-1,n-1, czero, czero, a(2,1), lda ) + if (n>1) call stdlib_wlaset( 'L', n-1,n-1, czero, czero, a(2,1), lda ) call stdlib_wgesvj( 'U','N','V', n, n, a, lda, sva, n, v, ldv,cwork, lwork, rwork, & lrwork, info ) scalem = rwork(1) @@ -6527,14 +6527,14 @@ module stdlib_linalg_lapack_w call stdlib_wcopy( n-p+1, a(p,p), lda, u(p,p), 1 ) call stdlib_wlacgv( n-p+1, u(p,p), 1 ) end do - call stdlib_wlaset( 'U', nr-1, nr-1, czero, czero, u(1,2), ldu ) + if (nr>1) call stdlib_wlaset( 'U', nr-1, nr-1, czero, czero, u(1,2), ldu ) call stdlib_wgeqrf( n, nr, u, ldu, cwork(n+1), cwork(2*n+1),lwork-2*n, ierr ) do p = 1, nr - 1 call stdlib_wcopy( nr-p, u(p,p+1), ldu, u(p+1,p), 1 ) call stdlib_wlacgv( n-p+1, u(p,p), 1 ) end do - call stdlib_wlaset( 'U', nr-1, nr-1, czero, czero, u(1,2), ldu ) + if (nr>1) call stdlib_wlaset( 'U', nr-1, nr-1, czero, czero, u(1,2), ldu ) call stdlib_wgesvj( 'L', 'U', 'N', nr,nr, u, ldu, sva, nr, a,lda, cwork(n+1), lwork-& n, rwork, lrwork, info ) scalem = rwork(1) @@ -6593,7 +6593,7 @@ module stdlib_linalg_lapack_w end do end do else - call stdlib_wlaset( 'U', nr-1, nr-1, czero, czero, v(1,2), ldv ) + if (nr>1) call stdlib_wlaset( 'U', nr-1, nr-1, czero, czero, v(1,2), ldv ) end if ! estimate the row scaled condition number of r1 ! (if r1 is rectangular, n > nr, then the condition number @@ -6675,7 +6675,7 @@ module stdlib_linalg_lapack_w end do end do else - call stdlib_wlaset( 'L',nr-1,nr-1,czero,czero,v(2,1),ldv ) + if (nr>1) call stdlib_wlaset( 'L',nr-1,nr-1,czero,czero,v(2,1),ldv ) end if ! now, compute r2 = l3 * q3, the lq factorization. call stdlib_wgelqf( nr, nr, v, ldv, cwork(2*n+n*nr+1),cwork(2*n+n*nr+nr+1), & @@ -6709,7 +6709,7 @@ module stdlib_linalg_lapack_w end do end do else - call stdlib_wlaset( 'U', nr-1,nr-1, czero,czero, v(1,2), ldv ) + if (nr>1) call stdlib_wlaset( 'U', nr-1,nr-1, czero,czero, v(1,2), ldv ) end if ! second preconditioning finished; continue with jacobi svd ! the input matrix is lower trinagular. @@ -6928,7 +6928,7 @@ module stdlib_linalg_lapack_w end do end do else - call stdlib_wlaset( 'U', nr-1, nr-1, czero, czero, v(1,2), ldv ) + if (nr>1) call stdlib_wlaset( 'U', nr-1, nr-1, czero, czero, v(1,2), ldv ) end if call stdlib_wgeqrf( n, nr, v, ldv, cwork(n+1), cwork(2*n+1),lwork-2*n, ierr ) @@ -6947,7 +6947,7 @@ module stdlib_linalg_lapack_w end do end do else - call stdlib_wlaset('U', nr-1, nr-1, czero, czero, u(1,2), ldu ) + if (nr>1) call stdlib_wlaset('U', nr-1, nr-1, czero, czero, u(1,2), ldu ) end if call stdlib_wgesvj( 'L', 'U', 'V', nr, nr, u, ldu, sva,n, v, ldv, cwork(2*n+n*nr+1),& lwork-2*n-n*nr,rwork, lrwork, info ) @@ -13396,7 +13396,7 @@ module stdlib_linalg_lapack_w call stdlib_wgelqf( m, n, a, lda, work( itau ), work( iwork ),lwork-iwork+1, & ierr ) ! zero out above l - call stdlib_wlaset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) + if (m>1) call stdlib_wlaset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) ie = 1 itauq = 1 itaup = itauq + m @@ -13552,7 +13552,7 @@ module stdlib_linalg_lapack_w 1, ierr ) ! copy l to u, zeroing about above it call stdlib_wlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_wlaset( 'U', m-1, m-1, czero, czero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_wlaset( 'U', m-1, m-1, czero, czero, u( 1, 2 ),ldu ) ! generate q in a ! (cworkspace: need m*m+2*m, prefer m*m+m+m*nb) ! (rworkspace: 0) @@ -13609,7 +13609,7 @@ module stdlib_linalg_lapack_w 1, ierr ) ! copy l to u, zeroing out above it call stdlib_wlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_wlaset( 'U', m-1, m-1, czero, czero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_wlaset( 'U', m-1, m-1, czero, czero, u( 1, 2 ),ldu ) ! generate q in a ! (cworkspace: need 2*m, prefer m+m*nb) ! (rworkspace: 0) @@ -13723,7 +13723,7 @@ module stdlib_linalg_lapack_w itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_wlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) + if (m>1) call stdlib_wlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) ! bidiagonalize l in a ! (cworkspace: need 3*m, prefer 2*m+2*m*nb) ! (rworkspace: need m) @@ -13843,7 +13843,7 @@ module stdlib_linalg_lapack_w itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_wlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) + if (m>1) call stdlib_wlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) ! bidiagonalize l in a ! (cworkspace: need 3*m, prefer 2*m+2*m*nb) ! (rworkspace: need m) @@ -13951,7 +13951,7 @@ module stdlib_linalg_lapack_w lwork-iwork+1, ierr ) ! copy l to u, zeroing out above it call stdlib_wlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_wlaset( 'U', m-1, m-1, czero, czero,u( 1, 2 ), ldu ) + if (m>1) call stdlib_wlaset( 'U', m-1, m-1, czero, czero,u( 1, 2 ), ldu ) ie = 1 itauq = itau itaup = itauq + m @@ -14064,7 +14064,7 @@ module stdlib_linalg_lapack_w itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_wlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) + if (m>1) call stdlib_wlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) ! bidiagonalize l in a ! (cworkspace: need 3*m, prefer 2*m+2*m*nb) ! (rworkspace: need m) @@ -14186,7 +14186,7 @@ module stdlib_linalg_lapack_w itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_wlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) + if (m>1) call stdlib_wlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) ! bidiagonalize l in a ! (cworkspace: need 3*m, prefer 2*m+2*m*nb) ! (rworkspace: need m) @@ -14297,7 +14297,7 @@ module stdlib_linalg_lapack_w lwork-iwork+1, ierr ) ! copy l to u, zeroing out above it call stdlib_wlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_wlaset( 'U', m-1, m-1, czero, czero,u( 1, 2 ), ldu ) + if (m>1) call stdlib_wlaset( 'U', m-1, m-1, czero, czero,u( 1, 2 ), ldu ) ie = 1 itauq = itau itaup = itauq + m @@ -15166,7 +15166,7 @@ module stdlib_linalg_lapack_w v(q,p) = conjg(u(p,nr+q)) end do end do - call stdlib_wlaset('U',nr-1,nr-1,czero,czero,v(1,2),ldv) + if (nr>1) call stdlib_wlaset('U',nr-1,nr-1,czero,czero,v(1,2),ldv) call stdlib_wgesvd( 'S', 'O', nr, nr, v, ldv, s, u, ldu,v,ldv, cwork(n+nr+& 1),lcwork-n-nr,rwork, info ) call stdlib_wlaset('A',n-nr,nr,czero,czero,v(nr+1,1),ldv) diff --git a/src/stdlib_linalg_lapack_z.fypp b/src/stdlib_linalg_lapack_z.fypp index eb37665dd..c7e5efecf 100644 --- a/src/stdlib_linalg_lapack_z.fypp +++ b/src/stdlib_linalg_lapack_z.fypp @@ -68805,7 +68805,7 @@ module stdlib_linalg_lapack_z call stdlib_zgelqf( m, n, a, lda, work( itau ), work( iwork ),lwork-iwork+1, & ierr ) ! zero out above l - call stdlib_zlaset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) + if (m>1) call stdlib_zlaset( 'U', m-1, m-1, czero, czero, a( 1, 2 ),lda ) ie = 1 itauq = 1 itaup = itauq + m @@ -68961,7 +68961,7 @@ module stdlib_linalg_lapack_z 1, ierr ) ! copy l to u, zeroing about above it call stdlib_zlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_zlaset( 'U', m-1, m-1, czero, czero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_zlaset( 'U', m-1, m-1, czero, czero, u( 1, 2 ),ldu ) ! generate q in a ! (cworkspace: need m*m+2*m, prefer m*m+m+m*nb) ! (rworkspace: 0) @@ -69018,7 +69018,7 @@ module stdlib_linalg_lapack_z 1, ierr ) ! copy l to u, zeroing out above it call stdlib_zlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_zlaset( 'U', m-1, m-1, czero, czero, u( 1, 2 ),ldu ) + if (m>1) call stdlib_zlaset( 'U', m-1, m-1, czero, czero, u( 1, 2 ),ldu ) ! generate q in a ! (cworkspace: need 2*m, prefer m+m*nb) ! (rworkspace: 0) @@ -69132,7 +69132,7 @@ module stdlib_linalg_lapack_z itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_zlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) + if (m>1) call stdlib_zlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) ! bidiagonalize l in a ! (cworkspace: need 3*m, prefer 2*m+2*m*nb) ! (rworkspace: need m) @@ -69252,7 +69252,7 @@ module stdlib_linalg_lapack_z itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_zlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) + if (m>1) call stdlib_zlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) ! bidiagonalize l in a ! (cworkspace: need 3*m, prefer 2*m+2*m*nb) ! (rworkspace: need m) @@ -69360,7 +69360,7 @@ module stdlib_linalg_lapack_z lwork-iwork+1, ierr ) ! copy l to u, zeroing out above it call stdlib_zlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_zlaset( 'U', m-1, m-1, czero, czero,u( 1, 2 ), ldu ) + if (m>1) call stdlib_zlaset( 'U', m-1, m-1, czero, czero,u( 1, 2 ), ldu ) ie = 1 itauq = itau itaup = itauq + m @@ -69473,7 +69473,7 @@ module stdlib_linalg_lapack_z itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_zlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) + if (m>1) call stdlib_zlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) ! bidiagonalize l in a ! (cworkspace: need 3*m, prefer 2*m+2*m*nb) ! (rworkspace: need m) @@ -69595,7 +69595,7 @@ module stdlib_linalg_lapack_z itaup = itauq + m iwork = itaup + m ! zero out above l in a - call stdlib_zlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) + if (m>1) call stdlib_zlaset( 'U', m-1, m-1, czero, czero,a( 1, 2 ), lda ) ! bidiagonalize l in a ! (cworkspace: need 3*m, prefer 2*m+2*m*nb) ! (rworkspace: need m) @@ -69706,7 +69706,7 @@ module stdlib_linalg_lapack_z lwork-iwork+1, ierr ) ! copy l to u, zeroing out above it call stdlib_zlacpy( 'L', m, m, a, lda, u, ldu ) - call stdlib_zlaset( 'U', m-1, m-1, czero, czero,u( 1, 2 ), ldu ) + if (m>1) call stdlib_zlaset( 'U', m-1, m-1, czero, czero,u( 1, 2 ), ldu ) ie = 1 itauq = itau itaup = itauq + m @@ -70575,7 +70575,7 @@ module stdlib_linalg_lapack_z v(q,p) = conjg(u(p,nr+q)) end do end do - call stdlib_zlaset('U',nr-1,nr-1,czero,czero,v(1,2),ldv) + if (nr>1) call stdlib_zlaset('U',nr-1,nr-1,czero,czero,v(1,2),ldv) call stdlib_zgesvd( 'S', 'O', nr, nr, v, ldv, s, u, ldu,v,ldv, cwork(n+nr+& 1),lcwork-n-nr,rwork, info ) call stdlib_zlaset('A',n-nr,nr,czero,czero,v(nr+1,1),ldv) @@ -75638,7 +75638,7 @@ module stdlib_linalg_lapack_z end do end do else - call stdlib_zlaset( 'U', nr-1,nr-1, czero,czero, a(1,2),lda ) + if (nr>1) call stdlib_zlaset( 'U', nr-1,nr-1, czero,czero, a(1,2),lda ) end if ! Second Preconditioning Using The Qr Factorization call stdlib_zgeqrf( n,nr, a,lda, cwork, cwork(n+1),lwork-n, ierr ) @@ -75663,7 +75663,7 @@ module stdlib_linalg_lapack_z end do end do else - call stdlib_zlaset( 'U', nr-1, nr-1, czero, czero, a(1,2), lda ) + if (nr>1) call stdlib_zlaset( 'U', nr-1, nr-1, czero, czero, a(1,2), lda ) end if ! .. and one-sided jacobi rotations are started on a lower ! triangular matrix (plus perturbation which is ignored in @@ -75681,7 +75681,7 @@ module stdlib_linalg_lapack_z call stdlib_zcopy( n-p+1, a(p,p), lda, v(p,p), 1 ) call stdlib_zlacgv( n-p+1, v(p,p), 1 ) end do - call stdlib_zlaset( 'U', nr-1,nr-1, czero, czero, v(1,2), ldv ) + if (nr>1) call stdlib_zlaset( 'U', nr-1,nr-1, czero, czero, v(1,2), ldv ) call stdlib_zgesvj( 'L','U','N', n, nr, v, ldv, sva, nr, a, lda,cwork, lwork, & rwork, lrwork, info ) scalem = rwork(1) @@ -75689,17 +75689,17 @@ module stdlib_linalg_lapack_z else ! .. two more qr factorizations ( one qrf is not enough, two require ! accumulated product of jacobi rotations, three are perfect ) - call stdlib_zlaset( 'L', nr-1,nr-1, czero, czero, a(2,1), lda ) + if (nr>1) call stdlib_zlaset( 'L', nr-1,nr-1, czero, czero, a(2,1), lda ) call stdlib_zgelqf( nr,n, a, lda, cwork, cwork(n+1), lwork-n, ierr) call stdlib_zlacpy( 'L', nr, nr, a, lda, v, ldv ) - call stdlib_zlaset( 'U', nr-1,nr-1, czero, czero, v(1,2), ldv ) + if (nr>1) call stdlib_zlaset( 'U', nr-1,nr-1, czero, czero, v(1,2), ldv ) call stdlib_zgeqrf( nr, nr, v, ldv, cwork(n+1), cwork(2*n+1),lwork-2*n, ierr ) do p = 1, nr call stdlib_zcopy( nr-p+1, v(p,p), ldv, v(p,p), 1 ) call stdlib_zlacgv( nr-p+1, v(p,p), 1 ) end do - call stdlib_zlaset('U', nr-1, nr-1, czero, czero, v(1,2), ldv) + if (nr>1) call stdlib_zlaset('U', nr-1, nr-1, czero, czero, v(1,2), ldv) call stdlib_zgesvj( 'L', 'U','N', nr, nr, v,ldv, sva, nr, u,ldu, cwork(n+1), & lwork-n, rwork, lrwork, info ) scalem = rwork(1) @@ -75722,7 +75722,7 @@ module stdlib_linalg_lapack_z call stdlib_zlacpy( 'A', n, n, v, ldv, u, ldu ) end if else if ( jracc .and. (.not. lsvec) .and. ( nr== n ) ) then - call stdlib_zlaset( 'L', n-1,n-1, czero, czero, a(2,1), lda ) + if (n>1) call stdlib_zlaset( 'L', n-1,n-1, czero, czero, a(2,1), lda ) call stdlib_zgesvj( 'U','N','V', n, n, a, lda, sva, n, v, ldv,cwork, lwork, rwork, & lrwork, info ) scalem = rwork(1) @@ -75736,14 +75736,14 @@ module stdlib_linalg_lapack_z call stdlib_zcopy( n-p+1, a(p,p), lda, u(p,p), 1 ) call stdlib_zlacgv( n-p+1, u(p,p), 1 ) end do - call stdlib_zlaset( 'U', nr-1, nr-1, czero, czero, u(1,2), ldu ) + if (nr>1) call stdlib_zlaset( 'U', nr-1, nr-1, czero, czero, u(1,2), ldu ) call stdlib_zgeqrf( n, nr, u, ldu, cwork(n+1), cwork(2*n+1),lwork-2*n, ierr ) do p = 1, nr - 1 call stdlib_zcopy( nr-p, u(p,p+1), ldu, u(p+1,p), 1 ) call stdlib_zlacgv( n-p+1, u(p,p), 1 ) end do - call stdlib_zlaset( 'U', nr-1, nr-1, czero, czero, u(1,2), ldu ) + if (nr>1) call stdlib_zlaset( 'U', nr-1, nr-1, czero, czero, u(1,2), ldu ) call stdlib_zgesvj( 'L', 'U', 'N', nr,nr, u, ldu, sva, nr, a,lda, cwork(n+1), lwork-& n, rwork, lrwork, info ) scalem = rwork(1) @@ -75802,7 +75802,7 @@ module stdlib_linalg_lapack_z end do end do else - call stdlib_zlaset( 'U', nr-1, nr-1, czero, czero, v(1,2), ldv ) + if (nr>1) call stdlib_zlaset( 'U', nr-1, nr-1, czero, czero, v(1,2), ldv ) end if ! estimate the row scaled condition number of r1 ! (if r1 is rectangular, n > nr, then the condition number @@ -75884,7 +75884,7 @@ module stdlib_linalg_lapack_z end do end do else - call stdlib_zlaset( 'L',nr-1,nr-1,czero,czero,v(2,1),ldv ) + if (nr>1) call stdlib_zlaset( 'L',nr-1,nr-1,czero,czero,v(2,1),ldv ) end if ! now, compute r2 = l3 * q3, the lq factorization. call stdlib_zgelqf( nr, nr, v, ldv, cwork(2*n+n*nr+1),cwork(2*n+n*nr+nr+1), & @@ -75918,7 +75918,7 @@ module stdlib_linalg_lapack_z end do end do else - call stdlib_zlaset( 'U', nr-1,nr-1, czero,czero, v(1,2), ldv ) + if (nr>1) call stdlib_zlaset( 'U', nr-1,nr-1, czero,czero, v(1,2), ldv ) end if ! second preconditioning finished; continue with jacobi svd ! the input matrix is lower trinagular. @@ -76137,7 +76137,7 @@ module stdlib_linalg_lapack_z end do end do else - call stdlib_zlaset( 'U', nr-1, nr-1, czero, czero, v(1,2), ldv ) + if (nr>1) call stdlib_zlaset( 'U', nr-1, nr-1, czero, czero, v(1,2), ldv ) end if call stdlib_zgeqrf( n, nr, v, ldv, cwork(n+1), cwork(2*n+1),lwork-2*n, ierr ) @@ -76156,7 +76156,7 @@ module stdlib_linalg_lapack_z end do end do else - call stdlib_zlaset('U', nr-1, nr-1, czero, czero, u(1,2), ldu ) + if (nr>1) call stdlib_zlaset('U', nr-1, nr-1, czero, czero, u(1,2), ldu ) end if call stdlib_zgesvj( 'L', 'U', 'V', nr, nr, u, ldu, sva,n, v, ldv, cwork(2*n+n*nr+1),& lwork-2*n-n*nr,rwork, lrwork, info ) From bde2f3cd30bcf8e85beb7e78d53ab94144c30ae1 Mon Sep 17 00:00:00 2001 From: Federico Perini Date: Sat, 15 Jun 2024 10:21:01 -0500 Subject: [PATCH 3/3] add test --- test/linalg/test_linalg_svd.fypp | 37 ++++++++++++++++++++++++++++++++ 1 file changed, 37 insertions(+) diff --git a/test/linalg/test_linalg_svd.fypp b/test/linalg/test_linalg_svd.fypp index 9fe8a889a..d5a01d123 100644 --- a/test/linalg/test_linalg_svd.fypp +++ b/test/linalg/test_linalg_svd.fypp @@ -31,6 +31,12 @@ module test_linalg_svd #:endif #:endfor + #:for rk,rt,ri in RC_KINDS_TYPES + #:if rk!="xdp" + tests = [tests,new_unittest("test_svd_row_${ri}$",test_svd_row_${ri}$)] + #:endif + #:endfor + end subroutine test_svd !> Real matrix svd @@ -240,6 +246,37 @@ module test_linalg_svd #:endif #:endfor + + #:for rk,rt,ri in RC_KINDS_TYPES + #:if rk!="xdp" + ! Issue #835: bounds checking triggers an error with 1-sized A matrix + subroutine test_svd_row_${ri}$(error) + type(error_type), allocatable, intent(out) :: error + + !> Reference solution + type(linalg_state_type) :: state + integer(ilp), parameter :: m = 1, n = 1 + real(${rk}$), parameter :: tol = sqrt(epsilon(0.0_${rk}$)) + real(${rk}$) :: Arand(m, n), S(n) + ${rt}$ :: A(m, n), U(m, m), Vt(n, n) + + ! Random matrix. + call random_number(Arand) + A = Arand + + call svd(A, S, U, Vt, err=state) + + call check(error,state%ok(),'1-row SVD: '//state%print()) + if (allocated(error)) return + call check(error, abs(S(1)-A(1,1))