Skip to content

2-by-1 CS decomposition fixes #405

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 3 commits into from
Feb 17, 2021
Merged

2-by-1 CS decomposition fixes #405

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
eb1921a
xUNCSD2BY1: improve LRWORK handling
christoph-conrads Apr 15, 2020
fe3cdc4
xORCSD2BY1, xUNCSD2BY1: fix U2 orthogonality
christoph-conrads Apr 21, 2020
0c32d4b
xUNCSD2BY1: allow workspace query with LRWORK=-1
christoph-conrads Apr 22, 2020
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
30 changes: 19 additions & 11 deletions SRC/cuncsd2by1.f
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -190,9 +190,10 @@
*> The dimension of the array WORK.
*>
*> If LWORK = -1, then a workspace query is assumed; the routine
*> only calculates the optimal size of the WORK array, returns
*> this value as the first entry of the work array, and no error
*> message related to LWORK is issued by XERBLA.
*> only calculates the optimal size of the WORK and RWORK
*> arrays, returns this value as the first entry of the WORK
*> and RWORK array, respectively, and no error message related
*> to LWORK or LRWORK is issued by XERBLA.
*> \endverbatim
*>
*> \param[out] RWORK
Expand All @@ -211,10 +212,11 @@
*> LRWORK is INTEGER
*> The dimension of the array RWORK.
*>
*> If LRWORK = -1, then a workspace query is assumed; the routine
*> only calculates the optimal size of the RWORK array, returns
*> this value as the first entry of the work array, and no error
*> message related to LRWORK is issued by XERBLA.
*> If LRWORK=-1, then a workspace query is assumed; the routine
*> only calculates the optimal size of the WORK and RWORK
*> arrays, returns this value as the first entry of the WORK
*> and RWORK array, respectively, and no error message related
*> to LWORK or LRWORK is issued by XERBLA.
*> \endverbatim
*
*> \param[out] IWORK
Expand Down Expand Up @@ -313,7 +315,7 @@ SUBROUTINE CUNCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
WANTU1 = LSAME( JOBU1, 'Y' )
WANTU2 = LSAME( JOBU2, 'Y' )
WANTV1T = LSAME( JOBV1T, 'Y' )
LQUERY = LWORK .EQ. -1
LQUERY = ( LWORK.EQ.-1 ) .OR. ( LRWORK.EQ.-1 )
*
IF( M .LT. 0 ) THEN
INFO = -4
Expand Down Expand Up @@ -513,6 +515,9 @@ SUBROUTINE CUNCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
IF( LWORK .LT. LWORKMIN .AND. .NOT.LQUERY ) THEN
INFO = -19
END IF
IF( LRWORK .LT. LRWORKMIN .AND. .NOT.LQUERY ) THEN
INFO = -21
END IF
END IF
IF( INFO .NE. 0 ) THEN
CALL XERBLA( 'CUNCSD2BY1', -INFO )
Expand Down Expand Up @@ -566,8 +571,8 @@ SUBROUTINE CUNCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
$ RWORK(IPHI), U1, LDU1, U2, LDU2, V1T, LDV1T, CDUM,
$ 1, RWORK(IB11D), RWORK(IB11E), RWORK(IB12D),
$ RWORK(IB12E), RWORK(IB21D), RWORK(IB21E),
$ RWORK(IB22D), RWORK(IB22E), RWORK(IBBCSD), LBBCSD,
$ CHILDINFO )
$ RWORK(IB22D), RWORK(IB22E), RWORK(IBBCSD),
$ LRWORK-IBBCSD+1, CHILDINFO )
*
* Permute rows and columns to place zero submatrices in
* preferred positions
Expand Down Expand Up @@ -708,6 +713,10 @@ SUBROUTINE CUNCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
*
* Accumulate Householder reflectors
*

