![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DLARFG 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download DLARFG + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlarfg.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlarfg.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlarfg.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE DLARFG( N, ALPHA, X, INCX, TAU ) 00022 * 00023 * .. Scalar Arguments .. 00024 * INTEGER INCX, N 00025 * DOUBLE PRECISION ALPHA, TAU 00026 * .. 00027 * .. Array Arguments .. 00028 * DOUBLE PRECISION X( * ) 00029 * .. 00030 * 00031 * 00032 *> \par Purpose: 00033 * ============= 00034 *> 00035 *> \verbatim 00036 *> 00037 *> DLARFG generates a real elementary reflector H of order n, such 00038 *> that 00039 *> 00040 *> H * ( alpha ) = ( beta ), H**T * H = I. 00041 *> ( x ) ( 0 ) 00042 *> 00043 *> where alpha and beta are scalars, and x is an (n-1)-element real 00044 *> vector. H is represented in the form 00045 *> 00046 *> H = I - tau * ( 1 ) * ( 1 v**T ) , 00047 *> ( v ) 00048 *> 00049 *> where tau is a real scalar and v is a real (n-1)-element 00050 *> vector. 00051 *> 00052 *> If the elements of x are all zero, then tau = 0 and H is taken to be 00053 *> the unit matrix. 00054 *> 00055 *> Otherwise 1 <= tau <= 2. 00056 *> \endverbatim 00057 * 00058 * Arguments: 00059 * ========== 00060 * 00061 *> \param[in] N 00062 *> \verbatim 00063 *> N is INTEGER 00064 *> The order of the elementary reflector. 00065 *> \endverbatim 00066 *> 00067 *> \param[in,out] ALPHA 00068 *> \verbatim 00069 *> ALPHA is DOUBLE PRECISION 00070 *> On entry, the value alpha. 00071 *> On exit, it is overwritten with the value beta. 00072 *> \endverbatim 00073 *> 00074 *> \param[in,out] X 00075 *> \verbatim 00076 *> X is DOUBLE PRECISION array, dimension 00077 *> (1+(N-2)*abs(INCX)) 00078 *> On entry, the vector x. 00079 *> On exit, it is overwritten with the vector v. 00080 *> \endverbatim 00081 *> 00082 *> \param[in] INCX 00083 *> \verbatim 00084 *> INCX is INTEGER 00085 *> The increment between elements of X. INCX > 0. 00086 *> \endverbatim 00087 *> 00088 *> \param[out] TAU 00089 *> \verbatim 00090 *> TAU is DOUBLE PRECISION 00091 *> The value tau. 00092 *> \endverbatim 00093 * 00094 * Authors: 00095 * ======== 00096 * 00097 *> \author Univ. of Tennessee 00098 *> \author Univ. of California Berkeley 00099 *> \author Univ. of Colorado Denver 00100 *> \author NAG Ltd. 00101 * 00102 *> \date November 2011 00103 * 00104 *> \ingroup doubleOTHERauxiliary 00105 * 00106 * ===================================================================== 00107 SUBROUTINE DLARFG( N, ALPHA, X, INCX, TAU ) 00108 * 00109 * -- LAPACK auxiliary routine (version 3.4.0) -- 00110 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00111 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00112 * November 2011 00113 * 00114 * .. Scalar Arguments .. 00115 INTEGER INCX, N 00116 DOUBLE PRECISION ALPHA, TAU 00117 * .. 00118 * .. Array Arguments .. 00119 DOUBLE PRECISION X( * ) 00120 * .. 00121 * 00122 * ===================================================================== 00123 * 00124 * .. Parameters .. 00125 DOUBLE PRECISION ONE, ZERO 00126 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 00127 * .. 00128 * .. Local Scalars .. 00129 INTEGER J, KNT 00130 DOUBLE PRECISION BETA, RSAFMN, SAFMIN, XNORM 00131 * .. 00132 * .. External Functions .. 00133 DOUBLE PRECISION DLAMCH, DLAPY2, DNRM2 00134 EXTERNAL DLAMCH, DLAPY2, DNRM2 00135 * .. 00136 * .. Intrinsic Functions .. 00137 INTRINSIC ABS, SIGN 00138 * .. 00139 * .. External Subroutines .. 00140 EXTERNAL DSCAL 00141 * .. 00142 * .. Executable Statements .. 00143 * 00144 IF( N.LE.1 ) THEN 00145 TAU = ZERO 00146 RETURN 00147 END IF 00148 * 00149 XNORM = DNRM2( N-1, X, INCX ) 00150 * 00151 IF( XNORM.EQ.ZERO ) THEN 00152 * 00153 * H = I 00154 * 00155 TAU = ZERO 00156 ELSE 00157 * 00158 * general case 00159 * 00160 BETA = -SIGN( DLAPY2( ALPHA, XNORM ), ALPHA ) 00161 SAFMIN = DLAMCH( 'S' ) / DLAMCH( 'E' ) 00162 KNT = 0 00163 IF( ABS( BETA ).LT.SAFMIN ) THEN 00164 * 00165 * XNORM, BETA may be inaccurate; scale X and recompute them 00166 * 00167 RSAFMN = ONE / SAFMIN 00168 10 CONTINUE 00169 KNT = KNT + 1 00170 CALL DSCAL( N-1, RSAFMN, X, INCX ) 00171 BETA = BETA*RSAFMN 00172 ALPHA = ALPHA*RSAFMN 00173 IF( ABS( BETA ).LT.SAFMIN ) 00174 $ GO TO 10 00175 * 00176 * New BETA is at most 1, at least SAFMIN 00177 * 00178 XNORM = DNRM2( N-1, X, INCX ) 00179 BETA = -SIGN( DLAPY2( ALPHA, XNORM ), ALPHA ) 00180 END IF 00181 TAU = ( BETA-ALPHA ) / BETA 00182 CALL DSCAL( N-1, ONE / ( ALPHA-BETA ), X, INCX ) 00183 * 00184 * If ALPHA is subnormal, it may lose relative accuracy 00185 * 00186 DO 20 J = 1, KNT 00187 BETA = BETA*SAFMIN 00188 20 CONTINUE 00189 ALPHA = BETA 00190 END IF 00191 * 00192 RETURN 00193 * 00194 * End of DLARFG 00195 * 00196 END