![]() ![]() | ![]() |
Mise à jour de lapack vers la version 3.3.0.
1: SUBROUTINE DLARFG( N, ALPHA, X, INCX, TAU ) 2: * 3: * -- LAPACK auxiliary routine (version 3.2) -- 4: * -- LAPACK is a software package provided by Univ. of Tennessee, -- 5: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 6: * November 2006 7: * 8: * .. Scalar Arguments .. 9: INTEGER INCX, N 10: DOUBLE PRECISION ALPHA, TAU 11: * .. 12: * .. Array Arguments .. 13: DOUBLE PRECISION X( * ) 14: * .. 15: * 16: * Purpose 17: * ======= 18: * 19: * DLARFG generates a real elementary reflector H of order n, such 20: * that 21: * 22: * H * ( alpha ) = ( beta ), H' * H = I. 23: * ( x ) ( 0 ) 24: * 25: * where alpha and beta are scalars, and x is an (n-1)-element real 26: * vector. H is represented in the form 27: * 28: * H = I - tau * ( 1 ) * ( 1 v' ) , 29: * ( v ) 30: * 31: * where tau is a real scalar and v is a real (n-1)-element 32: * vector. 33: * 34: * If the elements of x are all zero, then tau = 0 and H is taken to be 35: * the unit matrix. 36: * 37: * Otherwise 1 <= tau <= 2. 38: * 39: * Arguments 40: * ========= 41: * 42: * N (input) INTEGER 43: * The order of the elementary reflector. 44: * 45: * ALPHA (input/output) DOUBLE PRECISION 46: * On entry, the value alpha. 47: * On exit, it is overwritten with the value beta. 48: * 49: * X (input/output) DOUBLE PRECISION array, dimension 50: * (1+(N-2)*abs(INCX)) 51: * On entry, the vector x. 52: * On exit, it is overwritten with the vector v. 53: * 54: * INCX (input) INTEGER 55: * The increment between elements of X. INCX > 0. 56: * 57: * TAU (output) DOUBLE PRECISION 58: * The value tau. 59: * 60: * ===================================================================== 61: * 62: * .. Parameters .. 63: DOUBLE PRECISION ONE, ZERO 64: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 65: * .. 66: * .. Local Scalars .. 67: INTEGER J, KNT 68: DOUBLE PRECISION BETA, RSAFMN, SAFMIN, XNORM 69: * .. 70: * .. External Functions .. 71: DOUBLE PRECISION DLAMCH, DLAPY2, DNRM2 72: EXTERNAL DLAMCH, DLAPY2, DNRM2 73: * .. 74: * .. Intrinsic Functions .. 75: INTRINSIC ABS, SIGN 76: * .. 77: * .. External Subroutines .. 78: EXTERNAL DSCAL 79: * .. 80: * .. Executable Statements .. 81: * 82: IF( N.LE.1 ) THEN 83: TAU = ZERO 84: RETURN 85: END IF 86: * 87: XNORM = DNRM2( N-1, X, INCX ) 88: * 89: IF( XNORM.EQ.ZERO ) THEN 90: * 91: * H = I 92: * 93: TAU = ZERO 94: ELSE 95: * 96: * general case 97: * 98: BETA = -SIGN( DLAPY2( ALPHA, XNORM ), ALPHA ) 99: SAFMIN = DLAMCH( 'S' ) / DLAMCH( 'E' ) 100: KNT = 0 101: IF( ABS( BETA ).LT.SAFMIN ) THEN 102: * 103: * XNORM, BETA may be inaccurate; scale X and recompute them 104: * 105: RSAFMN = ONE / SAFMIN 106: 10 CONTINUE 107: KNT = KNT + 1 108: CALL DSCAL( N-1, RSAFMN, X, INCX ) 109: BETA = BETA*RSAFMN 110: ALPHA = ALPHA*RSAFMN 111: IF( ABS( BETA ).LT.SAFMIN ) 112: $ GO TO 10 113: * 114: * New BETA is at most 1, at least SAFMIN 115: * 116: XNORM = DNRM2( N-1, X, INCX ) 117: BETA = -SIGN( DLAPY2( ALPHA, XNORM ), ALPHA ) 118: END IF 119: TAU = ( BETA-ALPHA ) / BETA 120: CALL DSCAL( N-1, ONE / ( ALPHA-BETA ), X, INCX ) 121: * 122: * If ALPHA is subnormal, it may lose relative accuracy 123: * 124: DO 20 J = 1, KNT 125: BETA = BETA*SAFMIN 126: 20 CONTINUE 127: ALPHA = BETA 128: END IF 129: * 130: RETURN 131: * 132: * End of DLARFG 133: * 134: END