1 SUBROUTINE ZBDT02( M, N, B, LDB, C, LDC, U, LDU, WORK, RWORK,
2 $ RESID )
3 *
4 * -- LAPACK test routine (version 3.1) --
5 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
6 * November 2006
7 *
8 * .. Scalar Arguments ..
9 INTEGER LDB, LDC, LDU, M, N
10 DOUBLE PRECISION RESID
11 * ..
12 * .. Array Arguments ..
13 DOUBLE PRECISION RWORK( * )
14 COMPLEX*16 B( LDB, * ), C( LDC, * ), U( LDU, * ),
15 $ WORK( * )
16 * ..
17 *
18 * Purpose
19 * =======
20 *
21 * ZBDT02 tests the change of basis C = U' * B by computing the residual
22 *
23 * RESID = norm( B - U * C ) / ( max(m,n) * norm(B) * EPS ),
24 *
25 * where B and C are M by N matrices, U is an M by M orthogonal matrix,
26 * and EPS is the machine precision.
27 *
28 * Arguments
29 * =========
30 *
31 * M (input) INTEGER
32 * The number of rows of the matrices B and C and the order of
33 * the matrix Q.
34 *
35 * N (input) INTEGER
36 * The number of columns of the matrices B and C.
37 *
38 * B (input) COMPLEX*16 array, dimension (LDB,N)
39 * The m by n matrix B.
40 *
41 * LDB (input) INTEGER
42 * The leading dimension of the array B. LDB >= max(1,M).
43 *
44 * C (input) COMPLEX*16 array, dimension (LDC,N)
45 * The m by n matrix C, assumed to contain U' * B.
46 *
47 * LDC (input) INTEGER
48 * The leading dimension of the array C. LDC >= max(1,M).
49 *
50 * U (input) COMPLEX*16 array, dimension (LDU,M)
51 * The m by m orthogonal matrix U.
52 *
53 * LDU (input) INTEGER
54 * The leading dimension of the array U. LDU >= max(1,M).
55 *
56 * WORK (workspace) COMPLEX*16 array, dimension (M)
57 *
58 * RWORK (workspace) DOUBLE PRECISION array, dimension (M)
59 *
60 * RESID (output) DOUBLE PRECISION
61 * RESID = norm( B - U * C ) / ( max(m,n) * norm(B) * EPS ),
62 *
63 * ======================================================================
64 *
65 * .. Parameters ..
66 DOUBLE PRECISION ZERO, ONE
67 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 )
68 * ..
69 * .. Local Scalars ..
70 INTEGER J
71 DOUBLE PRECISION BNORM, EPS, REALMN
72 * ..
73 * .. External Functions ..
74 DOUBLE PRECISION DLAMCH, DZASUM, ZLANGE
75 EXTERNAL DLAMCH, DZASUM, ZLANGE
76 * ..
77 * .. External Subroutines ..
78 EXTERNAL ZCOPY, ZGEMV
79 * ..
80 * .. Intrinsic Functions ..
81 INTRINSIC DBLE, DCMPLX, MAX, MIN
82 * ..
83 * .. Executable Statements ..
84 *
85 * Quick return if possible
86 *
87 RESID = ZERO
88 IF( M.LE.0 .OR. N.LE.0 )
89 $ RETURN
90 REALMN = DBLE( MAX( M, N ) )
91 EPS = DLAMCH( 'Precision' )
92 *
93 * Compute norm( B - U * C )
94 *
95 DO 10 J = 1, N
96 CALL ZCOPY( M, B( 1, J ), 1, WORK, 1 )
97 CALL ZGEMV( 'No transpose', M, M, -DCMPLX( ONE ), U, LDU,
98 $ C( 1, J ), 1, DCMPLX( ONE ), WORK, 1 )
99 RESID = MAX( RESID, DZASUM( M, WORK, 1 ) )
100 10 CONTINUE
101 *
102 * Compute norm of B.
103 *
104 BNORM = ZLANGE( '1', M, N, B, LDB, RWORK )
105 *
106 IF( BNORM.LE.ZERO ) THEN
107 IF( RESID.NE.ZERO )
108 $ RESID = ONE / EPS
109 ELSE
110 IF( BNORM.GE.RESID ) THEN
111 RESID = ( RESID / BNORM ) / ( REALMN*EPS )
112 ELSE
113 IF( BNORM.LT.ONE ) THEN
114 RESID = ( MIN( RESID, REALMN*BNORM ) / BNORM ) /
115 $ ( REALMN*EPS )
116 ELSE
117 RESID = MIN( RESID / BNORM, REALMN ) / ( REALMN*EPS )
118 END IF
119 END IF
120 END IF
121 RETURN
122 *
123 * End of ZBDT02
124 *
125 END
2 $ RESID )
3 *
4 * -- LAPACK test routine (version 3.1) --
5 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
6 * November 2006
7 *
8 * .. Scalar Arguments ..
9 INTEGER LDB, LDC, LDU, M, N
10 DOUBLE PRECISION RESID
11 * ..
12 * .. Array Arguments ..
13 DOUBLE PRECISION RWORK( * )
14 COMPLEX*16 B( LDB, * ), C( LDC, * ), U( LDU, * ),
15 $ WORK( * )
16 * ..
17 *
18 * Purpose
19 * =======
20 *
21 * ZBDT02 tests the change of basis C = U' * B by computing the residual
22 *
23 * RESID = norm( B - U * C ) / ( max(m,n) * norm(B) * EPS ),
24 *
25 * where B and C are M by N matrices, U is an M by M orthogonal matrix,
26 * and EPS is the machine precision.
27 *
28 * Arguments
29 * =========
30 *
31 * M (input) INTEGER
32 * The number of rows of the matrices B and C and the order of
33 * the matrix Q.
34 *
35 * N (input) INTEGER
36 * The number of columns of the matrices B and C.
37 *
38 * B (input) COMPLEX*16 array, dimension (LDB,N)
39 * The m by n matrix B.
40 *
41 * LDB (input) INTEGER
42 * The leading dimension of the array B. LDB >= max(1,M).
43 *
44 * C (input) COMPLEX*16 array, dimension (LDC,N)
45 * The m by n matrix C, assumed to contain U' * B.
46 *
47 * LDC (input) INTEGER
48 * The leading dimension of the array C. LDC >= max(1,M).
49 *
50 * U (input) COMPLEX*16 array, dimension (LDU,M)
51 * The m by m orthogonal matrix U.
52 *
53 * LDU (input) INTEGER
54 * The leading dimension of the array U. LDU >= max(1,M).
55 *
56 * WORK (workspace) COMPLEX*16 array, dimension (M)
57 *
58 * RWORK (workspace) DOUBLE PRECISION array, dimension (M)
59 *
60 * RESID (output) DOUBLE PRECISION
61 * RESID = norm( B - U * C ) / ( max(m,n) * norm(B) * EPS ),
62 *
63 * ======================================================================
64 *
65 * .. Parameters ..
66 DOUBLE PRECISION ZERO, ONE
67 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 )
68 * ..
69 * .. Local Scalars ..
70 INTEGER J
71 DOUBLE PRECISION BNORM, EPS, REALMN
72 * ..
73 * .. External Functions ..
74 DOUBLE PRECISION DLAMCH, DZASUM, ZLANGE
75 EXTERNAL DLAMCH, DZASUM, ZLANGE
76 * ..
77 * .. External Subroutines ..
78 EXTERNAL ZCOPY, ZGEMV
79 * ..
80 * .. Intrinsic Functions ..
81 INTRINSIC DBLE, DCMPLX, MAX, MIN
82 * ..
83 * .. Executable Statements ..
84 *
85 * Quick return if possible
86 *
87 RESID = ZERO
88 IF( M.LE.0 .OR. N.LE.0 )
89 $ RETURN
90 REALMN = DBLE( MAX( M, N ) )
91 EPS = DLAMCH( 'Precision' )
92 *
93 * Compute norm( B - U * C )
94 *
95 DO 10 J = 1, N
96 CALL ZCOPY( M, B( 1, J ), 1, WORK, 1 )
97 CALL ZGEMV( 'No transpose', M, M, -DCMPLX( ONE ), U, LDU,
98 $ C( 1, J ), 1, DCMPLX( ONE ), WORK, 1 )
99 RESID = MAX( RESID, DZASUM( M, WORK, 1 ) )
100 10 CONTINUE
101 *
102 * Compute norm of B.
103 *
104 BNORM = ZLANGE( '1', M, N, B, LDB, RWORK )
105 *
106 IF( BNORM.LE.ZERO ) THEN
107 IF( RESID.NE.ZERO )
108 $ RESID = ONE / EPS
109 ELSE
110 IF( BNORM.GE.RESID ) THEN
111 RESID = ( RESID / BNORM ) / ( REALMN*EPS )
112 ELSE
113 IF( BNORM.LT.ONE ) THEN
114 RESID = ( MIN( RESID, REALMN*BNORM ) / BNORM ) /
115 $ ( REALMN*EPS )
116 ELSE
117 RESID = MIN( RESID / BNORM, REALMN ) / ( REALMN*EPS )
118 END IF
119 END IF
120 END IF
121 RETURN
122 *
123 * End of ZBDT02
124 *
125 END