![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DLAQP2 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download DLAQP2 + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlaqp2.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlaqp2.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlaqp2.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE DLAQP2( M, N, OFFSET, A, LDA, JPVT, TAU, VN1, VN2, 00022 * WORK ) 00023 * 00024 * .. Scalar Arguments .. 00025 * INTEGER LDA, M, N, OFFSET 00026 * .. 00027 * .. Array Arguments .. 00028 * INTEGER JPVT( * ) 00029 * DOUBLE PRECISION A( LDA, * ), TAU( * ), VN1( * ), VN2( * ), 00030 * $ WORK( * ) 00031 * .. 00032 * 00033 * 00034 *> \par Purpose: 00035 * ============= 00036 *> 00037 *> \verbatim 00038 *> 00039 *> DLAQP2 computes a QR factorization with column pivoting of 00040 *> the block A(OFFSET+1:M,1:N). 00041 *> The block A(1:OFFSET,1:N) is accordingly pivoted, but not factorized. 00042 *> \endverbatim 00043 * 00044 * Arguments: 00045 * ========== 00046 * 00047 *> \param[in] M 00048 *> \verbatim 00049 *> M is INTEGER 00050 *> The number of rows of the matrix A. M >= 0. 00051 *> \endverbatim 00052 *> 00053 *> \param[in] N 00054 *> \verbatim 00055 *> N is INTEGER 00056 *> The number of columns of the matrix A. N >= 0. 00057 *> \endverbatim 00058 *> 00059 *> \param[in] OFFSET 00060 *> \verbatim 00061 *> OFFSET is INTEGER 00062 *> The number of rows of the matrix A that must be pivoted 00063 *> but no factorized. OFFSET >= 0. 00064 *> \endverbatim 00065 *> 00066 *> \param[in,out] A 00067 *> \verbatim 00068 *> A is DOUBLE PRECISION array, dimension (LDA,N) 00069 *> On entry, the M-by-N matrix A. 00070 *> On exit, the upper triangle of block A(OFFSET+1:M,1:N) is 00071 *> the triangular factor obtained; the elements in block 00072 *> A(OFFSET+1:M,1:N) below the diagonal, together with the 00073 *> array TAU, represent the orthogonal matrix Q as a product of 00074 *> elementary reflectors. Block A(1:OFFSET,1:N) has been 00075 *> accordingly pivoted, but no factorized. 00076 *> \endverbatim 00077 *> 00078 *> \param[in] LDA 00079 *> \verbatim 00080 *> LDA is INTEGER 00081 *> The leading dimension of the array A. LDA >= max(1,M). 00082 *> \endverbatim 00083 *> 00084 *> \param[in,out] JPVT 00085 *> \verbatim 00086 *> JPVT is INTEGER array, dimension (N) 00087 *> On entry, if JPVT(i) .ne. 0, the i-th column of A is permuted 00088 *> to the front of A*P (a leading column); if JPVT(i) = 0, 00089 *> the i-th column of A is a free column. 00090 *> On exit, if JPVT(i) = k, then the i-th column of A*P 00091 *> was the k-th column of A. 00092 *> \endverbatim 00093 *> 00094 *> \param[out] TAU 00095 *> \verbatim 00096 *> TAU is DOUBLE PRECISION array, dimension (min(M,N)) 00097 *> The scalar factors of the elementary reflectors. 00098 *> \endverbatim 00099 *> 00100 *> \param[in,out] VN1 00101 *> \verbatim 00102 *> VN1 is DOUBLE PRECISION array, dimension (N) 00103 *> The vector with the partial column norms. 00104 *> \endverbatim 00105 *> 00106 *> \param[in,out] VN2 00107 *> \verbatim 00108 *> VN2 is DOUBLE PRECISION array, dimension (N) 00109 *> The vector with the exact column norms. 00110 *> \endverbatim 00111 *> 00112 *> \param[out] WORK 00113 *> \verbatim 00114 *> WORK is DOUBLE PRECISION array, dimension (N) 00115 *> \endverbatim 00116 * 00117 * Authors: 00118 * ======== 00119 * 00120 *> \author Univ. of Tennessee 00121 *> \author Univ. of California Berkeley 00122 *> \author Univ. of Colorado Denver 00123 *> \author NAG Ltd. 00124 * 00125 *> \date November 2011 00126 * 00127 *> \ingroup doubleOTHERauxiliary 00128 * 00129 *> \par Contributors: 00130 * ================== 00131 *> 00132 *> G. Quintana-Orti, Depto. de Informatica, Universidad Jaime I, Spain 00133 *> X. Sun, Computer Science Dept., Duke University, USA 00134 *> \n 00135 *> Partial column norm updating strategy modified on April 2011 00136 *> Z. Drmac and Z. Bujanovic, Dept. of Mathematics, 00137 *> University of Zagreb, Croatia. 00138 * 00139 *> \par References: 00140 * ================ 00141 *> 00142 *> LAPACK Working Note 176 00143 * 00144 *> \htmlonly 00145 *> <a href="http://www.netlib.org/lapack/lawnspdf/lawn176.pdf">[PDF]</a> 00146 *> \endhtmlonly 00147 * 00148 * ===================================================================== 00149 SUBROUTINE DLAQP2( M, N, OFFSET, A, LDA, JPVT, TAU, VN1, VN2, 00150 $ WORK ) 00151 * 00152 * -- LAPACK auxiliary routine (version 3.4.0) -- 00153 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00154 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00155 * November 2011 00156 * 00157 * .. Scalar Arguments .. 00158 INTEGER LDA, M, N, OFFSET 00159 * .. 00160 * .. Array Arguments .. 00161 INTEGER JPVT( * ) 00162 DOUBLE PRECISION A( LDA, * ), TAU( * ), VN1( * ), VN2( * ), 00163 $ WORK( * ) 00164 * .. 00165 * 00166 * ===================================================================== 00167 * 00168 * .. Parameters .. 00169 DOUBLE PRECISION ZERO, ONE 00170 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 ) 00171 * .. 00172 * .. Local Scalars .. 00173 INTEGER I, ITEMP, J, MN, OFFPI, PVT 00174 DOUBLE PRECISION AII, TEMP, TEMP2, TOL3Z 00175 * .. 00176 * .. External Subroutines .. 00177 EXTERNAL DLARF, DLARFG, DSWAP 00178 * .. 00179 * .. Intrinsic Functions .. 00180 INTRINSIC ABS, MAX, MIN, SQRT 00181 * .. 00182 * .. External Functions .. 00183 INTEGER IDAMAX 00184 DOUBLE PRECISION DLAMCH, DNRM2 00185 EXTERNAL IDAMAX, DLAMCH, DNRM2 00186 * .. 00187 * .. Executable Statements .. 00188 * 00189 MN = MIN( M-OFFSET, N ) 00190 TOL3Z = SQRT(DLAMCH('Epsilon')) 00191 * 00192 * Compute factorization. 00193 * 00194 DO 20 I = 1, MN 00195 * 00196 OFFPI = OFFSET + I 00197 * 00198 * Determine ith pivot column and swap if necessary. 00199 * 00200 PVT = ( I-1 ) + IDAMAX( N-I+1, VN1( I ), 1 ) 00201 * 00202 IF( PVT.NE.I ) THEN 00203 CALL DSWAP( M, A( 1, PVT ), 1, A( 1, I ), 1 ) 00204 ITEMP = JPVT( PVT ) 00205 JPVT( PVT ) = JPVT( I ) 00206 JPVT( I ) = ITEMP 00207 VN1( PVT ) = VN1( I ) 00208 VN2( PVT ) = VN2( I ) 00209 END IF 00210 * 00211 * Generate elementary reflector H(i). 00212 * 00213 IF( OFFPI.LT.M ) THEN 00214 CALL DLARFG( M-OFFPI+1, A( OFFPI, I ), A( OFFPI+1, I ), 1, 00215 $ TAU( I ) ) 00216 ELSE 00217 CALL DLARFG( 1, A( M, I ), A( M, I ), 1, TAU( I ) ) 00218 END IF 00219 * 00220 IF( I.LE.N ) THEN 00221 * 00222 * Apply H(i)**T to A(offset+i:m,i+1:n) from the left. 00223 * 00224 AII = A( OFFPI, I ) 00225 A( OFFPI, I ) = ONE 00226 CALL DLARF( 'Left', M-OFFPI+1, N-I, A( OFFPI, I ), 1, 00227 $ TAU( I ), A( OFFPI, I+1 ), LDA, WORK( 1 ) ) 00228 A( OFFPI, I ) = AII 00229 END IF 00230 * 00231 * Update partial column norms. 00232 * 00233 DO 10 J = I + 1, N 00234 IF( VN1( J ).NE.ZERO ) THEN 00235 * 00236 * NOTE: The following 4 lines follow from the analysis in 00237 * Lapack Working Note 176. 00238 * 00239 TEMP = ONE - ( ABS( A( OFFPI, J ) ) / VN1( J ) )**2 00240 TEMP = MAX( TEMP, ZERO ) 00241 TEMP2 = TEMP*( VN1( J ) / VN2( J ) )**2 00242 IF( TEMP2 .LE. TOL3Z ) THEN 00243 IF( OFFPI.LT.M ) THEN 00244 VN1( J ) = DNRM2( M-OFFPI, A( OFFPI+1, J ), 1 ) 00245 VN2( J ) = VN1( J ) 00246 ELSE 00247 VN1( J ) = ZERO 00248 VN2( J ) = ZERO 00249 END IF 00250 ELSE 00251 VN1( J ) = VN1( J )*SQRT( TEMP ) 00252 END IF 00253 END IF 00254 10 CONTINUE 00255 * 00256 20 CONTINUE 00257 * 00258 RETURN 00259 * 00260 * End of DLAQP2 00261 * 00262 END