![]() ![]() | ![]() |
Mise à jour de lapack vers la version 3.3.0.
1: SUBROUTINE DLARTGS( X, Y, SIGMA, CS, SN ) 2: IMPLICIT NONE 3: * 4: * -- LAPACK routine (version 3.3.0) -- 5: * 6: * -- Contributed by Brian Sutton of the Randolph-Macon College -- 7: * -- November 2010 8: * 9: * -- LAPACK is a software package provided by Univ. of Tennessee, -- 10: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 11: * 12: * .. Scalar Arguments .. 13: DOUBLE PRECISION CS, SIGMA, SN, X, Y 14: * .. 15: * 16: * Purpose 17: * ======= 18: * 19: * DLARTGS generates a plane rotation designed to introduce a bulge in 20: * Golub-Reinsch-style implicit QR iteration for the bidiagonal SVD 21: * problem. X and Y are the top-row entries, and SIGMA is the shift. 22: * The computed CS and SN define a plane rotation satisfying 23: * 24: * [ CS SN ] . [ X^2 - SIGMA ] = [ R ], 25: * [ -SN CS ] [ X * Y ] [ 0 ] 26: * 27: * with R nonnegative. If X^2 - SIGMA and X * Y are 0, then the 28: * rotation is by PI/2. 29: * 30: * Arguments 31: * ========= 32: * 33: * X (input) DOUBLE PRECISION 34: * The (1,1) entry of an upper bidiagonal matrix. 35: * 36: * Y (input) DOUBLE PRECISION 37: * The (1,2) entry of an upper bidiagonal matrix. 38: * 39: * SIGMA (input) DOUBLE PRECISION 40: * The shift. 41: * 42: * CS (output) DOUBLE PRECISION 43: * The cosine of the rotation. 44: * 45: * SN (output) DOUBLE PRECISION 46: * The sine of the rotation. 47: * 48: * =================================================================== 49: * 50: * .. Parameters .. 51: DOUBLE PRECISION NEGONE, ONE, ZERO 52: PARAMETER ( NEGONE = -1.0D0, ONE = 1.0D0, ZERO = 0.0D0 ) 53: * .. 54: * .. Local Scalars .. 55: DOUBLE PRECISION R, S, THRESH, W, Z 56: * .. 57: * .. External Functions .. 58: DOUBLE PRECISION DLAMCH 59: EXTERNAL DLAMCH 60: * .. Executable Statements .. 61: * 62: THRESH = DLAMCH('E') 63: * 64: * Compute the first column of B'*B - SIGMA^2*I, up to a scale 65: * factor. 66: * 67: IF( (SIGMA .EQ. ZERO .AND. ABS(X) .LT. THRESH) .OR. 68: $ (ABS(X) .EQ. SIGMA .AND. Y .EQ. ZERO) ) THEN 69: Z = ZERO 70: W = ZERO 71: ELSE IF( SIGMA .EQ. ZERO ) THEN 72: IF( X .GE. ZERO ) THEN 73: Z = X 74: W = Y 75: ELSE 76: Z = -X 77: W = -Y 78: END IF 79: ELSE IF( ABS(X) .LT. THRESH ) THEN 80: Z = -SIGMA*SIGMA 81: W = ZERO 82: ELSE 83: IF( X .GE. ZERO ) THEN 84: S = ONE 85: ELSE 86: S = NEGONE 87: END IF 88: Z = S * (ABS(X)-SIGMA) * (S+SIGMA/X) 89: W = S * Y 90: END IF 91: * 92: * Generate the rotation. 93: * CALL DLARTGP( Z, W, CS, SN, R ) might seem more natural; 94: * reordering the arguments ensures that if Z = 0 then the rotation 95: * is by PI/2. 96: * 97: CALL DLARTGP( W, Z, SN, CS, R ) 98: * 99: RETURN 100: * 101: * End DLARTGS 102: * 103: END 104: