From 1b0d1a35823660f0c4340a985a01bbb705b2be24 Mon Sep 17 00:00:00 2001 From: julie Date: Sun, 17 Nov 2013 00:01:25 +0000 Subject: Patch provided by Lawrence Mulholland from Nag on Nov 1st 2013 Email below: ============ I have been incorporating some routines into the NAG Library, which means some automatic code translation and writing some example and test programs. The routines I have been adding are: ?geqrt, ?gemqrt, ?tpqrt, ?tpmqrt, ?orcsd, ?uncsd At the end of this message I will give you my current svn status and svn diff for consideration and approval before I commit. In each case, when testing immediate exits, my tests failed because constraints were mutually exclusive for the immediate return case. I have already committed changes to the constraints for some of the above to allow immediate exit. I have completed this for the remainder of this set. Less importantly, there are things in the code that trip up a checking compiler: a) an IF ( clause1(i) .AND. clause2(array(i)) ) THEN where array(i) is either not initialized or is out of bounds if clause1(i) is .FALSE. This is wrong since a Fortran compiler is at liberty to test clause2 first. In my changes this has been split into two as best suits the case. b) an CALL SUB (i, array(N-i+2)) with i = 1 and array(N+1) either not initialized or out of bounds, but internally array(N+1) is not referenced. In this case I don't think the Fortran standard is clear, but it trips up the nagfor compiler with checking on. So in the NAG incorporated versions of Lapack routines such calls are protected and/or a special i=1 call is made. The changes I want to commit also do this. c) workspace queries passing zero instead of array references e.g. lwork = -1 call barf(n,m,0,0,0,0,0,-1,info) a checking compiler won't like this. I have changed cases like this to pass available arrays of sufficient size and the right shape in place of the zeros. --- SRC/sorbdb.f | 176 +++++++++++++++++++++++++++++++++++++++++++---------------- 1 file changed, 129 insertions(+), 47 deletions(-) (limited to 'SRC/sorbdb.f') diff --git a/SRC/sorbdb.f b/SRC/sorbdb.f index 433b3c03..9f406bfc 100644 --- a/SRC/sorbdb.f +++ b/SRC/sorbdb.f @@ -415,19 +415,36 @@ THETA(I) = ATAN2( SNRM2( M-P-I+1, X21(I,I), 1 ), $ SNRM2( P-I+1, X11(I,I), 1 ) ) * - CALL SLARFGP( P-I+1, X11(I,I), X11(I+1,I), 1, TAUP1(I) ) + IF( P .GT. I ) THEN + CALL SLARFGP( P-I+1, X11(I,I), X11(I+1,I), 1, TAUP1(I) ) + ELSE IF( P .EQ. I ) THEN + CALL SLARFGP( P-I+1, X11(I,I), X11(I,I), 1, TAUP1(I) ) + END IF X11(I,I) = ONE - CALL SLARFGP( M-P-I+1, X21(I,I), X21(I+1,I), 1, TAUP2(I) ) + IF ( M-P .GT. I ) THEN + CALL SLARFGP( M-P-I+1, X21(I,I), X21(I+1,I), 1, + $ TAUP2(I) ) + ELSE IF ( M-P .EQ. I ) THEN + CALL SLARFGP( M-P-I+1, X21(I,I), X21(I,I), 1, TAUP2(I) ) + END IF X21(I,I) = ONE * - CALL SLARF( 'L', P-I+1, Q-I, X11(I,I), 1, TAUP1(I), - $ X11(I,I+1), LDX11, WORK ) - CALL SLARF( 'L', P-I+1, M-Q-I+1, X11(I,I), 1, TAUP1(I), - $ X12(I,I), LDX12, WORK ) - CALL SLARF( 'L', M-P-I+1, Q-I, X21(I,I), 1, TAUP2(I), - $ X21(I,I+1), LDX21, WORK ) - CALL SLARF( 'L', M-P-I+1, M-Q-I+1, X21(I,I), 1, TAUP2(I), - $ X22(I,I), LDX22, WORK ) + IF ( Q .GT. I ) THEN + CALL SLARF( 'L', P-I+1, Q-I, X11(I,I), 1, TAUP1(I), + $ X11(I,I+1), LDX11, WORK ) + END IF + IF ( M-Q+1 .GT. I ) THEN + CALL SLARF( 'L', P-I+1, M-Q-I+1, X11(I,I), 1, TAUP1(I), + $ X12(I,I), LDX12, WORK ) + END IF + IF ( Q .GT. I ) THEN + CALL SLARF( 'L', M-P-I+1, Q-I, X21(I,I), 1, TAUP2(I), + $ X21(I,I+1), LDX21, WORK ) + END IF + IF ( M-Q+1 .GT. I ) THEN + CALL SLARF( 'L', M-P-I+1, M-Q-I+1, X21(I,I), 1, TAUP2(I), + $ X22(I,I), LDX22, WORK ) + END IF * IF( I .LT. Q ) THEN CALL SSCAL( Q-I, -Z1*Z3*SIN(THETA(I)), X11(I,I+1), @@ -444,12 +461,24 @@ $ SNRM2( M-Q-I+1, X12(I,I), LDX12 ) ) * IF( I .LT. Q ) THEN - CALL SLARFGP( Q-I, X11(I,I+1), X11(I,I+2), LDX11, - $ TAUQ1(I) ) + IF ( Q-I .EQ. 1 ) THEN + CALL SLARFGP( Q-I, X11(I,I+1), X11(I,I+1), LDX11, + $ TAUQ1(I) ) + ELSE + CALL SLARFGP( Q-I, X11(I,I+1), X11(I,I+2), LDX11, + $ TAUQ1(I) ) + END IF X11(I,I+1) = ONE END IF - CALL SLARFGP( M-Q-I+1, X12(I,I), X12(I,I+1), LDX12, - $ TAUQ2(I) ) + IF ( Q+I-1 .LT. M ) THEN + IF ( M-Q .EQ. I ) THEN + CALL SLARFGP( M-Q-I+1, X12(I,I), X12(I,I), LDX12, + $ TAUQ2(I) ) + ELSE + CALL SLARFGP( M-Q-I+1, X12(I,I), X12(I,I+1), LDX12, + $ TAUQ2(I) ) + END IF + END IF X12(I,I) = ONE * IF( I .LT. Q ) THEN @@ -458,10 +487,14 @@ CALL SLARF( 'R', M-P-I, Q-I, X11(I,I+1), LDX11, TAUQ1(I), $ X21(I+1,I+1), LDX21, WORK ) END IF - CALL SLARF( 'R', P-I, M-Q-I+1, X12(I,I), LDX12, TAUQ2(I), - $ X12(I+1,I), LDX12, WORK ) - CALL SLARF( 'R', M-P-I, M-Q-I+1, X12(I,I), LDX12, TAUQ2(I), - $ X22(I+1,I), LDX22, WORK ) + IF ( P .GT. I ) THEN + CALL SLARF( 'R', P-I, M-Q-I+1, X12(I,I), LDX12, TAUQ2(I), + $ X12(I+1,I), LDX12, WORK ) + END IF + IF ( M-P .GT. I ) THEN + CALL SLARF( 'R', M-P-I, M-Q-I+1, X12(I,I), LDX12, + $ TAUQ2(I), X22(I+1,I), LDX22, WORK ) + END IF * END DO * @@ -470,12 +503,19 @@ DO I = Q + 1, P * CALL SSCAL( M-Q-I+1, -Z1*Z4, X12(I,I), LDX12 ) - CALL SLARFGP( M-Q-I+1, X12(I,I), X12(I,I+1), LDX12, - $ TAUQ2(I) ) + IF ( I .GE. M-Q ) THEN + CALL SLARFGP( M-Q-I+1, X12(I,I), X12(I,I), LDX12, + $ TAUQ2(I) ) + ELSE + CALL SLARFGP( M-Q-I+1, X12(I,I), X12(I,I+1), LDX12, + $ TAUQ2(I) ) + END IF X12(I,I) = ONE * - CALL SLARF( 'R', P-I, M-Q-I+1, X12(I,I), LDX12, TAUQ2(I), - $ X12(I+1,I), LDX12, WORK ) + IF ( P. GT. I ) THEN + CALL SLARF( 'R', P-I, M-Q-I+1, X12(I,I), LDX12, TAUQ2(I), + $ X12(I+1,I), LDX12, WORK ) + END IF IF( M-P-Q .GE. 1 ) $ CALL SLARF( 'R', M-P-Q, M-Q-I+1, X12(I,I), LDX12, $ TAUQ2(I), X22(Q+1,I), LDX22, WORK ) @@ -487,11 +527,18 @@ DO I = 1, M - P - Q * CALL SSCAL( M-P-Q-I+1, Z2*Z4, X22(Q+I,P+I), LDX22 ) - CALL SLARFGP( M-P-Q-I+1, X22(Q+I,P+I), X22(Q+I,P+I+1), - $ LDX22, TAUQ2(P+I) ) + IF ( I .EQ. M-P-Q ) THEN + CALL SLARFGP( M-P-Q-I+1, X22(Q+I,P+I), X22(Q+I,P+I), + $ LDX22, TAUQ2(P+I) ) + ELSE + CALL SLARFGP( M-P-Q-I+1, X22(Q+I,P+I), X22(Q+I,P+I+1), + $ LDX22, TAUQ2(P+I) ) + END IF X22(Q+I,P+I) = ONE - CALL SLARF( 'R', M-P-Q-I, M-P-Q-I+1, X22(Q+I,P+I), LDX22, - $ TAUQ2(P+I), X22(Q+I+1,P+I), LDX22, WORK ) + IF ( I .LT. M-P-Q ) THEN + CALL SLARF( 'R', M-P-Q-I, M-P-Q-I+1, X22(Q+I,P+I), LDX22, + $ TAUQ2(P+I), X22(Q+I+1,P+I), LDX22, WORK ) + END IF * END DO * @@ -521,18 +568,31 @@ * CALL SLARFGP( P-I+1, X11(I,I), X11(I,I+1), LDX11, TAUP1(I) ) X11(I,I) = ONE - CALL SLARFGP( M-P-I+1, X21(I,I), X21(I,I+1), LDX21, - $ TAUP2(I) ) + IF ( I .EQ. M-P ) THEN + CALL SLARFGP( M-P-I+1, X21(I,I), X21(I,I), LDX21, + $ TAUP2(I) ) + ELSE + CALL SLARFGP( M-P-I+1, X21(I,I), X21(I,I+1), LDX21, + $ TAUP2(I) ) + END IF X21(I,I) = ONE * - CALL SLARF( 'R', Q-I, P-I+1, X11(I,I), LDX11, TAUP1(I), - $ X11(I+1,I), LDX11, WORK ) - CALL SLARF( 'R', M-Q-I+1, P-I+1, X11(I,I), LDX11, TAUP1(I), - $ X12(I,I), LDX12, WORK ) - CALL SLARF( 'R', Q-I, M-P-I+1, X21(I,I), LDX21, TAUP2(I), - $ X21(I+1,I), LDX21, WORK ) - CALL SLARF( 'R', M-Q-I+1, M-P-I+1, X21(I,I), LDX21, - $ TAUP2(I), X22(I,I), LDX22, WORK ) + IF ( Q .GT. I ) THEN + CALL SLARF( 'R', Q-I, P-I+1, X11(I,I), LDX11, TAUP1(I), + $ X11(I+1,I), LDX11, WORK ) + END IF + IF ( M-Q+1 .GT. I ) THEN + CALL SLARF( 'R', M-Q-I+1, P-I+1, X11(I,I), LDX11, + $ TAUP1(I), X12(I,I), LDX12, WORK ) + END IF + IF ( Q .GT. I ) THEN + CALL SLARF( 'R', Q-I, M-P-I+1, X21(I,I), LDX21, TAUP2(I), + $ X21(I+1,I), LDX21, WORK ) + END IF + IF ( M-Q+1 .GT. I ) THEN + CALL SLARF( 'R', M-Q-I+1, M-P-I+1, X21(I,I), LDX21, + $ TAUP2(I), X22(I,I), LDX22, WORK ) + END IF * IF( I .LT. Q ) THEN CALL SSCAL( Q-I, -Z1*Z3*SIN(THETA(I)), X11(I+1,I), 1 ) @@ -548,10 +608,22 @@ $ SNRM2( M-Q-I+1, X12(I,I), 1 ) ) * IF( I .LT. Q ) THEN - CALL SLARFGP( Q-I, X11(I+1,I), X11(I+2,I), 1, TAUQ1(I) ) + IF ( Q-I .EQ. 1) THEN + CALL SLARFGP( Q-I, X11(I+1,I), X11(I+1,I), 1, + $ TAUQ1(I) ) + ELSE + CALL SLARFGP( Q-I, X11(I+1,I), X11(I+2,I), 1, + $ TAUQ1(I) ) + END IF X11(I+1,I) = ONE END IF - CALL SLARFGP( M-Q-I+1, X12(I,I), X12(I+1,I), 1, TAUQ2(I) ) + IF ( M-Q .GT. I ) THEN + CALL SLARFGP( M-Q-I+1, X12(I,I), X12(I+1,I), 1, + $ TAUQ2(I) ) + ELSE + CALL SLARFGP( M-Q-I+1, X12(I,I), X12(I,I), 1, + $ TAUQ2(I) ) + END IF X12(I,I) = ONE * IF( I .LT. Q ) THEN @@ -562,8 +634,10 @@ END IF CALL SLARF( 'L', M-Q-I+1, P-I, X12(I,I), 1, TAUQ2(I), $ X12(I,I+1), LDX12, WORK ) - CALL SLARF( 'L', M-Q-I+1, M-P-I, X12(I,I), 1, TAUQ2(I), - $ X22(I,I+1), LDX22, WORK ) + IF ( M-P-I .GT. 0 ) THEN + CALL SLARF( 'L', M-Q-I+1, M-P-I, X12(I,I), 1, TAUQ2(I), + $ X22(I,I+1), LDX22, WORK ) + END IF * END DO * @@ -575,8 +649,10 @@ CALL SLARFGP( M-Q-I+1, X12(I,I), X12(I+1,I), 1, TAUQ2(I) ) X12(I,I) = ONE * - CALL SLARF( 'L', M-Q-I+1, P-I, X12(I,I), 1, TAUQ2(I), - $ X12(I,I+1), LDX12, WORK ) + IF ( P .GT. I ) THEN + CALL SLARF( 'L', M-Q-I+1, P-I, X12(I,I), 1, TAUQ2(I), + $ X12(I,I+1), LDX12, WORK ) + END IF IF( M-P-Q .GE. 1 ) $ CALL SLARF( 'L', M-Q-I+1, M-P-Q, X12(I,I), 1, TAUQ2(I), $ X22(I,Q+1), LDX22, WORK ) @@ -588,12 +664,18 @@ DO I = 1, M - P - Q * CALL SSCAL( M-P-Q-I+1, Z2*Z4, X22(P+I,Q+I), 1 ) - CALL SLARFGP( M-P-Q-I+1, X22(P+I,Q+I), X22(P+I+1,Q+I), 1, - $ TAUQ2(P+I) ) - X22(P+I,Q+I) = ONE + IF ( M-P-Q .EQ. I ) THEN + CALL SLARFGP( M-P-Q-I+1, X22(P+I,Q+I), X22(P+I,Q+I), 1, + $ TAUQ2(P+I) ) + X22(P+I,Q+I) = ONE + ELSE + CALL SLARFGP( M-P-Q-I+1, X22(P+I,Q+I), X22(P+I+1,Q+I), 1, + $ TAUQ2(P+I) ) + X22(P+I,Q+I) = ONE + CALL SLARF( 'L', M-P-Q-I+1, M-P-Q-I, X22(P+I,Q+I), 1, + $ TAUQ2(P+I), X22(P+I,Q+I+1), LDX22, WORK ) + END IF * - CALL SLARF( 'L', M-P-Q-I+1, M-P-Q-I, X22(P+I,Q+I), 1, - $ TAUQ2(P+I), X22(P+I,Q+I+1), LDX22, WORK ) * END DO * -- cgit v1.2.3