--- rpl/lapack/lapack/dgesvj.f 2010/12/21 13:53:26 1.5 +++ rpl/lapack/lapack/dgesvj.f 2011/07/22 07:38:05 1.6 @@ -1,11 +1,11 @@ SUBROUTINE DGESVJ( JOBA, JOBU, JOBV, M, N, A, LDA, SVA, MV, V, - + LDV, WORK, LWORK, INFO ) + $ LDV, WORK, LWORK, INFO ) * -* -- LAPACK routine (version 3.3.0) -- +* -- LAPACK routine (version 3.3.1) -- * * -- Contributed by Zlatko Drmac of the University of Zagreb and -- * -- Kresimir Veselic of the Fernuniversitaet Hagen -- -* November 2010 +* -- April 2011 -- * * -- LAPACK is a software package provided by Univ. of Tennessee, -- * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- @@ -23,7 +23,7 @@ * .. * .. Array Arguments .. DOUBLE PRECISION A( LDA, * ), SVA( N ), V( LDV, * ), - + WORK( LWORK ) + $ WORK( LWORK ) * .. * * Purpose @@ -133,7 +133,7 @@ * referenced * * M (input) INTEGER -* The number of rows of the input matrix A. M >= 0. +* The number of rows of the input matrix A. 1/DLAMCH('E') > M >= 0. * * N (input) INTEGER * The number of columns of the input matrix A. @@ -256,22 +256,22 @@ * .. Local Parameters .. DOUBLE PRECISION ZERO, HALF, ONE, TWO PARAMETER ( ZERO = 0.0D0, HALF = 0.5D0, ONE = 1.0D0, - + TWO = 2.0D0 ) + $ TWO = 2.0D0 ) INTEGER NSWEEP PARAMETER ( NSWEEP = 30 ) * .. * .. Local Scalars .. DOUBLE PRECISION AAPP, AAPP0, AAPQ, AAQQ, APOAQ, AQOAP, BIG, - + BIGTHETA, CS, CTOL, EPSLN, LARGE, MXAAPQ, - + MXSINJ, ROOTBIG, ROOTEPS, ROOTSFMIN, ROOTTOL, - + SKL, SFMIN, SMALL, SN, T, TEMP1, THETA, - + THSIGN, TOL + $ BIGTHETA, CS, CTOL, EPSLN, LARGE, MXAAPQ, + $ MXSINJ, ROOTBIG, ROOTEPS, ROOTSFMIN, ROOTTOL, + $ SKL, SFMIN, SMALL, SN, T, TEMP1, THETA, + $ THSIGN, TOL INTEGER BLSKIP, EMPTSW, i, ibr, IERR, igl, IJBLSK, ir1, - + ISWROT, jbc, jgl, KBL, LKAHEAD, MVL, N2, N34, - + N4, NBL, NOTROT, p, PSKIPPED, q, ROWSKIP, - + SWBAND + $ ISWROT, jbc, jgl, KBL, LKAHEAD, MVL, N2, N34, + $ N4, NBL, NOTROT, p, PSKIPPED, q, ROWSKIP, + $ SWBAND LOGICAL APPLV, GOSCALE, LOWER, LSVEC, NOSCALE, ROTOK, - + RSVEC, UCTOL, UPPER + $ RSVEC, UCTOL, UPPER * .. * .. Local Arrays .. DOUBLE PRECISION FASTR( 5 ) @@ -327,7 +327,7 @@ ELSE IF( MV.LT.0 ) THEN INFO = -9 ELSE IF( ( RSVEC .AND. ( LDV.LT.N ) ) .OR. - + ( APPLV .AND. ( LDV.LT.MV ) ) ) THEN + $ ( APPLV .AND. ( LDV.LT.MV ) ) ) THEN INFO = -11 ELSE IF( UCTOL .AND. ( WORK( 1 ).LE.ONE ) ) THEN INFO = -12 @@ -383,7 +383,7 @@ ROOTTOL = DSQRT( TOL ) * IF( DBLE( M )*EPSLN.GE.ONE ) THEN - INFO = -5 + INFO = -4 CALL XERBLA( 'DGESVJ', -INFO ) RETURN END IF @@ -518,7 +518,7 @@ * IF( N.EQ.1 ) THEN IF( LSVEC )CALL DLASCL( 'G', 0, 0, SVA( 1 ), SKL, M, 1, - + A( 1, 1 ), LDA, IERR ) + $ A( 1, 1 ), LDA, IERR ) WORK( 1 ) = ONE / SKL IF( SVA( 1 ).GE.SFMIN ) THEN WORK( 2 ) = ONE @@ -538,7 +538,7 @@ SN = DSQRT( SFMIN / EPSLN ) TEMP1 = DSQRT( BIG / DBLE( N ) ) IF( ( AAPP.LE.SN ) .OR. ( AAQQ.GE.TEMP1 ) .OR. - + ( ( SN.LE.AAQQ ) .AND. ( AAPP.LE.TEMP1 ) ) ) THEN + $ ( ( SN.LE.AAQQ ) .AND. ( AAPP.LE.TEMP1 ) ) ) THEN TEMP1 = DMIN1( BIG, TEMP1 / AAPP ) * AAQQ = AAQQ*TEMP1 * AAPP = AAPP*TEMP1 @@ -638,54 +638,54 @@ * [+ + x x] [x x]. [x x] * CALL DGSVJ0( JOBV, M-N34, N-N34, A( N34+1, N34+1 ), LDA, - + WORK( N34+1 ), SVA( N34+1 ), MVL, - + V( N34*q+1, N34+1 ), LDV, EPSLN, SFMIN, TOL, - + 2, WORK( N+1 ), LWORK-N, IERR ) + $ WORK( N34+1 ), SVA( N34+1 ), MVL, + $ V( N34*q+1, N34+1 ), LDV, EPSLN, SFMIN, TOL, + $ 2, WORK( N+1 ), LWORK-N, IERR ) * CALL DGSVJ0( JOBV, M-N2, N34-N2, A( N2+1, N2+1 ), LDA, - + WORK( N2+1 ), SVA( N2+1 ), MVL, - + V( N2*q+1, N2+1 ), LDV, EPSLN, SFMIN, TOL, 2, - + WORK( N+1 ), LWORK-N, IERR ) + $ WORK( N2+1 ), SVA( N2+1 ), MVL, + $ V( N2*q+1, N2+1 ), LDV, EPSLN, SFMIN, TOL, 2, + $ WORK( N+1 ), LWORK-N, IERR ) * CALL DGSVJ1( JOBV, M-N2, N-N2, N4, A( N2+1, N2+1 ), LDA, - + WORK( N2+1 ), SVA( N2+1 ), MVL, - + V( N2*q+1, N2+1 ), LDV, EPSLN, SFMIN, TOL, 1, - + WORK( N+1 ), LWORK-N, IERR ) + $ WORK( N2+1 ), SVA( N2+1 ), MVL, + $ V( N2*q+1, N2+1 ), LDV, EPSLN, SFMIN, TOL, 1, + $ WORK( N+1 ), LWORK-N, IERR ) * CALL DGSVJ0( JOBV, M-N4, N2-N4, A( N4+1, N4+1 ), LDA, - + WORK( N4+1 ), SVA( N4+1 ), MVL, - + V( N4*q+1, N4+1 ), LDV, EPSLN, SFMIN, TOL, 1, - + WORK( N+1 ), LWORK-N, IERR ) + $ WORK( N4+1 ), SVA( N4+1 ), MVL, + $ V( N4*q+1, N4+1 ), LDV, EPSLN, SFMIN, TOL, 1, + $ WORK( N+1 ), LWORK-N, IERR ) * CALL DGSVJ0( JOBV, M, N4, A, LDA, WORK, SVA, MVL, V, LDV, - + EPSLN, SFMIN, TOL, 1, WORK( N+1 ), LWORK-N, - + IERR ) + $ EPSLN, SFMIN, TOL, 1, WORK( N+1 ), LWORK-N, + $ IERR ) * CALL DGSVJ1( JOBV, M, N2, N4, A, LDA, WORK, SVA, MVL, V, - + LDV, EPSLN, SFMIN, TOL, 1, WORK( N+1 ), - + LWORK-N, IERR ) + $ LDV, EPSLN, SFMIN, TOL, 1, WORK( N+1 ), + $ LWORK-N, IERR ) * * ELSE IF( UPPER ) THEN * * CALL DGSVJ0( JOBV, N4, N4, A, LDA, WORK, SVA, MVL, V, LDV, - + EPSLN, SFMIN, TOL, 2, WORK( N+1 ), LWORK-N, - + IERR ) + $ EPSLN, SFMIN, TOL, 2, WORK( N+1 ), LWORK-N, + $ IERR ) * CALL DGSVJ0( JOBV, N2, N4, A( 1, N4+1 ), LDA, WORK( N4+1 ), - + SVA( N4+1 ), MVL, V( N4*q+1, N4+1 ), LDV, - + EPSLN, SFMIN, TOL, 1, WORK( N+1 ), LWORK-N, - + IERR ) + $ SVA( N4+1 ), MVL, V( N4*q+1, N4+1 ), LDV, + $ EPSLN, SFMIN, TOL, 1, WORK( N+1 ), LWORK-N, + $ IERR ) * CALL DGSVJ1( JOBV, N2, N2, N4, A, LDA, WORK, SVA, MVL, V, - + LDV, EPSLN, SFMIN, TOL, 1, WORK( N+1 ), - + LWORK-N, IERR ) + $ LDV, EPSLN, SFMIN, TOL, 1, WORK( N+1 ), + $ LWORK-N, IERR ) * CALL DGSVJ0( JOBV, N2+N4, N4, A( 1, N2+1 ), LDA, - + WORK( N2+1 ), SVA( N2+1 ), MVL, - + V( N2*q+1, N2+1 ), LDV, EPSLN, SFMIN, TOL, 1, - + WORK( N+1 ), LWORK-N, IERR ) + $ WORK( N2+1 ), SVA( N2+1 ), MVL, + $ V( N2*q+1, N2+1 ), LDV, EPSLN, SFMIN, TOL, 1, + $ WORK( N+1 ), LWORK-N, IERR ) END IF * @@ -725,7 +725,7 @@ IF( p.NE.q ) THEN CALL DSWAP( M, A( 1, p ), 1, A( 1, q ), 1 ) IF( RSVEC )CALL DSWAP( MVL, V( 1, p ), 1, - + V( 1, q ), 1 ) + $ V( 1, q ), 1 ) TEMP1 = SVA( p ) SVA( p ) = SVA( q ) SVA( q ) = TEMP1 @@ -749,7 +749,7 @@ * below should read "AAPP = DNRM2( M, A(1,p), 1 ) * WORK(p)". * IF( ( SVA( p ).LT.ROOTBIG ) .AND. - + ( SVA( p ).GT.ROOTSFMIN ) ) THEN + $ ( SVA( p ).GT.ROOTSFMIN ) ) THEN SVA( p ) = DNRM2( M, A( 1, p ), 1 )*WORK( p ) ELSE TEMP1 = ZERO @@ -777,31 +777,31 @@ ROTOK = ( SMALL*AAPP ).LE.AAQQ IF( AAPP.LT.( BIG / AAQQ ) ) THEN AAPQ = ( DDOT( M, A( 1, p ), 1, A( 1, - + q ), 1 )*WORK( p )*WORK( q ) / - + AAQQ ) / AAPP + $ q ), 1 )*WORK( p )*WORK( q ) / + $ AAQQ ) / AAPP ELSE CALL DCOPY( M, A( 1, p ), 1, - + WORK( N+1 ), 1 ) + $ WORK( N+1 ), 1 ) CALL DLASCL( 'G', 0, 0, AAPP, - + WORK( p ), M, 1, - + WORK( N+1 ), LDA, IERR ) + $ WORK( p ), M, 1, + $ WORK( N+1 ), LDA, IERR ) AAPQ = DDOT( M, WORK( N+1 ), 1, - + A( 1, q ), 1 )*WORK( q ) / AAQQ + $ A( 1, q ), 1 )*WORK( q ) / AAQQ END IF ELSE ROTOK = AAPP.LE.( AAQQ / SMALL ) IF( AAPP.GT.( SMALL / AAQQ ) ) THEN AAPQ = ( DDOT( M, A( 1, p ), 1, A( 1, - + q ), 1 )*WORK( p )*WORK( q ) / - + AAQQ ) / AAPP + $ q ), 1 )*WORK( p )*WORK( q ) / + $ AAQQ ) / AAPP ELSE CALL DCOPY( M, A( 1, q ), 1, - + WORK( N+1 ), 1 ) + $ WORK( N+1 ), 1 ) CALL DLASCL( 'G', 0, 0, AAQQ, - + WORK( q ), M, 1, - + WORK( N+1 ), LDA, IERR ) + $ WORK( q ), M, 1, + $ WORK( N+1 ), LDA, IERR ) AAPQ = DDOT( M, WORK( N+1 ), 1, - + A( 1, p ), 1 )*WORK( p ) / AAPP + $ A( 1, p ), 1 )*WORK( p ) / AAPP END IF END IF * @@ -824,25 +824,24 @@ * AQOAP = AAQQ / AAPP APOAQ = AAPP / AAQQ - THETA = -HALF*DABS( AQOAP-APOAQ ) / - + AAPQ + THETA = -HALF*DABS(AQOAP-APOAQ)/AAPQ * IF( DABS( THETA ).GT.BIGTHETA ) THEN * T = HALF / THETA FASTR( 3 ) = T*WORK( p ) / WORK( q ) FASTR( 4 ) = -T*WORK( q ) / - + WORK( p ) + $ WORK( p ) CALL DROTM( M, A( 1, p ), 1, - + A( 1, q ), 1, FASTR ) + $ A( 1, q ), 1, FASTR ) IF( RSVEC )CALL DROTM( MVL, - + V( 1, p ), 1, - + V( 1, q ), 1, - + FASTR ) + $ V( 1, p ), 1, + $ V( 1, q ), 1, + $ FASTR ) SVA( q ) = AAQQ*DSQRT( DMAX1( ZERO, - + ONE+T*APOAQ*AAPQ ) ) + $ ONE+T*APOAQ*AAPQ ) ) AAPP = AAPP*DSQRT( DMAX1( ZERO, - + ONE-T*AQOAP*AAPQ ) ) + $ ONE-T*AQOAP*AAPQ ) ) MXSINJ = DMAX1( MXSINJ, DABS( T ) ) * ELSE @@ -851,15 +850,15 @@ * THSIGN = -DSIGN( ONE, AAPQ ) T = ONE / ( THETA+THSIGN* - + DSQRT( ONE+THETA*THETA ) ) + $ DSQRT( ONE+THETA*THETA ) ) CS = DSQRT( ONE / ( ONE+T*T ) ) SN = T*CS * MXSINJ = DMAX1( MXSINJ, DABS( SN ) ) SVA( q ) = AAQQ*DSQRT( DMAX1( ZERO, - + ONE+T*APOAQ*AAPQ ) ) + $ ONE+T*APOAQ*AAPQ ) ) AAPP = AAPP*DSQRT( DMAX1( ZERO, - + ONE-T*AQOAP*AAPQ ) ) + $ ONE-T*AQOAP*AAPQ ) ) * APOAQ = WORK( p ) / WORK( q ) AQOAP = WORK( q ) / WORK( p ) @@ -870,88 +869,88 @@ WORK( p ) = WORK( p )*CS WORK( q ) = WORK( q )*CS CALL DROTM( M, A( 1, p ), 1, - + A( 1, q ), 1, - + FASTR ) + $ A( 1, q ), 1, + $ FASTR ) IF( RSVEC )CALL DROTM( MVL, - + V( 1, p ), 1, V( 1, q ), - + 1, FASTR ) + $ V( 1, p ), 1, V( 1, q ), + $ 1, FASTR ) ELSE CALL DAXPY( M, -T*AQOAP, - + A( 1, q ), 1, - + A( 1, p ), 1 ) + $ A( 1, q ), 1, + $ A( 1, p ), 1 ) CALL DAXPY( M, CS*SN*APOAQ, - + A( 1, p ), 1, - + A( 1, q ), 1 ) + $ A( 1, p ), 1, + $ A( 1, q ), 1 ) WORK( p ) = WORK( p )*CS WORK( q ) = WORK( q ) / CS IF( RSVEC ) THEN CALL DAXPY( MVL, -T*AQOAP, - + V( 1, q ), 1, - + V( 1, p ), 1 ) + $ V( 1, q ), 1, + $ V( 1, p ), 1 ) CALL DAXPY( MVL, - + CS*SN*APOAQ, - + V( 1, p ), 1, - + V( 1, q ), 1 ) + $ CS*SN*APOAQ, + $ V( 1, p ), 1, + $ V( 1, q ), 1 ) END IF END IF ELSE IF( WORK( q ).GE.ONE ) THEN CALL DAXPY( M, T*APOAQ, - + A( 1, p ), 1, - + A( 1, q ), 1 ) + $ A( 1, p ), 1, + $ A( 1, q ), 1 ) CALL DAXPY( M, -CS*SN*AQOAP, - + A( 1, q ), 1, - + A( 1, p ), 1 ) + $ A( 1, q ), 1, + $ A( 1, p ), 1 ) WORK( p ) = WORK( p ) / CS WORK( q ) = WORK( q )*CS IF( RSVEC ) THEN CALL DAXPY( MVL, T*APOAQ, - + V( 1, p ), 1, - + V( 1, q ), 1 ) + $ V( 1, p ), 1, + $ V( 1, q ), 1 ) CALL DAXPY( MVL, - + -CS*SN*AQOAP, - + V( 1, q ), 1, - + V( 1, p ), 1 ) + $ -CS*SN*AQOAP, + $ V( 1, q ), 1, + $ V( 1, p ), 1 ) END IF ELSE IF( WORK( p ).GE.WORK( q ) ) - + THEN + $ THEN CALL DAXPY( M, -T*AQOAP, - + A( 1, q ), 1, - + A( 1, p ), 1 ) + $ A( 1, q ), 1, + $ A( 1, p ), 1 ) CALL DAXPY( M, CS*SN*APOAQ, - + A( 1, p ), 1, - + A( 1, q ), 1 ) + $ A( 1, p ), 1, + $ A( 1, q ), 1 ) WORK( p ) = WORK( p )*CS WORK( q ) = WORK( q ) / CS IF( RSVEC ) THEN CALL DAXPY( MVL, - + -T*AQOAP, - + V( 1, q ), 1, - + V( 1, p ), 1 ) + $ -T*AQOAP, + $ V( 1, q ), 1, + $ V( 1, p ), 1 ) CALL DAXPY( MVL, - + CS*SN*APOAQ, - + V( 1, p ), 1, - + V( 1, q ), 1 ) + $ CS*SN*APOAQ, + $ V( 1, p ), 1, + $ V( 1, q ), 1 ) END IF ELSE CALL DAXPY( M, T*APOAQ, - + A( 1, p ), 1, - + A( 1, q ), 1 ) + $ A( 1, p ), 1, + $ A( 1, q ), 1 ) CALL DAXPY( M, - + -CS*SN*AQOAP, - + A( 1, q ), 1, - + A( 1, p ), 1 ) + $ -CS*SN*AQOAP, + $ A( 1, q ), 1, + $ A( 1, p ), 1 ) WORK( p ) = WORK( p ) / CS WORK( q ) = WORK( q )*CS IF( RSVEC ) THEN CALL DAXPY( MVL, - + T*APOAQ, V( 1, p ), - + 1, V( 1, q ), 1 ) + $ T*APOAQ, V( 1, p ), + $ 1, V( 1, q ), 1 ) CALL DAXPY( MVL, - + -CS*SN*AQOAP, - + V( 1, q ), 1, - + V( 1, p ), 1 ) + $ -CS*SN*AQOAP, + $ V( 1, q ), 1, + $ V( 1, p ), 1 ) END IF END IF END IF @@ -961,19 +960,19 @@ ELSE * .. have to use modified Gram-Schmidt like transformation CALL DCOPY( M, A( 1, p ), 1, - + WORK( N+1 ), 1 ) + $ WORK( N+1 ), 1 ) CALL DLASCL( 'G', 0, 0, AAPP, ONE, M, - + 1, WORK( N+1 ), LDA, - + IERR ) + $ 1, WORK( N+1 ), LDA, + $ IERR ) CALL DLASCL( 'G', 0, 0, AAQQ, ONE, M, - + 1, A( 1, q ), LDA, IERR ) + $ 1, A( 1, q ), LDA, IERR ) TEMP1 = -AAPQ*WORK( p ) / WORK( q ) CALL DAXPY( M, TEMP1, WORK( N+1 ), 1, - + A( 1, q ), 1 ) + $ A( 1, q ), 1 ) CALL DLASCL( 'G', 0, 0, ONE, AAQQ, M, - + 1, A( 1, q ), LDA, IERR ) + $ 1, A( 1, q ), LDA, IERR ) SVA( q ) = AAQQ*DSQRT( DMAX1( ZERO, - + ONE-AAPQ*AAPQ ) ) + $ ONE-AAPQ*AAPQ ) ) MXSINJ = DMAX1( MXSINJ, SFMIN ) END IF * END IF ROTOK THEN ... ELSE @@ -982,29 +981,29 @@ * recompute SVA(q), SVA(p). * IF( ( SVA( q ) / AAQQ )**2.LE.ROOTEPS ) - + THEN + $ THEN IF( ( AAQQ.LT.ROOTBIG ) .AND. - + ( AAQQ.GT.ROOTSFMIN ) ) THEN + $ ( AAQQ.GT.ROOTSFMIN ) ) THEN SVA( q ) = DNRM2( M, A( 1, q ), 1 )* - + WORK( q ) + $ WORK( q ) ELSE T = ZERO AAQQ = ONE CALL DLASSQ( M, A( 1, q ), 1, T, - + AAQQ ) + $ AAQQ ) SVA( q ) = T*DSQRT( AAQQ )*WORK( q ) END IF END IF IF( ( AAPP / AAPP0 ).LE.ROOTEPS ) THEN IF( ( AAPP.LT.ROOTBIG ) .AND. - + ( AAPP.GT.ROOTSFMIN ) ) THEN + $ ( AAPP.GT.ROOTSFMIN ) ) THEN AAPP = DNRM2( M, A( 1, p ), 1 )* - + WORK( p ) + $ WORK( p ) ELSE T = ZERO AAPP = ONE CALL DLASSQ( M, A( 1, p ), 1, T, - + AAPP ) + $ AAPP ) AAPP = T*DSQRT( AAPP )*WORK( p ) END IF SVA( p ) = AAPP @@ -1023,7 +1022,7 @@ END IF * IF( ( i.LE.SWBAND ) .AND. - + ( PSKIPPED.GT.ROWSKIP ) ) THEN + $ ( PSKIPPED.GT.ROWSKIP ) ) THEN IF( ir1.EQ.0 )AAPP = -AAPP NOTROT = 0 GO TO 2103 @@ -1040,7 +1039,7 @@ ELSE SVA( p ) = AAPP IF( ( ir1.EQ.0 ) .AND. ( AAPP.EQ.ZERO ) ) - + NOTROT = NOTROT + MIN0( igl+KBL-1, N ) - p + $ NOTROT = NOTROT + MIN0( igl+KBL-1, N ) - p END IF * 2001 CONTINUE @@ -1085,16 +1084,16 @@ END IF IF( AAPP.LT.( BIG / AAQQ ) ) THEN AAPQ = ( DDOT( M, A( 1, p ), 1, A( 1, - + q ), 1 )*WORK( p )*WORK( q ) / - + AAQQ ) / AAPP + $ q ), 1 )*WORK( p )*WORK( q ) / + $ AAQQ ) / AAPP ELSE CALL DCOPY( M, A( 1, p ), 1, - + WORK( N+1 ), 1 ) + $ WORK( N+1 ), 1 ) CALL DLASCL( 'G', 0, 0, AAPP, - + WORK( p ), M, 1, - + WORK( N+1 ), LDA, IERR ) + $ WORK( p ), M, 1, + $ WORK( N+1 ), LDA, IERR ) AAPQ = DDOT( M, WORK( N+1 ), 1, - + A( 1, q ), 1 )*WORK( q ) / AAQQ + $ A( 1, q ), 1 )*WORK( q ) / AAQQ END IF ELSE IF( AAPP.GE.AAQQ ) THEN @@ -1104,16 +1103,16 @@ END IF IF( AAPP.GT.( SMALL / AAQQ ) ) THEN AAPQ = ( DDOT( M, A( 1, p ), 1, A( 1, - + q ), 1 )*WORK( p )*WORK( q ) / - + AAQQ ) / AAPP + $ q ), 1 )*WORK( p )*WORK( q ) / + $ AAQQ ) / AAPP ELSE CALL DCOPY( M, A( 1, q ), 1, - + WORK( N+1 ), 1 ) + $ WORK( N+1 ), 1 ) CALL DLASCL( 'G', 0, 0, AAQQ, - + WORK( q ), M, 1, - + WORK( N+1 ), LDA, IERR ) + $ WORK( q ), M, 1, + $ WORK( N+1 ), LDA, IERR ) AAPQ = DDOT( M, WORK( N+1 ), 1, - + A( 1, p ), 1 )*WORK( p ) / AAPP + $ A( 1, p ), 1 )*WORK( p ) / AAPP END IF END IF * @@ -1131,25 +1130,24 @@ * AQOAP = AAQQ / AAPP APOAQ = AAPP / AAQQ - THETA = -HALF*DABS( AQOAP-APOAQ ) / - + AAPQ + THETA = -HALF*DABS(AQOAP-APOAQ)/AAPQ IF( AAQQ.GT.AAPP0 )THETA = -THETA * IF( DABS( THETA ).GT.BIGTHETA ) THEN T = HALF / THETA FASTR( 3 ) = T*WORK( p ) / WORK( q ) FASTR( 4 ) = -T*WORK( q ) / - + WORK( p ) + $ WORK( p ) CALL DROTM( M, A( 1, p ), 1, - + A( 1, q ), 1, FASTR ) + $ A( 1, q ), 1, FASTR ) IF( RSVEC )CALL DROTM( MVL, - + V( 1, p ), 1, - + V( 1, q ), 1, - + FASTR ) + $ V( 1, p ), 1, + $ V( 1, q ), 1, + $ FASTR ) SVA( q ) = AAQQ*DSQRT( DMAX1( ZERO, - + ONE+T*APOAQ*AAPQ ) ) + $ ONE+T*APOAQ*AAPQ ) ) AAPP = AAPP*DSQRT( DMAX1( ZERO, - + ONE-T*AQOAP*AAPQ ) ) + $ ONE-T*AQOAP*AAPQ ) ) MXSINJ = DMAX1( MXSINJ, DABS( T ) ) ELSE * @@ -1158,14 +1156,14 @@ THSIGN = -DSIGN( ONE, AAPQ ) IF( AAQQ.GT.AAPP0 )THSIGN = -THSIGN T = ONE / ( THETA+THSIGN* - + DSQRT( ONE+THETA*THETA ) ) + $ DSQRT( ONE+THETA*THETA ) ) CS = DSQRT( ONE / ( ONE+T*T ) ) SN = T*CS MXSINJ = DMAX1( MXSINJ, DABS( SN ) ) SVA( q ) = AAQQ*DSQRT( DMAX1( ZERO, - + ONE+T*APOAQ*AAPQ ) ) + $ ONE+T*APOAQ*AAPQ ) ) AAPP = AAPP*DSQRT( DMAX1( ZERO, - + ONE-T*AQOAP*AAPQ ) ) + $ ONE-T*AQOAP*AAPQ ) ) * APOAQ = WORK( p ) / WORK( q ) AQOAP = WORK( q ) / WORK( p ) @@ -1177,26 +1175,26 @@ WORK( p ) = WORK( p )*CS WORK( q ) = WORK( q )*CS CALL DROTM( M, A( 1, p ), 1, - + A( 1, q ), 1, - + FASTR ) + $ A( 1, q ), 1, + $ FASTR ) IF( RSVEC )CALL DROTM( MVL, - + V( 1, p ), 1, V( 1, q ), - + 1, FASTR ) + $ V( 1, p ), 1, V( 1, q ), + $ 1, FASTR ) ELSE CALL DAXPY( M, -T*AQOAP, - + A( 1, q ), 1, - + A( 1, p ), 1 ) + $ A( 1, q ), 1, + $ A( 1, p ), 1 ) CALL DAXPY( M, CS*SN*APOAQ, - + A( 1, p ), 1, - + A( 1, q ), 1 ) + $ A( 1, p ), 1, + $ A( 1, q ), 1 ) IF( RSVEC ) THEN CALL DAXPY( MVL, -T*AQOAP, - + V( 1, q ), 1, - + V( 1, p ), 1 ) + $ V( 1, q ), 1, + $ V( 1, p ), 1 ) CALL DAXPY( MVL, - + CS*SN*APOAQ, - + V( 1, p ), 1, - + V( 1, q ), 1 ) + $ CS*SN*APOAQ, + $ V( 1, p ), 1, + $ V( 1, q ), 1 ) END IF WORK( p ) = WORK( p )*CS WORK( q ) = WORK( q ) / CS @@ -1204,61 +1202,61 @@ ELSE IF( WORK( q ).GE.ONE ) THEN CALL DAXPY( M, T*APOAQ, - + A( 1, p ), 1, - + A( 1, q ), 1 ) + $ A( 1, p ), 1, + $ A( 1, q ), 1 ) CALL DAXPY( M, -CS*SN*AQOAP, - + A( 1, q ), 1, - + A( 1, p ), 1 ) + $ A( 1, q ), 1, + $ A( 1, p ), 1 ) IF( RSVEC ) THEN CALL DAXPY( MVL, T*APOAQ, - + V( 1, p ), 1, - + V( 1, q ), 1 ) + $ V( 1, p ), 1, + $ V( 1, q ), 1 ) CALL DAXPY( MVL, - + -CS*SN*AQOAP, - + V( 1, q ), 1, - + V( 1, p ), 1 ) + $ -CS*SN*AQOAP, + $ V( 1, q ), 1, + $ V( 1, p ), 1 ) END IF WORK( p ) = WORK( p ) / CS WORK( q ) = WORK( q )*CS ELSE IF( WORK( p ).GE.WORK( q ) ) - + THEN + $ THEN CALL DAXPY( M, -T*AQOAP, - + A( 1, q ), 1, - + A( 1, p ), 1 ) + $ A( 1, q ), 1, + $ A( 1, p ), 1 ) CALL DAXPY( M, CS*SN*APOAQ, - + A( 1, p ), 1, - + A( 1, q ), 1 ) + $ A( 1, p ), 1, + $ A( 1, q ), 1 ) WORK( p ) = WORK( p )*CS WORK( q ) = WORK( q ) / CS IF( RSVEC ) THEN CALL DAXPY( MVL, - + -T*AQOAP, - + V( 1, q ), 1, - + V( 1, p ), 1 ) + $ -T*AQOAP, + $ V( 1, q ), 1, + $ V( 1, p ), 1 ) CALL DAXPY( MVL, - + CS*SN*APOAQ, - + V( 1, p ), 1, - + V( 1, q ), 1 ) + $ CS*SN*APOAQ, + $ V( 1, p ), 1, + $ V( 1, q ), 1 ) END IF ELSE CALL DAXPY( M, T*APOAQ, - + A( 1, p ), 1, - + A( 1, q ), 1 ) + $ A( 1, p ), 1, + $ A( 1, q ), 1 ) CALL DAXPY( M, - + -CS*SN*AQOAP, - + A( 1, q ), 1, - + A( 1, p ), 1 ) + $ -CS*SN*AQOAP, + $ A( 1, q ), 1, + $ A( 1, p ), 1 ) WORK( p ) = WORK( p ) / CS WORK( q ) = WORK( q )*CS IF( RSVEC ) THEN CALL DAXPY( MVL, - + T*APOAQ, V( 1, p ), - + 1, V( 1, q ), 1 ) + $ T*APOAQ, V( 1, p ), + $ 1, V( 1, q ), 1 ) CALL DAXPY( MVL, - + -CS*SN*AQOAP, - + V( 1, q ), 1, - + V( 1, p ), 1 ) + $ -CS*SN*AQOAP, + $ V( 1, q ), 1, + $ V( 1, p ), 1 ) END IF END IF END IF @@ -1268,39 +1266,39 @@ ELSE IF( AAPP.GT.AAQQ ) THEN CALL DCOPY( M, A( 1, p ), 1, - + WORK( N+1 ), 1 ) + $ WORK( N+1 ), 1 ) CALL DLASCL( 'G', 0, 0, AAPP, ONE, - + M, 1, WORK( N+1 ), LDA, - + IERR ) + $ M, 1, WORK( N+1 ), LDA, + $ IERR ) CALL DLASCL( 'G', 0, 0, AAQQ, ONE, - + M, 1, A( 1, q ), LDA, - + IERR ) + $ M, 1, A( 1, q ), LDA, + $ IERR ) TEMP1 = -AAPQ*WORK( p ) / WORK( q ) CALL DAXPY( M, TEMP1, WORK( N+1 ), - + 1, A( 1, q ), 1 ) + $ 1, A( 1, q ), 1 ) CALL DLASCL( 'G', 0, 0, ONE, AAQQ, - + M, 1, A( 1, q ), LDA, - + IERR ) + $ M, 1, A( 1, q ), LDA, + $ IERR ) SVA( q ) = AAQQ*DSQRT( DMAX1( ZERO, - + ONE-AAPQ*AAPQ ) ) + $ ONE-AAPQ*AAPQ ) ) MXSINJ = DMAX1( MXSINJ, SFMIN ) ELSE CALL DCOPY( M, A( 1, q ), 1, - + WORK( N+1 ), 1 ) + $ WORK( N+1 ), 1 ) CALL DLASCL( 'G', 0, 0, AAQQ, ONE, - + M, 1, WORK( N+1 ), LDA, - + IERR ) + $ M, 1, WORK( N+1 ), LDA, + $ IERR ) CALL DLASCL( 'G', 0, 0, AAPP, ONE, - + M, 1, A( 1, p ), LDA, - + IERR ) + $ M, 1, A( 1, p ), LDA, + $ IERR ) TEMP1 = -AAPQ*WORK( q ) / WORK( p ) CALL DAXPY( M, TEMP1, WORK( N+1 ), - + 1, A( 1, p ), 1 ) + $ 1, A( 1, p ), 1 ) CALL DLASCL( 'G', 0, 0, ONE, AAPP, - + M, 1, A( 1, p ), LDA, - + IERR ) + $ M, 1, A( 1, p ), LDA, + $ IERR ) SVA( p ) = AAPP*DSQRT( DMAX1( ZERO, - + ONE-AAPQ*AAPQ ) ) + $ ONE-AAPQ*AAPQ ) ) MXSINJ = DMAX1( MXSINJ, SFMIN ) END IF END IF @@ -1309,29 +1307,29 @@ * In the case of cancellation in updating SVA(q) * .. recompute SVA(q) IF( ( SVA( q ) / AAQQ )**2.LE.ROOTEPS ) - + THEN + $ THEN IF( ( AAQQ.LT.ROOTBIG ) .AND. - + ( AAQQ.GT.ROOTSFMIN ) ) THEN + $ ( AAQQ.GT.ROOTSFMIN ) ) THEN SVA( q ) = DNRM2( M, A( 1, q ), 1 )* - + WORK( q ) + $ WORK( q ) ELSE T = ZERO AAQQ = ONE CALL DLASSQ( M, A( 1, q ), 1, T, - + AAQQ ) + $ AAQQ ) SVA( q ) = T*DSQRT( AAQQ )*WORK( q ) END IF END IF IF( ( AAPP / AAPP0 )**2.LE.ROOTEPS ) THEN IF( ( AAPP.LT.ROOTBIG ) .AND. - + ( AAPP.GT.ROOTSFMIN ) ) THEN + $ ( AAPP.GT.ROOTSFMIN ) ) THEN AAPP = DNRM2( M, A( 1, p ), 1 )* - + WORK( p ) + $ WORK( p ) ELSE T = ZERO AAPP = ONE CALL DLASSQ( M, A( 1, p ), 1, T, - + AAPP ) + $ AAPP ) AAPP = T*DSQRT( AAPP )*WORK( p ) END IF SVA( p ) = AAPP @@ -1350,13 +1348,13 @@ END IF * IF( ( i.LE.SWBAND ) .AND. ( IJBLSK.GE.BLSKIP ) ) - + THEN + $ THEN SVA( p ) = AAPP NOTROT = 0 GO TO 2011 END IF IF( ( i.LE.SWBAND ) .AND. - + ( PSKIPPED.GT.ROWSKIP ) ) THEN + $ ( PSKIPPED.GT.ROWSKIP ) ) THEN AAPP = -AAPP NOTROT = 0 GO TO 2203 @@ -1371,7 +1369,7 @@ ELSE * IF( AAPP.EQ.ZERO )NOTROT = NOTROT + - + MIN0( jgl+KBL-1, N ) - jgl + 1 + $ MIN0( jgl+KBL-1, N ) - jgl + 1 IF( AAPP.LT.ZERO )NOTROT = 0 * END IF @@ -1391,7 +1389,7 @@ * * .. update SVA(N) IF( ( SVA( N ).LT.ROOTBIG ) .AND. ( SVA( N ).GT.ROOTSFMIN ) ) - + THEN + $ THEN SVA( N ) = DNRM2( M, A( 1, N ), 1 )*WORK( N ) ELSE T = ZERO @@ -1403,10 +1401,10 @@ * Additional steering devices * IF( ( i.LT.SWBAND ) .AND. ( ( MXAAPQ.LE.ROOTTOL ) .OR. - + ( ISWROT.LE.N ) ) )SWBAND = i + $ ( ISWROT.LE.N ) ) )SWBAND = i * IF( ( i.GT.SWBAND+1 ) .AND. ( MXAAPQ.LT.DSQRT( DBLE( N ) )* - + TOL ) .AND. ( DBLE( N )*MXAAPQ*MXSINJ.LT.TOL ) ) THEN + $ TOL ) .AND. ( DBLE( N )*MXAAPQ*MXSINJ.LT.TOL ) ) THEN GO TO 1994 END IF * @@ -1479,8 +1477,8 @@ * * Undo scaling, if necessary (and possible). IF( ( ( SKL.GT.ONE ) .AND. ( SVA( 1 ).LT.( BIG / - + SKL) ) ) .OR. ( ( SKL.LT.ONE ) .AND. ( SVA( N2 ).GT. - + ( SFMIN / SKL) ) ) ) THEN + $ SKL) ) ) .OR. ( ( SKL.LT.ONE ) .AND. ( SVA( N2 ).GT. + $ ( SFMIN / SKL) ) ) ) THEN DO 2400 p = 1, N SVA( p ) = SKL*SVA( p ) 2400 CONTINUE