![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DLARZT 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download DLARZT + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlarzt.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlarzt.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlarzt.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE DLARZT( DIRECT, STOREV, N, K, V, LDV, TAU, T, LDT ) 00022 * 00023 * .. Scalar Arguments .. 00024 * CHARACTER DIRECT, STOREV 00025 * INTEGER K, LDT, LDV, N 00026 * .. 00027 * .. Array Arguments .. 00028 * DOUBLE PRECISION T( LDT, * ), TAU( * ), V( LDV, * ) 00029 * .. 00030 * 00031 * 00032 *> \par Purpose: 00033 * ============= 00034 *> 00035 *> \verbatim 00036 *> 00037 *> DLARZT forms the triangular factor T of a real block reflector 00038 *> H of order > n, which is defined as a product of k elementary 00039 *> reflectors. 00040 *> 00041 *> If DIRECT = 'F', H = H(1) H(2) . . . H(k) and T is upper triangular; 00042 *> 00043 *> If DIRECT = 'B', H = H(k) . . . H(2) H(1) and T is lower triangular. 00044 *> 00045 *> If STOREV = 'C', the vector which defines the elementary reflector 00046 *> H(i) is stored in the i-th column of the array V, and 00047 *> 00048 *> H = I - V * T * V**T 00049 *> 00050 *> If STOREV = 'R', the vector which defines the elementary reflector 00051 *> H(i) is stored in the i-th row of the array V, and 00052 *> 00053 *> H = I - V**T * T * V 00054 *> 00055 *> Currently, only STOREV = 'R' and DIRECT = 'B' are supported. 00056 *> \endverbatim 00057 * 00058 * Arguments: 00059 * ========== 00060 * 00061 *> \param[in] DIRECT 00062 *> \verbatim 00063 *> DIRECT is CHARACTER*1 00064 *> Specifies the order in which the elementary reflectors are 00065 *> multiplied to form the block reflector: 00066 *> = 'F': H = H(1) H(2) . . . H(k) (Forward, not supported yet) 00067 *> = 'B': H = H(k) . . . H(2) H(1) (Backward) 00068 *> \endverbatim 00069 *> 00070 *> \param[in] STOREV 00071 *> \verbatim 00072 *> STOREV is CHARACTER*1 00073 *> Specifies how the vectors which define the elementary 00074 *> reflectors are stored (see also Further Details): 00075 *> = 'C': columnwise (not supported yet) 00076 *> = 'R': rowwise 00077 *> \endverbatim 00078 *> 00079 *> \param[in] N 00080 *> \verbatim 00081 *> N is INTEGER 00082 *> The order of the block reflector H. N >= 0. 00083 *> \endverbatim 00084 *> 00085 *> \param[in] K 00086 *> \verbatim 00087 *> K is INTEGER 00088 *> The order of the triangular factor T (= the number of 00089 *> elementary reflectors). K >= 1. 00090 *> \endverbatim 00091 *> 00092 *> \param[in,out] V 00093 *> \verbatim 00094 *> V is DOUBLE PRECISION array, dimension 00095 *> (LDV,K) if STOREV = 'C' 00096 *> (LDV,N) if STOREV = 'R' 00097 *> The matrix V. See further details. 00098 *> \endverbatim 00099 *> 00100 *> \param[in] LDV 00101 *> \verbatim 00102 *> LDV is INTEGER 00103 *> The leading dimension of the array V. 00104 *> If STOREV = 'C', LDV >= max(1,N); if STOREV = 'R', LDV >= K. 00105 *> \endverbatim 00106 *> 00107 *> \param[in] TAU 00108 *> \verbatim 00109 *> TAU is DOUBLE PRECISION array, dimension (K) 00110 *> TAU(i) must contain the scalar factor of the elementary 00111 *> reflector H(i). 00112 *> \endverbatim 00113 *> 00114 *> \param[out] T 00115 *> \verbatim 00116 *> T is DOUBLE PRECISION array, dimension (LDT,K) 00117 *> The k by k triangular factor T of the block reflector. 00118 *> If DIRECT = 'F', T is upper triangular; if DIRECT = 'B', T is 00119 *> lower triangular. The rest of the array is not used. 00120 *> \endverbatim 00121 *> 00122 *> \param[in] LDT 00123 *> \verbatim 00124 *> LDT is INTEGER 00125 *> The leading dimension of the array T. LDT >= K. 00126 *> \endverbatim 00127 * 00128 * Authors: 00129 * ======== 00130 * 00131 *> \author Univ. of Tennessee 00132 *> \author Univ. of California Berkeley 00133 *> \author Univ. of Colorado Denver 00134 *> \author NAG Ltd. 00135 * 00136 *> \date November 2011 00137 * 00138 *> \ingroup doubleOTHERcomputational 00139 * 00140 *> \par Contributors: 00141 * ================== 00142 *> 00143 *> A. Petitet, Computer Science Dept., Univ. of Tenn., Knoxville, USA 00144 * 00145 *> \par Further Details: 00146 * ===================== 00147 *> 00148 *> \verbatim 00149 *> 00150 *> The shape of the matrix V and the storage of the vectors which define 00151 *> the H(i) is best illustrated by the following example with n = 5 and 00152 *> k = 3. The elements equal to 1 are not stored; the corresponding 00153 *> array elements are modified but restored on exit. The rest of the 00154 *> array is not used. 00155 *> 00156 *> DIRECT = 'F' and STOREV = 'C': DIRECT = 'F' and STOREV = 'R': 00157 *> 00158 *> ______V_____ 00159 *> ( v1 v2 v3 ) / \ 00160 *> ( v1 v2 v3 ) ( v1 v1 v1 v1 v1 . . . . 1 ) 00161 *> V = ( v1 v2 v3 ) ( v2 v2 v2 v2 v2 . . . 1 ) 00162 *> ( v1 v2 v3 ) ( v3 v3 v3 v3 v3 . . 1 ) 00163 *> ( v1 v2 v3 ) 00164 *> . . . 00165 *> . . . 00166 *> 1 . . 00167 *> 1 . 00168 *> 1 00169 *> 00170 *> DIRECT = 'B' and STOREV = 'C': DIRECT = 'B' and STOREV = 'R': 00171 *> 00172 *> ______V_____ 00173 *> 1 / \ 00174 *> . 1 ( 1 . . . . v1 v1 v1 v1 v1 ) 00175 *> . . 1 ( . 1 . . . v2 v2 v2 v2 v2 ) 00176 *> . . . ( . . 1 . . v3 v3 v3 v3 v3 ) 00177 *> . . . 00178 *> ( v1 v2 v3 ) 00179 *> ( v1 v2 v3 ) 00180 *> V = ( v1 v2 v3 ) 00181 *> ( v1 v2 v3 ) 00182 *> ( v1 v2 v3 ) 00183 *> \endverbatim 00184 *> 00185 * ===================================================================== 00186 SUBROUTINE DLARZT( DIRECT, STOREV, N, K, V, LDV, TAU, T, LDT ) 00187 * 00188 * -- LAPACK computational routine (version 3.4.0) -- 00189 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00190 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00191 * November 2011 00192 * 00193 * .. Scalar Arguments .. 00194 CHARACTER DIRECT, STOREV 00195 INTEGER K, LDT, LDV, N 00196 * .. 00197 * .. Array Arguments .. 00198 DOUBLE PRECISION T( LDT, * ), TAU( * ), V( LDV, * ) 00199 * .. 00200 * 00201 * ===================================================================== 00202 * 00203 * .. Parameters .. 00204 DOUBLE PRECISION ZERO 00205 PARAMETER ( ZERO = 0.0D+0 ) 00206 * .. 00207 * .. Local Scalars .. 00208 INTEGER I, INFO, J 00209 * .. 00210 * .. External Subroutines .. 00211 EXTERNAL DGEMV, DTRMV, XERBLA 00212 * .. 00213 * .. External Functions .. 00214 LOGICAL LSAME 00215 EXTERNAL LSAME 00216 * .. 00217 * .. Executable Statements .. 00218 * 00219 * Check for currently supported options 00220 * 00221 INFO = 0 00222 IF( .NOT.LSAME( DIRECT, 'B' ) ) THEN 00223 INFO = -1 00224 ELSE IF( .NOT.LSAME( STOREV, 'R' ) ) THEN 00225 INFO = -2 00226 END IF 00227 IF( INFO.NE.0 ) THEN 00228 CALL XERBLA( 'DLARZT', -INFO ) 00229 RETURN 00230 END IF 00231 * 00232 DO 20 I = K, 1, -1 00233 IF( TAU( I ).EQ.ZERO ) THEN 00234 * 00235 * H(i) = I 00236 * 00237 DO 10 J = I, K 00238 T( J, I ) = ZERO 00239 10 CONTINUE 00240 ELSE 00241 * 00242 * general case 00243 * 00244 IF( I.LT.K ) THEN 00245 * 00246 * T(i+1:k,i) = - tau(i) * V(i+1:k,1:n) * V(i,1:n)**T 00247 * 00248 CALL DGEMV( 'No transpose', K-I, N, -TAU( I ), 00249 $ V( I+1, 1 ), LDV, V( I, 1 ), LDV, ZERO, 00250 $ T( I+1, I ), 1 ) 00251 * 00252 * T(i+1:k,i) = T(i+1:k,i+1:k) * T(i+1:k,i) 00253 * 00254 CALL DTRMV( 'Lower', 'No transpose', 'Non-unit', K-I, 00255 $ T( I+1, I+1 ), LDT, T( I+1, I ), 1 ) 00256 END IF 00257 T( I, I ) = TAU( I ) 00258 END IF 00259 20 CONTINUE 00260 RETURN 00261 * 00262 * End of DLARZT 00263 * 00264 END