![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DQPT01 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 * Definition: 00009 * =========== 00010 * 00011 * DOUBLE PRECISION FUNCTION DQPT01( M, N, K, A, AF, LDA, TAU, JPVT, 00012 * WORK, LWORK ) 00013 * 00014 * .. Scalar Arguments .. 00015 * INTEGER K, LDA, LWORK, M, N 00016 * .. 00017 * .. Array Arguments .. 00018 * INTEGER JPVT( * ) 00019 * DOUBLE PRECISION A( LDA, * ), AF( LDA, * ), TAU( * ), 00020 * $ WORK( LWORK ) 00021 * .. 00022 * 00023 * 00024 *> \par Purpose: 00025 * ============= 00026 *> 00027 *> \verbatim 00028 *> 00029 *> DQPT01 tests the QR-factorization with pivoting of a matrix A. The 00030 *> array AF contains the (possibly partial) QR-factorization of A, where 00031 *> the upper triangle of AF(1:k,1:k) is a partial triangular factor, 00032 *> the entries below the diagonal in the first k columns are the 00033 *> Householder vectors, and the rest of AF contains a partially updated 00034 *> matrix. 00035 *> 00036 *> This function returns ||A*P - Q*R||/(||norm(A)||*eps*M) 00037 *> \endverbatim 00038 * 00039 * Arguments: 00040 * ========== 00041 * 00042 *> \param[in] M 00043 *> \verbatim 00044 *> M is INTEGER 00045 *> The number of rows of the matrices A and AF. 00046 *> \endverbatim 00047 *> 00048 *> \param[in] N 00049 *> \verbatim 00050 *> N is INTEGER 00051 *> The number of columns of the matrices A and AF. 00052 *> \endverbatim 00053 *> 00054 *> \param[in] K 00055 *> \verbatim 00056 *> K is INTEGER 00057 *> The number of columns of AF that have been reduced 00058 *> to upper triangular form. 00059 *> \endverbatim 00060 *> 00061 *> \param[in] A 00062 *> \verbatim 00063 *> A is DOUBLE PRECISION array, dimension (LDA, N) 00064 *> The original matrix A. 00065 *> \endverbatim 00066 *> 00067 *> \param[in] AF 00068 *> \verbatim 00069 *> AF is DOUBLE PRECISION array, dimension (LDA,N) 00070 *> The (possibly partial) output of DGEQPF. The upper triangle 00071 *> of AF(1:k,1:k) is a partial triangular factor, the entries 00072 *> below the diagonal in the first k columns are the Householder 00073 *> vectors, and the rest of AF contains a partially updated 00074 *> matrix. 00075 *> \endverbatim 00076 *> 00077 *> \param[in] LDA 00078 *> \verbatim 00079 *> LDA is INTEGER 00080 *> The leading dimension of the arrays A and AF. 00081 *> \endverbatim 00082 *> 00083 *> \param[in] TAU 00084 *> \verbatim 00085 *> TAU is DOUBLE PRECISION array, dimension (K) 00086 *> Details of the Householder transformations as returned by 00087 *> DGEQPF. 00088 *> \endverbatim 00089 *> 00090 *> \param[in] JPVT 00091 *> \verbatim 00092 *> JPVT is INTEGER array, dimension (N) 00093 *> Pivot information as returned by DGEQPF. 00094 *> \endverbatim 00095 *> 00096 *> \param[out] WORK 00097 *> \verbatim 00098 *> WORK is DOUBLE PRECISION array, dimension (LWORK) 00099 *> \endverbatim 00100 *> 00101 *> \param[in] LWORK 00102 *> \verbatim 00103 *> LWORK is INTEGER 00104 *> The length of the array WORK. LWORK >= M*N+N. 00105 *> \endverbatim 00106 * 00107 * Authors: 00108 * ======== 00109 * 00110 *> \author Univ. of Tennessee 00111 *> \author Univ. of California Berkeley 00112 *> \author Univ. of Colorado Denver 00113 *> \author NAG Ltd. 00114 * 00115 *> \date November 2011 00116 * 00117 *> \ingroup double_lin 00118 * 00119 * ===================================================================== 00120 DOUBLE PRECISION FUNCTION DQPT01( M, N, K, A, AF, LDA, TAU, JPVT, 00121 $ WORK, LWORK ) 00122 * 00123 * -- LAPACK test routine (version 3.4.0) -- 00124 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00125 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00126 * November 2011 00127 * 00128 * .. Scalar Arguments .. 00129 INTEGER K, LDA, LWORK, M, N 00130 * .. 00131 * .. Array Arguments .. 00132 INTEGER JPVT( * ) 00133 DOUBLE PRECISION A( LDA, * ), AF( LDA, * ), TAU( * ), 00134 $ WORK( LWORK ) 00135 * .. 00136 * 00137 * ===================================================================== 00138 * 00139 * .. Parameters .. 00140 DOUBLE PRECISION ZERO, ONE 00141 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 ) 00142 * .. 00143 * .. Local Scalars .. 00144 INTEGER I, INFO, J 00145 DOUBLE PRECISION NORMA 00146 * .. 00147 * .. Local Arrays .. 00148 DOUBLE PRECISION RWORK( 1 ) 00149 * .. 00150 * .. External Functions .. 00151 DOUBLE PRECISION DLAMCH, DLANGE 00152 EXTERNAL DLAMCH, DLANGE 00153 * .. 00154 * .. External Subroutines .. 00155 EXTERNAL DAXPY, DCOPY, DORMQR, XERBLA 00156 * .. 00157 * .. Intrinsic Functions .. 00158 INTRINSIC DBLE, MAX, MIN 00159 * .. 00160 * .. Executable Statements .. 00161 * 00162 DQPT01 = ZERO 00163 * 00164 * Test if there is enough workspace 00165 * 00166 IF( LWORK.LT.M*N+N ) THEN 00167 CALL XERBLA( 'DQPT01', 10 ) 00168 RETURN 00169 END IF 00170 * 00171 * Quick return if possible 00172 * 00173 IF( M.LE.0 .OR. N.LE.0 ) 00174 $ RETURN 00175 * 00176 NORMA = DLANGE( 'One-norm', M, N, A, LDA, RWORK ) 00177 * 00178 DO 30 J = 1, K 00179 DO 10 I = 1, MIN( J, M ) 00180 WORK( ( J-1 )*M+I ) = AF( I, J ) 00181 10 CONTINUE 00182 DO 20 I = J + 1, M 00183 WORK( ( J-1 )*M+I ) = ZERO 00184 20 CONTINUE 00185 30 CONTINUE 00186 DO 40 J = K + 1, N 00187 CALL DCOPY( M, AF( 1, J ), 1, WORK( ( J-1 )*M+1 ), 1 ) 00188 40 CONTINUE 00189 * 00190 CALL DORMQR( 'Left', 'No transpose', M, N, K, AF, LDA, TAU, WORK, 00191 $ M, WORK( M*N+1 ), LWORK-M*N, INFO ) 00192 * 00193 DO 50 J = 1, N 00194 * 00195 * Compare i-th column of QR and jpvt(i)-th column of A 00196 * 00197 CALL DAXPY( M, -ONE, A( 1, JPVT( J ) ), 1, WORK( ( J-1 )*M+1 ), 00198 $ 1 ) 00199 50 CONTINUE 00200 * 00201 DQPT01 = DLANGE( 'One-norm', M, N, WORK, M, RWORK ) / 00202 $ ( DBLE( MAX( M, N ) )*DLAMCH( 'Epsilon' ) ) 00203 IF( NORMA.NE.ZERO ) 00204 $ DQPT01 = DQPT01 / NORMA 00205 * 00206 RETURN 00207 * 00208 * End of DQPT01 00209 * 00210 END