LAPACK 3.3.0
|
00001 SUBROUTINE SLARTGS( X, Y, SIGMA, CS, SN ) 00002 IMPLICIT NONE 00003 * 00004 * -- LAPACK routine (version 3.3.0) -- 00005 * 00006 * -- Contributed by Brian Sutton of the Randolph-Macon College -- 00007 * -- November 2010 00008 * 00009 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00010 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00011 * 00012 * .. Scalar Arguments .. 00013 REAL CS, SIGMA, SN, X, Y 00014 * .. 00015 * 00016 * Purpose 00017 * ======= 00018 * 00019 * SLARTGS generates a plane rotation designed to introduce a bulge in 00020 * Golub-Reinsch-style implicit QR iteration for the bidiagonal SVD 00021 * problem. X and Y are the top-row entries, and SIGMA is the shift. 00022 * The computed CS and SN define a plane rotation satisfying 00023 * 00024 * [ CS SN ] . [ X^2 - SIGMA ] = [ R ], 00025 * [ -SN CS ] [ X * Y ] [ 0 ] 00026 * 00027 * with R nonnegative. If X^2 - SIGMA and X * Y are 0, then the 00028 * rotation is by PI/2. 00029 * 00030 * Arguments 00031 * ========= 00032 * 00033 * X (input) REAL 00034 * The (1,1) entry of an upper bidiagonal matrix. 00035 * 00036 * Y (input) REAL 00037 * The (1,2) entry of an upper bidiagonal matrix. 00038 * 00039 * SIGMA (input) REAL 00040 * The shift. 00041 * 00042 * CS (output) REAL 00043 * The cosine of the rotation. 00044 * 00045 * SN (output) REAL 00046 * The sine of the rotation. 00047 * 00048 * =================================================================== 00049 * 00050 * .. Parameters .. 00051 REAL NEGONE, ONE, ZERO 00052 PARAMETER ( NEGONE = -1.0E0, ONE = 1.0E0, ZERO = 0.0E0 ) 00053 * .. 00054 * .. Local Scalars .. 00055 REAL R, S, THRESH, W, Z 00056 * .. 00057 * .. External Functions .. 00058 REAL SLAMCH 00059 EXTERNAL SLAMCH 00060 * .. Executable Statements .. 00061 * 00062 THRESH = SLAMCH('E') 00063 * 00064 * Compute the first column of B'*B - SIGMA^2*I, up to a scale 00065 * factor. 00066 * 00067 IF( (SIGMA .EQ. ZERO .AND. ABS(X) .LT. THRESH) .OR. 00068 $ (ABS(X) .EQ. SIGMA .AND. Y .EQ. ZERO) ) THEN 00069 Z = ZERO 00070 W = ZERO 00071 ELSE IF( SIGMA .EQ. ZERO ) THEN 00072 IF( X .GE. ZERO ) THEN 00073 Z = X 00074 W = Y 00075 ELSE 00076 Z = -X 00077 W = -Y 00078 END IF 00079 ELSE IF( ABS(X) .LT. THRESH ) THEN 00080 Z = -SIGMA*SIGMA 00081 W = ZERO 00082 ELSE 00083 IF( X .GE. ZERO ) THEN 00084 S = ONE 00085 ELSE 00086 S = NEGONE 00087 END IF 00088 Z = S * (ABS(X)-SIGMA) * (S+SIGMA/X) 00089 W = S * Y 00090 END IF 00091 * 00092 * Generate the rotation. 00093 * CALL SLARTGP( Z, W, CS, SN, R ) might seem more natural; 00094 * reordering the arguments ensures that if Z = 0 then the rotation 00095 * is by PI/2. 00096 * 00097 CALL SLARTGP( W, Z, SN, CS, R ) 00098 * 00099 RETURN 00100 * 00101 * End SLARTGS 00102 * 00103 END 00104