IF( WANTU2 .AND. M-P .GT. 0 ) THEN
CALL CCOPY( M-P, WORK(IORBDB+P), 1, U2, 1 )
END IF
IF( WANTU1 .AND. P .GT. 0 ) THEN
CALL CCOPY( P, WORK(IORBDB), 1, U1, 1 )
DO J = 2, P
Expand All @@ -719,7 +728,6 @@ SUBROUTINE CUNCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
$ WORK(IORGQR), LORGQR, CHILDINFO )
END IF
IF( WANTU2 .AND. M-P .GT. 0 ) THEN
CALL CCOPY( M-P, WORK(IORBDB+P), 1, U2, 1 )
DO J = 2, M-P
U2(1,J) = ZERO
END DO
Expand Down
4 changes: 3 additions & 1 deletion SRC/dorcsd2by1.f
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -674,6 +674,9 @@ SUBROUTINE DORCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
*
* Accumulate Householder reflectors
*
IF( WANTU2 .AND. M-P .GT. 0 ) THEN
CALL DCOPY( M-P, WORK(IORBDB+P), 1, U2, 1 )
END IF
IF( WANTU1 .AND. P .GT. 0 ) THEN
CALL DCOPY( P, WORK(IORBDB), 1, U1, 1 )
DO J = 2, P
Expand All @@ -685,7 +688,6 @@ SUBROUTINE DORCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
$ WORK(IORGQR), LORGQR, CHILDINFO )
END IF
IF( WANTU2 .AND. M-P .GT. 0 ) THEN
CALL DCOPY( M-P, WORK(IORBDB+P), 1, U2, 1 )
DO J = 2, M-P
U2(1,J) = ZERO
END DO
Expand Down
4 changes: 3 additions & 1 deletion SRC/sorcsd2by1.f
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -674,6 +674,9 @@ SUBROUTINE SORCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
*
* Accumulate Householder reflectors
*
IF( WANTU2 .AND. M-P .GT. 0 ) THEN
CALL SCOPY( M-P, WORK(IORBDB+P), 1, U2, 1 )
END IF
IF( WANTU1 .AND. P .GT. 0 ) THEN
CALL SCOPY( P, WORK(IORBDB), 1, U1, 1 )
DO J = 2, P
Expand All @@ -685,7 +688,6 @@ SUBROUTINE SORCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
$ WORK(IORGQR), LORGQR, CHILDINFO )
END IF
IF( WANTU2 .AND. M-P .GT. 0 ) THEN
CALL SCOPY( M-P, WORK(IORBDB+P), 1, U2, 1 )
DO J = 2, M-P
U2(1,J) = ZERO
END DO
Expand Down
29 changes: 18 additions & 11 deletions SRC/zuncsd2by1.f
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -189,9 +189,10 @@
*> The dimension of the array WORK.
*>
*> If LWORK = -1, then a workspace query is assumed; the routine
*> only calculates the optimal size of the WORK array, returns
*> this value as the first entry of the work array, and no error
*> message related to LWORK is issued by XERBLA.
*> only calculates the optimal size of the WORK and RWORK
*> arrays, returns this value as the first entry of the WORK
*> and RWORK array, respectively, and no error message related
*> to LWORK or LRWORK is issued by XERBLA.
*> \endverbatim
*>
*> \param[out] RWORK
Expand All @@ -210,10 +211,11 @@
*> LRWORK is INTEGER
*> The dimension of the array RWORK.
*>
*> If LRWORK = -1, then a workspace query is assumed; the routine
*> only calculates the optimal size of the RWORK array, returns
*> this value as the first entry of the work array, and no error
*> message related to LRWORK is issued by XERBLA.
*> If LRWORK=-1, then a workspace query is assumed; the routine
*> only calculates the optimal size of the WORK and RWORK
*> arrays, returns this value as the first entry of the WORK
*> and RWORK array, respectively, and no error message related
*> to LWORK or LRWORK is issued by XERBLA.
*> \endverbatim
*
*> \param[out] IWORK
Expand Down Expand Up @@ -312,7 +314,7 @@ SUBROUTINE ZUNCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
WANTU1 = LSAME( JOBU1, 'Y' )
WANTU2 = LSAME( JOBU2, 'Y' )
WANTV1T = LSAME( JOBV1T, 'Y' )
LQUERY = LWORK .EQ. -1
LQUERY = ( LWORK.EQ.-1 ) .OR. ( LRWORK.EQ.-1 )
*
IF( M .LT. 0 ) THEN
INFO = -4
Expand Down Expand Up @@ -511,6 +513,9 @@ SUBROUTINE ZUNCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
IF( LWORK .LT. LWORKMIN .AND. .NOT.LQUERY ) THEN
INFO = -19
END IF
IF( LRWORK .LT. LRWORKMIN .AND. .NOT.LQUERY ) THEN
INFO = -21
END IF
END IF
IF( INFO .NE. 0 ) THEN
CALL XERBLA( 'ZUNCSD2BY1', -INFO )
Expand Down Expand Up @@ -564,8 +569,8 @@ SUBROUTINE ZUNCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
$ RWORK(IPHI), U1, LDU1, U2, LDU2, V1T, LDV1T, CDUM,
$ 1, RWORK(IB11D), RWORK(IB11E), RWORK(IB12D),
$ RWORK(IB12E), RWORK(IB21D), RWORK(IB21E),
$ RWORK(IB22D), RWORK(IB22E), RWORK(IBBCSD), LBBCSD,
$ CHILDINFO )
$ RWORK(IB22D), RWORK(IB22E), RWORK(IBBCSD),
$ LRWORK-IBBCSD+1, CHILDINFO )
*
* Permute rows and columns to place zero submatrices in
* preferred positions
Expand Down Expand Up @@ -706,6 +711,9 @@ SUBROUTINE ZUNCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
*
* Accumulate Householder reflectors
*
IF( WANTU2 .AND. M-P .GT. 0 ) THEN
CALL ZCOPY( M-P, WORK(IORBDB+P), 1, U2, 1 )
END IF
IF( WANTU1 .AND. P .GT. 0 ) THEN
CALL ZCOPY( P, WORK(IORBDB), 1, U1, 1 )
DO J = 2, P
Expand All @@ -717,7 +725,6 @@ SUBROUTINE ZUNCSD2BY1( JOBU1, JOBU2, JOBV1T, M, P, Q, X11, LDX11,
$ WORK(IORGQR), LORGQR, CHILDINFO )
END IF
IF( WANTU2 .AND. M-P .GT. 0 ) THEN
CALL ZCOPY( M-P, WORK(IORBDB+P), 1, U2, 1 )
DO J = 2, M-P
U2(1,J) = ZERO
END DO
Expand Down