--- rpl/lapack/lapack/dgsvj1.f 2010/12/21 13:53:27 1.5 +++ rpl/lapack/lapack/dgsvj1.f 2011/07/22 07:38:05 1.6 @@ -1,11 +1,11 @@ SUBROUTINE DGSVJ1( JOBV, M, N, N1, A, LDA, D, SVA, MV, V, LDV, - + EPS, SFMIN, TOL, NSWEEP, WORK, LWORK, INFO ) + $ EPS, SFMIN, TOL, NSWEEP, 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..-- @@ -24,7 +24,7 @@ * .. * .. Array Arguments .. DOUBLE PRECISION A( LDA, * ), D( N ), SVA( N ), V( LDV, * ), - + WORK( LWORK ) + $ WORK( LWORK ) * .. * * Purpose @@ -162,16 +162,16 @@ * .. Local Parameters .. DOUBLE PRECISION ZERO, HALF, ONE, TWO PARAMETER ( ZERO = 0.0D0, HALF = 0.5D0, ONE = 1.0D0, - + TWO = 2.0D0 ) + $ TWO = 2.0D0 ) * .. * .. Local Scalars .. DOUBLE PRECISION AAPP, AAPP0, AAPQ, AAQQ, APOAQ, AQOAP, BIG, - + BIGTHETA, CS, LARGE, MXAAPQ, MXSINJ, ROOTBIG, - + ROOTEPS, ROOTSFMIN, ROOTTOL, SMALL, SN, T, - + TEMP1, THETA, THSIGN + $ BIGTHETA, CS, LARGE, MXAAPQ, MXSINJ, ROOTBIG, + $ ROOTEPS, ROOTSFMIN, ROOTTOL, SMALL, SN, T, + $ TEMP1, THETA, THSIGN INTEGER BLSKIP, EMPTSW, i, ibr, igl, IERR, IJBLSK, - + ISWROT, jbc, jgl, KBL, MVL, NOTROT, nblc, nblr, - + p, PSKIPPED, q, ROWSKIP, SWBAND + $ ISWROT, jbc, jgl, KBL, MVL, NOTROT, nblc, nblr, + $ p, PSKIPPED, q, ROWSKIP, SWBAND LOGICAL APPLV, ROTOK, RSVEC * .. * .. Local Arrays .. @@ -208,7 +208,7 @@ ELSE IF( ( RSVEC.OR.APPLV ) .AND. ( 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( TOL.LE.EPS ) THEN INFO = -14 @@ -333,14 +333,14 @@ END IF IF( AAPP.LT.( BIG / AAQQ ) ) THEN AAPQ = ( DDOT( M, A( 1, p ), 1, A( 1, - + q ), 1 )*D( p )*D( q ) / AAQQ ) - + / AAPP + $ q ), 1 )*D( p )*D( q ) / AAQQ ) + $ / AAPP ELSE CALL DCOPY( M, A( 1, p ), 1, WORK, 1 ) CALL DLASCL( 'G', 0, 0, AAPP, D( p ), - + M, 1, WORK, LDA, IERR ) + $ M, 1, WORK, LDA, IERR ) AAPQ = DDOT( M, WORK, 1, A( 1, q ), - + 1 )*D( q ) / AAQQ + $ 1 )*D( q ) / AAQQ END IF ELSE IF( AAPP.GE.AAQQ ) THEN @@ -350,14 +350,14 @@ END IF IF( AAPP.GT.( SMALL / AAQQ ) ) THEN AAPQ = ( DDOT( M, A( 1, p ), 1, A( 1, - + q ), 1 )*D( p )*D( q ) / AAQQ ) - + / AAPP + $ q ), 1 )*D( p )*D( q ) / AAQQ ) + $ / AAPP ELSE CALL DCOPY( M, A( 1, q ), 1, WORK, 1 ) CALL DLASCL( 'G', 0, 0, AAQQ, D( q ), - + M, 1, WORK, LDA, IERR ) + $ M, 1, WORK, LDA, IERR ) AAPQ = DDOT( M, WORK, 1, A( 1, p ), - + 1 )*D( p ) / AAPP + $ 1 )*D( p ) / AAPP END IF END IF @@ -375,8 +375,7 @@ * 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 @@ -384,15 +383,15 @@ FASTR( 3 ) = T*D( p ) / D( q ) FASTR( 4 ) = -T*D( q ) / D( 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 * @@ -401,14 +400,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 = D( p ) / D( q ) AQOAP = D( q ) / D( p ) @@ -420,26 +419,26 @@ D( p ) = D( p )*CS D( q ) = D( 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 D( p ) = D( p )*CS D( q ) = D( q ) / CS @@ -447,60 +446,60 @@ ELSE IF( D( 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 D( p ) = D( p ) / CS D( q ) = D( q )*CS ELSE IF( D( p ).GE.D( q ) ) 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 ) D( p ) = D( p )*CS D( q ) = D( 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 ) D( p ) = D( p ) / CS D( q ) = D( 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 @@ -510,37 +509,37 @@ ELSE IF( AAPP.GT.AAQQ ) THEN CALL DCOPY( M, A( 1, p ), 1, WORK, - + 1 ) + $ 1 ) CALL DLASCL( 'G', 0, 0, AAPP, ONE, - + M, 1, WORK, LDA, IERR ) + $ M, 1, WORK, 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*D( p ) / D( q ) CALL DAXPY( M, TEMP1, WORK, 1, - + A( 1, q ), 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, - + 1 ) + $ 1 ) CALL DLASCL( 'G', 0, 0, AAQQ, ONE, - + M, 1, WORK, LDA, IERR ) + $ M, 1, WORK, 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*D( q ) / D( p ) CALL DAXPY( M, TEMP1, WORK, 1, - + A( 1, p ), 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 @@ -549,29 +548,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 )* - + D( q ) + $ D( q ) ELSE T = ZERO AAQQ = ONE CALL DLASSQ( M, A( 1, q ), 1, T, - + AAQQ ) + $ AAQQ ) SVA( q ) = T*DSQRT( AAQQ )*D( 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 )* - + D( p ) + $ D( p ) ELSE T = ZERO AAPP = ONE CALL DLASSQ( M, A( 1, p ), 1, T, - + AAPP ) + $ AAPP ) AAPP = T*DSQRT( AAPP )*D( p ) END IF SVA( p ) = AAPP @@ -591,13 +590,13 @@ * IF ( NOTROT .GE. EMPTSW ) GO TO 2011 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 @@ -612,7 +611,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 *** IF ( NOTROT .GE. EMPTSW ) GO TO 2011 END IF @@ -632,7 +631,7 @@ * * .. update SVA(N) IF( ( SVA( N ).LT.ROOTBIG ) .AND. ( SVA( N ).GT.ROOTSFMIN ) ) - + THEN + $ THEN SVA( N ) = DNRM2( M, A( 1, N ), 1 )*D( N ) ELSE T = ZERO @@ -644,10 +643,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.DBLE( N )*TOL ) .AND. - + ( DBLE( N )*MXAAPQ*MXSINJ.LT.TOL ) ) THEN + $ ( DBLE( N )*MXAAPQ*MXSINJ.LT.TOL ) ) THEN GO TO 1994 END IF