1 SUBROUTINE SLASQ6( I0, N0, Z, PP, DMIN, DMIN1, DMIN2, DN, 2 $ DNM1, DNM2 ) 3* 4* -- LAPACK auxiliary routine (version 3.0) -- 5* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., 6* Courant Institute, Argonne National Lab, and Rice University 7* October 31, 1999 8* 9* .. Scalar Arguments .. 10 INTEGER I0, N0, PP 11 REAL DMIN, DMIN1, DMIN2, DN, DNM1, DNM2 12* .. 13* .. Array Arguments .. 14 REAL Z( * ) 15* .. 16* 17* Purpose 18* ======= 19* 20* SLASQ6 computes one dqd (shift equal to zero) transform in 21* ping-pong form, with protection against underflow and overflow. 22* 23* Arguments 24* ========= 25* 26* I0 (input) INTEGER 27* First index. 28* 29* N0 (input) INTEGER 30* Last index. 31* 32* Z (input) REAL array, dimension ( 4*N ) 33* Z holds the qd array. EMIN is stored in Z(4*N0) to avoid 34* an extra argument. 35* 36* PP (input) INTEGER 37* PP=0 for ping, PP=1 for pong. 38* 39* DMIN (output) REAL 40* Minimum value of d. 41* 42* DMIN1 (output) REAL 43* Minimum value of d, excluding D( N0 ). 44* 45* DMIN2 (output) REAL 46* Minimum value of d, excluding D( N0 ) and D( N0-1 ). 47* 48* DN (output) REAL 49* d(N0), the last value of d. 50* 51* DNM1 (output) REAL 52* d(N0-1). 53* 54* DNM2 (output) REAL 55* d(N0-2). 56* 57* ===================================================================== 58* 59* .. Parameter .. 60 REAL ZERO 61 PARAMETER ( ZERO = 0.0E0 ) 62* .. 63* .. Local Scalars .. 64 INTEGER J4, J4P2 65 REAL D, EMIN, SAFMIN, TEMP 66* .. 67* .. External Function .. 68 REAL SLAMCH 69 EXTERNAL SLAMCH 70* .. 71* .. Intrinsic Functions .. 72 INTRINSIC MIN 73* .. 74* .. Executable Statements .. 75* 76 IF( ( N0-I0-1 ).LE.0 ) 77 $ RETURN 78* 79 SAFMIN = SLAMCH( 'Safe minimum' ) 80 J4 = 4*I0 + PP - 3 81 EMIN = Z( J4+4 ) 82 D = Z( J4 ) 83 DMIN = D 84* 85 IF( PP.EQ.0 ) THEN 86 DO 10 J4 = 4*I0, 4*( N0-3 ), 4 87 Z( J4-2 ) = D + Z( J4-1 ) 88 IF( Z( J4-2 ).EQ.ZERO ) THEN 89 Z( J4 ) = ZERO 90 D = Z( J4+1 ) 91 DMIN = D 92 EMIN = ZERO 93 ELSE IF( SAFMIN*Z( J4+1 ).LT.Z( J4-2 ) .AND. 94 $ SAFMIN*Z( J4-2 ).LT.Z( J4+1 ) ) THEN 95 TEMP = Z( J4+1 ) / Z( J4-2 ) 96 Z( J4 ) = Z( J4-1 )*TEMP 97 D = D*TEMP 98 ELSE 99 Z( J4 ) = Z( J4+1 )*( Z( J4-1 ) / Z( J4-2 ) ) 100 D = Z( J4+1 )*( D / Z( J4-2 ) ) 101 END IF 102 DMIN = MIN( DMIN, D ) 103 EMIN = MIN( EMIN, Z( J4 ) ) 104 10 CONTINUE 105 ELSE 106 DO 20 J4 = 4*I0, 4*( N0-3 ), 4 107 Z( J4-3 ) = D + Z( J4 ) 108 IF( Z( J4-3 ).EQ.ZERO ) THEN 109 Z( J4-1 ) = ZERO 110 D = Z( J4+2 ) 111 DMIN = D 112 EMIN = ZERO 113 ELSE IF( SAFMIN*Z( J4+2 ).LT.Z( J4-3 ) .AND. 114 $ SAFMIN*Z( J4-3 ).LT.Z( J4+2 ) ) THEN 115 TEMP = Z( J4+2 ) / Z( J4-3 ) 116 Z( J4-1 ) = Z( J4 )*TEMP 117 D = D*TEMP 118 ELSE 119 Z( J4-1 ) = Z( J4+2 )*( Z( J4 ) / Z( J4-3 ) ) 120 D = Z( J4+2 )*( D / Z( J4-3 ) ) 121 END IF 122 DMIN = MIN( DMIN, D ) 123 EMIN = MIN( EMIN, Z( J4-1 ) ) 124 20 CONTINUE 125 END IF 126* 127* Unroll last two steps. 128* 129 DNM2 = D 130 DMIN2 = DMIN 131 J4 = 4*( N0-2 ) - PP 132 J4P2 = J4 + 2*PP - 1 133 Z( J4-2 ) = DNM2 + Z( J4P2 ) 134 IF( Z( J4-2 ).EQ.ZERO ) THEN 135 Z( J4 ) = ZERO 136 DNM1 = Z( J4P2+2 ) 137 DMIN = DNM1 138 EMIN = ZERO 139 ELSE IF( SAFMIN*Z( J4P2+2 ).LT.Z( J4-2 ) .AND. 140 $ SAFMIN*Z( J4-2 ).LT.Z( J4P2+2 ) ) THEN 141 TEMP = Z( J4P2+2 ) / Z( J4-2 ) 142 Z( J4 ) = Z( J4P2 )*TEMP 143 DNM1 = DNM2*TEMP 144 ELSE 145 Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) ) 146 DNM1 = Z( J4P2+2 )*( DNM2 / Z( J4-2 ) ) 147 END IF 148 DMIN = MIN( DMIN, DNM1 ) 149* 150 DMIN1 = DMIN 151 J4 = J4 + 4 152 J4P2 = J4 + 2*PP - 1 153 Z( J4-2 ) = DNM1 + Z( J4P2 ) 154 IF( Z( J4-2 ).EQ.ZERO ) THEN 155 Z( J4 ) = ZERO 156 DN = Z( J4P2+2 ) 157 DMIN = DN 158 EMIN = ZERO 159 ELSE IF( SAFMIN*Z( J4P2+2 ).LT.Z( J4-2 ) .AND. 160 $ SAFMIN*Z( J4-2 ).LT.Z( J4P2+2 ) ) THEN 161 TEMP = Z( J4P2+2 ) / Z( J4-2 ) 162 Z( J4 ) = Z( J4P2 )*TEMP 163 DN = DNM1*TEMP 164 ELSE 165 Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) ) 166 DN = Z( J4P2+2 )*( DNM1 / Z( J4-2 ) ) 167 END IF 168 DMIN = MIN( DMIN, DN ) 169* 170 Z( J4+2 ) = DN 171 Z( 4*N0-PP ) = EMIN 172 RETURN 173* 174* End of SLASQ6 175* 176 END 177