![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b SLASQ6 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download SLASQ6 + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/slasq6.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/slasq6.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/slasq6.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE SLASQ6( I0, N0, Z, PP, DMIN, DMIN1, DMIN2, DN, 00022 * DNM1, DNM2 ) 00023 * 00024 * .. Scalar Arguments .. 00025 * INTEGER I0, N0, PP 00026 * REAL DMIN, DMIN1, DMIN2, DN, DNM1, DNM2 00027 * .. 00028 * .. Array Arguments .. 00029 * REAL Z( * ) 00030 * .. 00031 * 00032 * 00033 *> \par Purpose: 00034 * ============= 00035 *> 00036 *> \verbatim 00037 *> 00038 *> SLASQ6 computes one dqd (shift equal to zero) transform in 00039 *> ping-pong form, with protection against underflow and overflow. 00040 *> \endverbatim 00041 * 00042 * Arguments: 00043 * ========== 00044 * 00045 *> \param[in] I0 00046 *> \verbatim 00047 *> I0 is INTEGER 00048 *> First index. 00049 *> \endverbatim 00050 *> 00051 *> \param[in] N0 00052 *> \verbatim 00053 *> N0 is INTEGER 00054 *> Last index. 00055 *> \endverbatim 00056 *> 00057 *> \param[in] Z 00058 *> \verbatim 00059 *> Z is REAL array, dimension ( 4*N ) 00060 *> Z holds the qd array. EMIN is stored in Z(4*N0) to avoid 00061 *> an extra argument. 00062 *> \endverbatim 00063 *> 00064 *> \param[in] PP 00065 *> \verbatim 00066 *> PP is INTEGER 00067 *> PP=0 for ping, PP=1 for pong. 00068 *> \endverbatim 00069 *> 00070 *> \param[out] DMIN 00071 *> \verbatim 00072 *> DMIN is REAL 00073 *> Minimum value of d. 00074 *> \endverbatim 00075 *> 00076 *> \param[out] DMIN1 00077 *> \verbatim 00078 *> DMIN1 is REAL 00079 *> Minimum value of d, excluding D( N0 ). 00080 *> \endverbatim 00081 *> 00082 *> \param[out] DMIN2 00083 *> \verbatim 00084 *> DMIN2 is REAL 00085 *> Minimum value of d, excluding D( N0 ) and D( N0-1 ). 00086 *> \endverbatim 00087 *> 00088 *> \param[out] DN 00089 *> \verbatim 00090 *> DN is REAL 00091 *> d(N0), the last value of d. 00092 *> \endverbatim 00093 *> 00094 *> \param[out] DNM1 00095 *> \verbatim 00096 *> DNM1 is REAL 00097 *> d(N0-1). 00098 *> \endverbatim 00099 *> 00100 *> \param[out] DNM2 00101 *> \verbatim 00102 *> DNM2 is REAL 00103 *> d(N0-2). 00104 *> \endverbatim 00105 * 00106 * Authors: 00107 * ======== 00108 * 00109 *> \author Univ. of Tennessee 00110 *> \author Univ. of California Berkeley 00111 *> \author Univ. of Colorado Denver 00112 *> \author NAG Ltd. 00113 * 00114 *> \date November 2011 00115 * 00116 *> \ingroup auxOTHERcomputational 00117 * 00118 * ===================================================================== 00119 SUBROUTINE SLASQ6( I0, N0, Z, PP, DMIN, DMIN1, DMIN2, DN, 00120 $ DNM1, DNM2 ) 00121 * 00122 * -- LAPACK computational routine (version 3.4.0) -- 00123 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00124 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00125 * November 2011 00126 * 00127 * .. Scalar Arguments .. 00128 INTEGER I0, N0, PP 00129 REAL DMIN, DMIN1, DMIN2, DN, DNM1, DNM2 00130 * .. 00131 * .. Array Arguments .. 00132 REAL Z( * ) 00133 * .. 00134 * 00135 * ===================================================================== 00136 * 00137 * .. Parameter .. 00138 REAL ZERO 00139 PARAMETER ( ZERO = 0.0E0 ) 00140 * .. 00141 * .. Local Scalars .. 00142 INTEGER J4, J4P2 00143 REAL D, EMIN, SAFMIN, TEMP 00144 * .. 00145 * .. External Function .. 00146 REAL SLAMCH 00147 EXTERNAL SLAMCH 00148 * .. 00149 * .. Intrinsic Functions .. 00150 INTRINSIC MIN 00151 * .. 00152 * .. Executable Statements .. 00153 * 00154 IF( ( N0-I0-1 ).LE.0 ) 00155 $ RETURN 00156 * 00157 SAFMIN = SLAMCH( 'Safe minimum' ) 00158 J4 = 4*I0 + PP - 3 00159 EMIN = Z( J4+4 ) 00160 D = Z( J4 ) 00161 DMIN = D 00162 * 00163 IF( PP.EQ.0 ) THEN 00164 DO 10 J4 = 4*I0, 4*( N0-3 ), 4 00165 Z( J4-2 ) = D + Z( J4-1 ) 00166 IF( Z( J4-2 ).EQ.ZERO ) THEN 00167 Z( J4 ) = ZERO 00168 D = Z( J4+1 ) 00169 DMIN = D 00170 EMIN = ZERO 00171 ELSE IF( SAFMIN*Z( J4+1 ).LT.Z( J4-2 ) .AND. 00172 $ SAFMIN*Z( J4-2 ).LT.Z( J4+1 ) ) THEN 00173 TEMP = Z( J4+1 ) / Z( J4-2 ) 00174 Z( J4 ) = Z( J4-1 )*TEMP 00175 D = D*TEMP 00176 ELSE 00177 Z( J4 ) = Z( J4+1 )*( Z( J4-1 ) / Z( J4-2 ) ) 00178 D = Z( J4+1 )*( D / Z( J4-2 ) ) 00179 END IF 00180 DMIN = MIN( DMIN, D ) 00181 EMIN = MIN( EMIN, Z( J4 ) ) 00182 10 CONTINUE 00183 ELSE 00184 DO 20 J4 = 4*I0, 4*( N0-3 ), 4 00185 Z( J4-3 ) = D + Z( J4 ) 00186 IF( Z( J4-3 ).EQ.ZERO ) THEN 00187 Z( J4-1 ) = ZERO 00188 D = Z( J4+2 ) 00189 DMIN = D 00190 EMIN = ZERO 00191 ELSE IF( SAFMIN*Z( J4+2 ).LT.Z( J4-3 ) .AND. 00192 $ SAFMIN*Z( J4-3 ).LT.Z( J4+2 ) ) THEN 00193 TEMP = Z( J4+2 ) / Z( J4-3 ) 00194 Z( J4-1 ) = Z( J4 )*TEMP 00195 D = D*TEMP 00196 ELSE 00197 Z( J4-1 ) = Z( J4+2 )*( Z( J4 ) / Z( J4-3 ) ) 00198 D = Z( J4+2 )*( D / Z( J4-3 ) ) 00199 END IF 00200 DMIN = MIN( DMIN, D ) 00201 EMIN = MIN( EMIN, Z( J4-1 ) ) 00202 20 CONTINUE 00203 END IF 00204 * 00205 * Unroll last two steps. 00206 * 00207 DNM2 = D 00208 DMIN2 = DMIN 00209 J4 = 4*( N0-2 ) - PP 00210 J4P2 = J4 + 2*PP - 1 00211 Z( J4-2 ) = DNM2 + Z( J4P2 ) 00212 IF( Z( J4-2 ).EQ.ZERO ) THEN 00213 Z( J4 ) = ZERO 00214 DNM1 = Z( J4P2+2 ) 00215 DMIN = DNM1 00216 EMIN = ZERO 00217 ELSE IF( SAFMIN*Z( J4P2+2 ).LT.Z( J4-2 ) .AND. 00218 $ SAFMIN*Z( J4-2 ).LT.Z( J4P2+2 ) ) THEN 00219 TEMP = Z( J4P2+2 ) / Z( J4-2 ) 00220 Z( J4 ) = Z( J4P2 )*TEMP 00221 DNM1 = DNM2*TEMP 00222 ELSE 00223 Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) ) 00224 DNM1 = Z( J4P2+2 )*( DNM2 / Z( J4-2 ) ) 00225 END IF 00226 DMIN = MIN( DMIN, DNM1 ) 00227 * 00228 DMIN1 = DMIN 00229 J4 = J4 + 4 00230 J4P2 = J4 + 2*PP - 1 00231 Z( J4-2 ) = DNM1 + Z( J4P2 ) 00232 IF( Z( J4-2 ).EQ.ZERO ) THEN 00233 Z( J4 ) = ZERO 00234 DN = Z( J4P2+2 ) 00235 DMIN = DN 00236 EMIN = ZERO 00237 ELSE IF( SAFMIN*Z( J4P2+2 ).LT.Z( J4-2 ) .AND. 00238 $ SAFMIN*Z( J4-2 ).LT.Z( J4P2+2 ) ) THEN 00239 TEMP = Z( J4P2+2 ) / Z( J4-2 ) 00240 Z( J4 ) = Z( J4P2 )*TEMP 00241 DN = DNM1*TEMP 00242 ELSE 00243 Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) ) 00244 DN = Z( J4P2+2 )*( DNM1 / Z( J4-2 ) ) 00245 END IF 00246 DMIN = MIN( DMIN, DN ) 00247 * 00248 Z( J4+2 ) = DN 00249 Z( 4*N0-PP ) = EMIN 00250 RETURN 00251 * 00252 * End of SLASQ6 00253 * 00254 END