1 SUBROUTINE DDRVRF3( NOUT, NN, NVAL, THRESH, A, LDA, ARF, B1, B2,
2 + D_WORK_DLANGE, D_WORK_DGEQRF, TAU )
3 *
4 * -- LAPACK test routine (version 3.2.0) --
5 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
6 * November 2008
7 *
8 * .. Scalar Arguments ..
9 INTEGER LDA, NN, NOUT
10 DOUBLE PRECISION THRESH
11 * ..
12 * .. Array Arguments ..
13 INTEGER NVAL( NN )
14 DOUBLE PRECISION A( LDA, * ), ARF( * ), B1( LDA, * ),
15 + B2( LDA, * ), D_WORK_DGEQRF( * ),
16 + D_WORK_DLANGE( * ), TAU( * )
17 * ..
18 *
19 * Purpose
20 * =======
21 *
22 * DDRVRF3 tests the LAPACK RFP routines:
23 * DTFSM
24 *
25 * Arguments
26 * =========
27 *
28 * NOUT (input) INTEGER
29 * The unit number for output.
30 *
31 * NN (input) INTEGER
32 * The number of values of N contained in the vector NVAL.
33 *
34 * NVAL (input) INTEGER array, dimension (NN)
35 * The values of the matrix dimension N.
36 *
37 * THRESH (input) DOUBLE PRECISION
38 * The threshold value for the test ratios. A result is
39 * included in the output file if RESULT >= THRESH. To have
40 * every test ratio printed, use THRESH = 0.
41 *
42 * A (workspace) DOUBLE PRECISION array, dimension (LDA,NMAX)
43 *
44 * LDA (input) INTEGER
45 * The leading dimension of the array A. LDA >= max(1,NMAX).
46 *
47 * ARF (workspace) DOUBLE PRECISION array, dimension ((NMAX*(NMAX+1))/2).
48 *
49 * B1 (workspace) DOUBLE PRECISION array, dimension (LDA,NMAX)
50 *
51 * B2 (workspace) DOUBLE PRECISION array, dimension (LDA,NMAX)
52 *
53 * D_WORK_DLANGE (workspace) DOUBLE PRECISION array, dimension (NMAX)
54 *
55 * D_WORK_DGEQRF (workspace) DOUBLE PRECISION array, dimension (NMAX)
56 *
57 * TAU (workspace) DOUBLE PRECISION array, dimension (NMAX)
58 *
59 * =====================================================================
60 * ..
61 * .. Parameters ..
62 DOUBLE PRECISION ZERO, ONE
63 PARAMETER ( ZERO = ( 0.0D+0, 0.0D+0 ) ,
64 + ONE = ( 1.0D+0, 0.0D+0 ) )
65 INTEGER NTESTS
66 PARAMETER ( NTESTS = 1 )
67 * ..
68 * .. Local Scalars ..
69 CHARACTER UPLO, CFORM, DIAG, TRANS, SIDE
70 INTEGER I, IFORM, IIM, IIN, INFO, IUPLO, J, M, N, NA,
71 + NFAIL, NRUN, ISIDE, IDIAG, IALPHA, ITRANS
72 DOUBLE PRECISION EPS, ALPHA
73 * ..
74 * .. Local Arrays ..
75 CHARACTER UPLOS( 2 ), FORMS( 2 ), TRANSS( 2 ),
76 + DIAGS( 2 ), SIDES( 2 )
77 INTEGER ISEED( 4 ), ISEEDY( 4 )
78 DOUBLE PRECISION RESULT( NTESTS )
79 * ..
80 * .. External Functions ..
81 DOUBLE PRECISION DLAMCH, DLANGE, DLARND
82 EXTERNAL DLAMCH, DLANGE, DLARND
83 * ..
84 * .. External Subroutines ..
85 EXTERNAL DTRTTF, DGEQRF, DGEQLF, DTFSM, DTRSM
86 * ..
87 * .. Intrinsic Functions ..
88 INTRINSIC MAX, SQRT
89 * ..
90 * .. Scalars in Common ..
91 CHARACTER*32 SRNAMT
92 * ..
93 * .. Common blocks ..
94 COMMON / SRNAMC / SRNAMT
95 * ..
96 * .. Data statements ..
97 DATA ISEEDY / 1988, 1989, 1990, 1991 /
98 DATA UPLOS / 'U', 'L' /
99 DATA FORMS / 'N', 'T' /
100 DATA SIDES / 'L', 'R' /
101 DATA TRANSS / 'N', 'T' /
102 DATA DIAGS / 'N', 'U' /
103 * ..
104 * .. Executable Statements ..
105 *
106 * Initialize constants and the random number seed.
107 *
108 NRUN = 0
109 NFAIL = 0
110 INFO = 0
111 DO 10 I = 1, 4
112 ISEED( I ) = ISEEDY( I )
113 10 CONTINUE
114 EPS = DLAMCH( 'Precision' )
115 *
116 DO 170 IIM = 1, NN
117 *
118 M = NVAL( IIM )
119 *
120 DO 160 IIN = 1, NN
121 *
122 N = NVAL( IIN )
123 *
124 DO 150 IFORM = 1, 2
125 *
126 CFORM = FORMS( IFORM )
127 *
128 DO 140 IUPLO = 1, 2
129 *
130 UPLO = UPLOS( IUPLO )
131 *
132 DO 130 ISIDE = 1, 2
133 *
134 SIDE = SIDES( ISIDE )
135 *
136 DO 120 ITRANS = 1, 2
137 *
138 TRANS = TRANSS( ITRANS )
139 *
140 DO 110 IDIAG = 1, 2
141 *
142 DIAG = DIAGS( IDIAG )
143 *
144 DO 100 IALPHA = 1, 3
145 *
146 IF ( IALPHA.EQ. 1) THEN
147 ALPHA = ZERO
148 ELSE IF ( IALPHA.EQ. 1) THEN
149 ALPHA = ONE
150 ELSE
151 ALPHA = DLARND( 2, ISEED )
152 END IF
153 *
154 * All the parameters are set:
155 * CFORM, SIDE, UPLO, TRANS, DIAG, M, N,
156 * and ALPHA
157 * READY TO TEST!
158 *
159 NRUN = NRUN + 1
160 *
161 IF ( ISIDE.EQ.1 ) THEN
162 *
163 * The case ISIDE.EQ.1 is when SIDE.EQ.'L'
164 * -> A is M-by-M ( B is M-by-N )
165 *
166 NA = M
167 *
168 ELSE
169 *
170 * The case ISIDE.EQ.2 is when SIDE.EQ.'R'
171 * -> A is N-by-N ( B is M-by-N )
172 *
173 NA = N
174 *
175 END IF
176 *
177 * Generate A our NA--by--NA triangular
178 * matrix.
179 * Our test is based on forward error so we
180 * do want A to be well conditionned! To get
181 * a well-conditionned triangular matrix, we
182 * take the R factor of the QR/LQ factorization
183 * of a random matrix.
184 *
185 DO J = 1, NA
186 DO I = 1, NA
187 A( I, J) = DLARND( 2, ISEED )
188 END DO
189 END DO
190 *
191 IF ( IUPLO.EQ.1 ) THEN
192 *
193 * The case IUPLO.EQ.1 is when SIDE.EQ.'U'
194 * -> QR factorization.
195 *
196 SRNAMT = 'DGEQRF'
197 CALL DGEQRF( NA, NA, A, LDA, TAU,
198 + D_WORK_DGEQRF, LDA,
199 + INFO )
200 ELSE
201 *
202 * The case IUPLO.EQ.2 is when SIDE.EQ.'L'
203 * -> QL factorization.
204 *
205 SRNAMT = 'DGELQF'
206 CALL DGELQF( NA, NA, A, LDA, TAU,
207 + D_WORK_DGEQRF, LDA,
208 + INFO )
209 END IF
210 *
211 * Store a copy of A in RFP format (in ARF).
212 *
213 SRNAMT = 'DTRTTF'
214 CALL DTRTTF( CFORM, UPLO, NA, A, LDA, ARF,
215 + INFO )
216 *
217 * Generate B1 our M--by--N right-hand side
218 * and store a copy in B2.
219 *
220 DO J = 1, N
221 DO I = 1, M
222 B1( I, J) = DLARND( 2, ISEED )
223 B2( I, J) = B1( I, J)
224 END DO
225 END DO
226 *
227 * Solve op( A ) X = B or X op( A ) = B
228 * with DTRSM
229 *
230 SRNAMT = 'DTRSM'
231 CALL DTRSM( SIDE, UPLO, TRANS, DIAG, M, N,
232 + ALPHA, A, LDA, B1, LDA )
233 *
234 * Solve op( A ) X = B or X op( A ) = B
235 * with DTFSM
236 *
237 SRNAMT = 'DTFSM'
238 CALL DTFSM( CFORM, SIDE, UPLO, TRANS,
239 + DIAG, M, N, ALPHA, ARF, B2,
240 + LDA )
241 *
242 * Check that the result agrees.
243 *
244 DO J = 1, N
245 DO I = 1, M
246 B1( I, J) = B2( I, J ) - B1( I, J )
247 END DO
248 END DO
249 *
250 RESULT(1) = DLANGE( 'I', M, N, B1, LDA,
251 + D_WORK_DLANGE )
252 *
253 RESULT(1) = RESULT(1) / SQRT( EPS )
254 + / MAX ( MAX( M, N), 1 )
255 *
256 IF( RESULT(1).GE.THRESH ) THEN
257 IF( NFAIL.EQ.0 ) THEN
258 WRITE( NOUT, * )
259 WRITE( NOUT, FMT = 9999 )
260 END IF
261 WRITE( NOUT, FMT = 9997 ) 'DTFSM',
262 + CFORM, SIDE, UPLO, TRANS, DIAG, M,
263 + N, RESULT(1)
264 NFAIL = NFAIL + 1
265 END IF
266 *
267 100 CONTINUE
268 110 CONTINUE
269 120 CONTINUE
270 130 CONTINUE
271 140 CONTINUE
272 150 CONTINUE
273 160 CONTINUE
274 170 CONTINUE
275 *
276 * Print a summary of the results.
277 *
278 IF ( NFAIL.EQ.0 ) THEN
279 WRITE( NOUT, FMT = 9996 ) 'DTFSM', NRUN
280 ELSE
281 WRITE( NOUT, FMT = 9995 ) 'DTFSM', NFAIL, NRUN
282 END IF
283 *
284 9999 FORMAT( 1X, ' *** Error(s) or Failure(s) while testing DTFSM
285 + ***')
286 9997 FORMAT( 1X, ' Failure in ',A5,', CFORM=''',A1,''',',
287 + ' SIDE=''',A1,''',',' UPLO=''',A1,''',',' TRANS=''',A1,''',',
288 + ' DIAG=''',A1,''',',' M=',I3,', N =', I3,', test=',G12.5)
289 9996 FORMAT( 1X, 'All tests for ',A5,' auxiliary routine passed the ',
290 + 'threshold (',I5,' tests run)')
291 9995 FORMAT( 1X, A6, ' auxiliary routine:',I5,' out of ',I5,
292 + ' tests failed to pass the threshold')
293 *
294 RETURN
295 *
296 * End of DDRVRF3
297 *
298 END
2 + D_WORK_DLANGE, D_WORK_DGEQRF, TAU )
3 *
4 * -- LAPACK test routine (version 3.2.0) --
5 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
6 * November 2008
7 *
8 * .. Scalar Arguments ..
9 INTEGER LDA, NN, NOUT
10 DOUBLE PRECISION THRESH
11 * ..
12 * .. Array Arguments ..
13 INTEGER NVAL( NN )
14 DOUBLE PRECISION A( LDA, * ), ARF( * ), B1( LDA, * ),
15 + B2( LDA, * ), D_WORK_DGEQRF( * ),
16 + D_WORK_DLANGE( * ), TAU( * )
17 * ..
18 *
19 * Purpose
20 * =======
21 *
22 * DDRVRF3 tests the LAPACK RFP routines:
23 * DTFSM
24 *
25 * Arguments
26 * =========
27 *
28 * NOUT (input) INTEGER
29 * The unit number for output.
30 *
31 * NN (input) INTEGER
32 * The number of values of N contained in the vector NVAL.
33 *
34 * NVAL (input) INTEGER array, dimension (NN)
35 * The values of the matrix dimension N.
36 *
37 * THRESH (input) DOUBLE PRECISION
38 * The threshold value for the test ratios. A result is
39 * included in the output file if RESULT >= THRESH. To have
40 * every test ratio printed, use THRESH = 0.
41 *
42 * A (workspace) DOUBLE PRECISION array, dimension (LDA,NMAX)
43 *
44 * LDA (input) INTEGER
45 * The leading dimension of the array A. LDA >= max(1,NMAX).
46 *
47 * ARF (workspace) DOUBLE PRECISION array, dimension ((NMAX*(NMAX+1))/2).
48 *
49 * B1 (workspace) DOUBLE PRECISION array, dimension (LDA,NMAX)
50 *
51 * B2 (workspace) DOUBLE PRECISION array, dimension (LDA,NMAX)
52 *
53 * D_WORK_DLANGE (workspace) DOUBLE PRECISION array, dimension (NMAX)
54 *
55 * D_WORK_DGEQRF (workspace) DOUBLE PRECISION array, dimension (NMAX)
56 *
57 * TAU (workspace) DOUBLE PRECISION array, dimension (NMAX)
58 *
59 * =====================================================================
60 * ..
61 * .. Parameters ..
62 DOUBLE PRECISION ZERO, ONE
63 PARAMETER ( ZERO = ( 0.0D+0, 0.0D+0 ) ,
64 + ONE = ( 1.0D+0, 0.0D+0 ) )
65 INTEGER NTESTS
66 PARAMETER ( NTESTS = 1 )
67 * ..
68 * .. Local Scalars ..
69 CHARACTER UPLO, CFORM, DIAG, TRANS, SIDE
70 INTEGER I, IFORM, IIM, IIN, INFO, IUPLO, J, M, N, NA,
71 + NFAIL, NRUN, ISIDE, IDIAG, IALPHA, ITRANS
72 DOUBLE PRECISION EPS, ALPHA
73 * ..
74 * .. Local Arrays ..
75 CHARACTER UPLOS( 2 ), FORMS( 2 ), TRANSS( 2 ),
76 + DIAGS( 2 ), SIDES( 2 )
77 INTEGER ISEED( 4 ), ISEEDY( 4 )
78 DOUBLE PRECISION RESULT( NTESTS )
79 * ..
80 * .. External Functions ..
81 DOUBLE PRECISION DLAMCH, DLANGE, DLARND
82 EXTERNAL DLAMCH, DLANGE, DLARND
83 * ..
84 * .. External Subroutines ..
85 EXTERNAL DTRTTF, DGEQRF, DGEQLF, DTFSM, DTRSM
86 * ..
87 * .. Intrinsic Functions ..
88 INTRINSIC MAX, SQRT
89 * ..
90 * .. Scalars in Common ..
91 CHARACTER*32 SRNAMT
92 * ..
93 * .. Common blocks ..
94 COMMON / SRNAMC / SRNAMT
95 * ..
96 * .. Data statements ..
97 DATA ISEEDY / 1988, 1989, 1990, 1991 /
98 DATA UPLOS / 'U', 'L' /
99 DATA FORMS / 'N', 'T' /
100 DATA SIDES / 'L', 'R' /
101 DATA TRANSS / 'N', 'T' /
102 DATA DIAGS / 'N', 'U' /
103 * ..
104 * .. Executable Statements ..
105 *
106 * Initialize constants and the random number seed.
107 *
108 NRUN = 0
109 NFAIL = 0
110 INFO = 0
111 DO 10 I = 1, 4
112 ISEED( I ) = ISEEDY( I )
113 10 CONTINUE
114 EPS = DLAMCH( 'Precision' )
115 *
116 DO 170 IIM = 1, NN
117 *
118 M = NVAL( IIM )
119 *
120 DO 160 IIN = 1, NN
121 *
122 N = NVAL( IIN )
123 *
124 DO 150 IFORM = 1, 2
125 *
126 CFORM = FORMS( IFORM )
127 *
128 DO 140 IUPLO = 1, 2
129 *
130 UPLO = UPLOS( IUPLO )
131 *
132 DO 130 ISIDE = 1, 2
133 *
134 SIDE = SIDES( ISIDE )
135 *
136 DO 120 ITRANS = 1, 2
137 *
138 TRANS = TRANSS( ITRANS )
139 *
140 DO 110 IDIAG = 1, 2
141 *
142 DIAG = DIAGS( IDIAG )
143 *
144 DO 100 IALPHA = 1, 3
145 *
146 IF ( IALPHA.EQ. 1) THEN
147 ALPHA = ZERO
148 ELSE IF ( IALPHA.EQ. 1) THEN
149 ALPHA = ONE
150 ELSE
151 ALPHA = DLARND( 2, ISEED )
152 END IF
153 *
154 * All the parameters are set:
155 * CFORM, SIDE, UPLO, TRANS, DIAG, M, N,
156 * and ALPHA
157 * READY TO TEST!
158 *
159 NRUN = NRUN + 1
160 *
161 IF ( ISIDE.EQ.1 ) THEN
162 *
163 * The case ISIDE.EQ.1 is when SIDE.EQ.'L'
164 * -> A is M-by-M ( B is M-by-N )
165 *
166 NA = M
167 *
168 ELSE
169 *
170 * The case ISIDE.EQ.2 is when SIDE.EQ.'R'
171 * -> A is N-by-N ( B is M-by-N )
172 *
173 NA = N
174 *
175 END IF
176 *
177 * Generate A our NA--by--NA triangular
178 * matrix.
179 * Our test is based on forward error so we
180 * do want A to be well conditionned! To get
181 * a well-conditionned triangular matrix, we
182 * take the R factor of the QR/LQ factorization
183 * of a random matrix.
184 *
185 DO J = 1, NA
186 DO I = 1, NA
187 A( I, J) = DLARND( 2, ISEED )
188 END DO
189 END DO
190 *
191 IF ( IUPLO.EQ.1 ) THEN
192 *
193 * The case IUPLO.EQ.1 is when SIDE.EQ.'U'
194 * -> QR factorization.
195 *
196 SRNAMT = 'DGEQRF'
197 CALL DGEQRF( NA, NA, A, LDA, TAU,
198 + D_WORK_DGEQRF, LDA,
199 + INFO )
200 ELSE
201 *
202 * The case IUPLO.EQ.2 is when SIDE.EQ.'L'
203 * -> QL factorization.
204 *
205 SRNAMT = 'DGELQF'
206 CALL DGELQF( NA, NA, A, LDA, TAU,
207 + D_WORK_DGEQRF, LDA,
208 + INFO )
209 END IF
210 *
211 * Store a copy of A in RFP format (in ARF).
212 *
213 SRNAMT = 'DTRTTF'
214 CALL DTRTTF( CFORM, UPLO, NA, A, LDA, ARF,
215 + INFO )
216 *
217 * Generate B1 our M--by--N right-hand side
218 * and store a copy in B2.
219 *
220 DO J = 1, N
221 DO I = 1, M
222 B1( I, J) = DLARND( 2, ISEED )
223 B2( I, J) = B1( I, J)
224 END DO
225 END DO
226 *
227 * Solve op( A ) X = B or X op( A ) = B
228 * with DTRSM
229 *
230 SRNAMT = 'DTRSM'
231 CALL DTRSM( SIDE, UPLO, TRANS, DIAG, M, N,
232 + ALPHA, A, LDA, B1, LDA )
233 *
234 * Solve op( A ) X = B or X op( A ) = B
235 * with DTFSM
236 *
237 SRNAMT = 'DTFSM'
238 CALL DTFSM( CFORM, SIDE, UPLO, TRANS,
239 + DIAG, M, N, ALPHA, ARF, B2,
240 + LDA )
241 *
242 * Check that the result agrees.
243 *
244 DO J = 1, N
245 DO I = 1, M
246 B1( I, J) = B2( I, J ) - B1( I, J )
247 END DO
248 END DO
249 *
250 RESULT(1) = DLANGE( 'I', M, N, B1, LDA,
251 + D_WORK_DLANGE )
252 *
253 RESULT(1) = RESULT(1) / SQRT( EPS )
254 + / MAX ( MAX( M, N), 1 )
255 *
256 IF( RESULT(1).GE.THRESH ) THEN
257 IF( NFAIL.EQ.0 ) THEN
258 WRITE( NOUT, * )
259 WRITE( NOUT, FMT = 9999 )
260 END IF
261 WRITE( NOUT, FMT = 9997 ) 'DTFSM',
262 + CFORM, SIDE, UPLO, TRANS, DIAG, M,
263 + N, RESULT(1)
264 NFAIL = NFAIL + 1
265 END IF
266 *
267 100 CONTINUE
268 110 CONTINUE
269 120 CONTINUE
270 130 CONTINUE
271 140 CONTINUE
272 150 CONTINUE
273 160 CONTINUE
274 170 CONTINUE
275 *
276 * Print a summary of the results.
277 *
278 IF ( NFAIL.EQ.0 ) THEN
279 WRITE( NOUT, FMT = 9996 ) 'DTFSM', NRUN
280 ELSE
281 WRITE( NOUT, FMT = 9995 ) 'DTFSM', NFAIL, NRUN
282 END IF
283 *
284 9999 FORMAT( 1X, ' *** Error(s) or Failure(s) while testing DTFSM
285 + ***')
286 9997 FORMAT( 1X, ' Failure in ',A5,', CFORM=''',A1,''',',
287 + ' SIDE=''',A1,''',',' UPLO=''',A1,''',',' TRANS=''',A1,''',',
288 + ' DIAG=''',A1,''',',' M=',I3,', N =', I3,', test=',G12.5)
289 9996 FORMAT( 1X, 'All tests for ',A5,' auxiliary routine passed the ',
290 + 'threshold (',I5,' tests run)')
291 9995 FORMAT( 1X, A6, ' auxiliary routine:',I5,' out of ',I5,
292 + ' tests failed to pass the threshold')
293 *
294 RETURN
295 *
296 * End of DDRVRF3
297 *
298 END