![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b ZLQT03 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 * Definition: 00009 * =========== 00010 * 00011 * SUBROUTINE ZLQT03( M, N, K, AF, C, CC, Q, LDA, TAU, WORK, LWORK, 00012 * RWORK, RESULT ) 00013 * 00014 * .. Scalar Arguments .. 00015 * INTEGER K, LDA, LWORK, M, N 00016 * .. 00017 * .. Array Arguments .. 00018 * DOUBLE PRECISION RESULT( * ), RWORK( * ) 00019 * COMPLEX*16 AF( LDA, * ), C( LDA, * ), CC( LDA, * ), 00020 * $ Q( LDA, * ), TAU( * ), WORK( LWORK ) 00021 * .. 00022 * 00023 * 00024 *> \par Purpose: 00025 * ============= 00026 *> 00027 *> \verbatim 00028 *> 00029 *> ZLQT03 tests ZUNMLQ, which computes Q*C, Q'*C, C*Q or C*Q'. 00030 *> 00031 *> ZLQT03 compares the results of a call to ZUNMLQ with the results of 00032 *> forming Q explicitly by a call to ZUNGLQ and then performing matrix 00033 *> multiplication by a call to ZGEMM. 00034 *> \endverbatim 00035 * 00036 * Arguments: 00037 * ========== 00038 * 00039 *> \param[in] M 00040 *> \verbatim 00041 *> M is INTEGER 00042 *> The number of rows or columns of the matrix C; C is n-by-m if 00043 *> Q is applied from the left, or m-by-n if Q is applied from 00044 *> the right. M >= 0. 00045 *> \endverbatim 00046 *> 00047 *> \param[in] N 00048 *> \verbatim 00049 *> N is INTEGER 00050 *> The order of the orthogonal matrix Q. N >= 0. 00051 *> \endverbatim 00052 *> 00053 *> \param[in] K 00054 *> \verbatim 00055 *> K is INTEGER 00056 *> The number of elementary reflectors whose product defines the 00057 *> orthogonal matrix Q. N >= K >= 0. 00058 *> \endverbatim 00059 *> 00060 *> \param[in] AF 00061 *> \verbatim 00062 *> AF is COMPLEX*16 array, dimension (LDA,N) 00063 *> Details of the LQ factorization of an m-by-n matrix, as 00064 *> returned by ZGELQF. See CGELQF for further details. 00065 *> \endverbatim 00066 *> 00067 *> \param[out] C 00068 *> \verbatim 00069 *> C is COMPLEX*16 array, dimension (LDA,N) 00070 *> \endverbatim 00071 *> 00072 *> \param[out] CC 00073 *> \verbatim 00074 *> CC is COMPLEX*16 array, dimension (LDA,N) 00075 *> \endverbatim 00076 *> 00077 *> \param[out] Q 00078 *> \verbatim 00079 *> Q is COMPLEX*16 array, dimension (LDA,N) 00080 *> \endverbatim 00081 *> 00082 *> \param[in] LDA 00083 *> \verbatim 00084 *> LDA is INTEGER 00085 *> The leading dimension of the arrays AF, C, CC, and Q. 00086 *> \endverbatim 00087 *> 00088 *> \param[in] TAU 00089 *> \verbatim 00090 *> TAU is COMPLEX*16 array, dimension (min(M,N)) 00091 *> The scalar factors of the elementary reflectors corresponding 00092 *> to the LQ factorization in AF. 00093 *> \endverbatim 00094 *> 00095 *> \param[out] WORK 00096 *> \verbatim 00097 *> WORK is COMPLEX*16 array, dimension (LWORK) 00098 *> \endverbatim 00099 *> 00100 *> \param[in] LWORK 00101 *> \verbatim 00102 *> LWORK is INTEGER 00103 *> The length of WORK. LWORK must be at least M, and should be 00104 *> M*NB, where NB is the blocksize for this environment. 00105 *> \endverbatim 00106 *> 00107 *> \param[out] RWORK 00108 *> \verbatim 00109 *> RWORK is DOUBLE PRECISION array, dimension (M) 00110 *> \endverbatim 00111 *> 00112 *> \param[out] RESULT 00113 *> \verbatim 00114 *> RESULT is DOUBLE PRECISION array, dimension (4) 00115 *> The test ratios compare two techniques for multiplying a 00116 *> random matrix C by an n-by-n orthogonal matrix Q. 00117 *> RESULT(1) = norm( Q*C - Q*C ) / ( N * norm(C) * EPS ) 00118 *> RESULT(2) = norm( C*Q - C*Q ) / ( N * norm(C) * EPS ) 00119 *> RESULT(3) = norm( Q'*C - Q'*C )/ ( N * norm(C) * EPS ) 00120 *> RESULT(4) = norm( C*Q' - C*Q' )/ ( N * norm(C) * EPS ) 00121 *> \endverbatim 00122 * 00123 * Authors: 00124 * ======== 00125 * 00126 *> \author Univ. of Tennessee 00127 *> \author Univ. of California Berkeley 00128 *> \author Univ. of Colorado Denver 00129 *> \author NAG Ltd. 00130 * 00131 *> \date November 2011 00132 * 00133 *> \ingroup complex16_lin 00134 * 00135 * ===================================================================== 00136 SUBROUTINE ZLQT03( M, N, K, AF, C, CC, Q, LDA, TAU, WORK, LWORK, 00137 $ RWORK, RESULT ) 00138 * 00139 * -- LAPACK test routine (version 3.4.0) -- 00140 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00141 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00142 * November 2011 00143 * 00144 * .. Scalar Arguments .. 00145 INTEGER K, LDA, LWORK, M, N 00146 * .. 00147 * .. Array Arguments .. 00148 DOUBLE PRECISION RESULT( * ), RWORK( * ) 00149 COMPLEX*16 AF( LDA, * ), C( LDA, * ), CC( LDA, * ), 00150 $ Q( LDA, * ), TAU( * ), WORK( LWORK ) 00151 * .. 00152 * 00153 * ===================================================================== 00154 * 00155 * .. Parameters .. 00156 DOUBLE PRECISION ZERO, ONE 00157 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 ) 00158 COMPLEX*16 ROGUE 00159 PARAMETER ( ROGUE = ( -1.0D+10, -1.0D+10 ) ) 00160 * .. 00161 * .. Local Scalars .. 00162 CHARACTER SIDE, TRANS 00163 INTEGER INFO, ISIDE, ITRANS, J, MC, NC 00164 DOUBLE PRECISION CNORM, EPS, RESID 00165 * .. 00166 * .. External Functions .. 00167 LOGICAL LSAME 00168 DOUBLE PRECISION DLAMCH, ZLANGE 00169 EXTERNAL LSAME, DLAMCH, ZLANGE 00170 * .. 00171 * .. External Subroutines .. 00172 EXTERNAL ZGEMM, ZLACPY, ZLARNV, ZLASET, ZUNGLQ, ZUNMLQ 00173 * .. 00174 * .. Local Arrays .. 00175 INTEGER ISEED( 4 ) 00176 * .. 00177 * .. Intrinsic Functions .. 00178 INTRINSIC DBLE, DCMPLX, MAX 00179 * .. 00180 * .. Scalars in Common .. 00181 CHARACTER*32 SRNAMT 00182 * .. 00183 * .. Common blocks .. 00184 COMMON / SRNAMC / SRNAMT 00185 * .. 00186 * .. Data statements .. 00187 DATA ISEED / 1988, 1989, 1990, 1991 / 00188 * .. 00189 * .. Executable Statements .. 00190 * 00191 EPS = DLAMCH( 'Epsilon' ) 00192 * 00193 * Copy the first k rows of the factorization to the array Q 00194 * 00195 CALL ZLASET( 'Full', N, N, ROGUE, ROGUE, Q, LDA ) 00196 CALL ZLACPY( 'Upper', K, N-1, AF( 1, 2 ), LDA, Q( 1, 2 ), LDA ) 00197 * 00198 * Generate the n-by-n matrix Q 00199 * 00200 SRNAMT = 'ZUNGLQ' 00201 CALL ZUNGLQ( N, N, K, Q, LDA, TAU, WORK, LWORK, INFO ) 00202 * 00203 DO 30 ISIDE = 1, 2 00204 IF( ISIDE.EQ.1 ) THEN 00205 SIDE = 'L' 00206 MC = N 00207 NC = M 00208 ELSE 00209 SIDE = 'R' 00210 MC = M 00211 NC = N 00212 END IF 00213 * 00214 * Generate MC by NC matrix C 00215 * 00216 DO 10 J = 1, NC 00217 CALL ZLARNV( 2, ISEED, MC, C( 1, J ) ) 00218 10 CONTINUE 00219 CNORM = ZLANGE( '1', MC, NC, C, LDA, RWORK ) 00220 IF( CNORM.EQ.ZERO ) 00221 $ CNORM = ONE 00222 * 00223 DO 20 ITRANS = 1, 2 00224 IF( ITRANS.EQ.1 ) THEN 00225 TRANS = 'N' 00226 ELSE 00227 TRANS = 'C' 00228 END IF 00229 * 00230 * Copy C 00231 * 00232 CALL ZLACPY( 'Full', MC, NC, C, LDA, CC, LDA ) 00233 * 00234 * Apply Q or Q' to C 00235 * 00236 SRNAMT = 'ZUNMLQ' 00237 CALL ZUNMLQ( SIDE, TRANS, MC, NC, K, AF, LDA, TAU, CC, LDA, 00238 $ WORK, LWORK, INFO ) 00239 * 00240 * Form explicit product and subtract 00241 * 00242 IF( LSAME( SIDE, 'L' ) ) THEN 00243 CALL ZGEMM( TRANS, 'No transpose', MC, NC, MC, 00244 $ DCMPLX( -ONE ), Q, LDA, C, LDA, 00245 $ DCMPLX( ONE ), CC, LDA ) 00246 ELSE 00247 CALL ZGEMM( 'No transpose', TRANS, MC, NC, NC, 00248 $ DCMPLX( -ONE ), C, LDA, Q, LDA, 00249 $ DCMPLX( ONE ), CC, LDA ) 00250 END IF 00251 * 00252 * Compute error in the difference 00253 * 00254 RESID = ZLANGE( '1', MC, NC, CC, LDA, RWORK ) 00255 RESULT( ( ISIDE-1 )*2+ITRANS ) = RESID / 00256 $ ( DBLE( MAX( 1, N ) )*CNORM*EPS ) 00257 * 00258 20 CONTINUE 00259 30 CONTINUE 00260 * 00261 RETURN 00262 * 00263 * End of ZLQT03 00264 * 00265 END