1 DOUBLE PRECISION FUNCTION ZQRT17( TRANS, IRESID, M, N, NRHS, A,
2 $ LDA, X, LDX, B, LDB, C, WORK, LWORK )
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 CHARACTER TRANS
10 INTEGER IRESID, LDA, LDB, LDX, LWORK, M, N, NRHS
11 * ..
12 * .. Array Arguments ..
13 COMPLEX*16 A( LDA, * ), B( LDB, * ), C( LDB, * ),
14 $ WORK( LWORK ), X( LDX, * )
15 * ..
16 *
17 * Purpose
18 * =======
19 *
20 * ZQRT17 computes the ratio
21 *
22 * || R'*op(A) ||/(||A||*alpha*max(M,N,NRHS)*eps)
23 *
24 * where R = op(A)*X - B, op(A) is A or A', and
25 *
26 * alpha = ||B|| if IRESID = 1 (zero-residual problem)
27 * alpha = ||R|| if IRESID = 2 (otherwise).
28 *
29 * Arguments
30 * =========
31 *
32 * TRANS (input) CHARACTER*1
33 * Specifies whether or not the transpose of A is used.
34 * = 'N': No transpose, op(A) = A.
35 * = 'C': Conjugate transpose, op(A) = A'.
36 *
37 * IRESID (input) INTEGER
38 * IRESID = 1 indicates zero-residual problem.
39 * IRESID = 2 indicates non-zero residual.
40 *
41 * M (input) INTEGER
42 * The number of rows of the matrix A.
43 * If TRANS = 'N', the number of rows of the matrix B.
44 * If TRANS = 'C', the number of rows of the matrix X.
45 *
46 * N (input) INTEGER
47 * The number of columns of the matrix A.
48 * If TRANS = 'N', the number of rows of the matrix X.
49 * If TRANS = 'C', the number of rows of the matrix B.
50 *
51 * NRHS (input) INTEGER
52 * The number of columns of the matrices X and B.
53 *
54 * A (input) COMPLEX*16 array, dimension (LDA,N)
55 * The m-by-n matrix A.
56 *
57 * LDA (input) INTEGER
58 * The leading dimension of the array A. LDA >= M.
59 *
60 * X (input) COMPLEX*16 array, dimension (LDX,NRHS)
61 * If TRANS = 'N', the n-by-nrhs matrix X.
62 * If TRANS = 'C', the m-by-nrhs matrix X.
63 *
64 * LDX (input) INTEGER
65 * The leading dimension of the array X.
66 * If TRANS = 'N', LDX >= N.
67 * If TRANS = 'C', LDX >= M.
68 *
69 * B (input) COMPLEX*16 array, dimension (LDB,NRHS)
70 * If TRANS = 'N', the m-by-nrhs matrix B.
71 * If TRANS = 'C', the n-by-nrhs matrix B.
72 *
73 * LDB (input) INTEGER
74 * The leading dimension of the array B.
75 * If TRANS = 'N', LDB >= M.
76 * If TRANS = 'C', LDB >= N.
77 *
78 * C (workspace) COMPLEX*16 array, dimension (LDB,NRHS)
79 *
80 * WORK (workspace) COMPLEX*16 array, dimension (LWORK)
81 *
82 * LWORK (input) INTEGER
83 * The length of the array WORK. LWORK >= NRHS*(M+N).
84 *
85 * =====================================================================
86 *
87 * .. Parameters ..
88 DOUBLE PRECISION ZERO, ONE
89 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
90 * ..
91 * .. Local Scalars ..
92 INTEGER INFO, ISCL, NCOLS, NROWS
93 DOUBLE PRECISION BIGNUM, ERR, NORMA, NORMB, NORMRS, NORMX,
94 $ SMLNUM
95 * ..
96 * .. Local Arrays ..
97 DOUBLE PRECISION RWORK( 1 )
98 * ..
99 * .. External Functions ..
100 LOGICAL LSAME
101 DOUBLE PRECISION DLAMCH, ZLANGE
102 EXTERNAL LSAME, DLAMCH, ZLANGE
103 * ..
104 * .. External Subroutines ..
105 EXTERNAL XERBLA, ZGEMM, ZLACPY, ZLASCL
106 * ..
107 * .. Intrinsic Functions ..
108 INTRINSIC DBLE, DCMPLX, MAX
109 * ..
110 * .. Executable Statements ..
111 *
112 ZQRT17 = ZERO
113 *
114 IF( LSAME( TRANS, 'N' ) ) THEN
115 NROWS = M
116 NCOLS = N
117 ELSE IF( LSAME( TRANS, 'C' ) ) THEN
118 NROWS = N
119 NCOLS = M
120 ELSE
121 CALL XERBLA( 'ZQRT17', 1 )
122 RETURN
123 END IF
124 *
125 IF( LWORK.LT.NCOLS*NRHS ) THEN
126 CALL XERBLA( 'ZQRT17', 13 )
127 RETURN
128 END IF
129 *
130 IF( M.LE.0 .OR. N.LE.0 .OR. NRHS.LE.0 )
131 $ RETURN
132 *
133 NORMA = ZLANGE( 'One-norm', M, N, A, LDA, RWORK )
134 SMLNUM = DLAMCH( 'Safe minimum' ) / DLAMCH( 'Precision' )
135 BIGNUM = ONE / SMLNUM
136 ISCL = 0
137 *
138 * compute residual and scale it
139 *
140 CALL ZLACPY( 'All', NROWS, NRHS, B, LDB, C, LDB )
141 CALL ZGEMM( TRANS, 'No transpose', NROWS, NRHS, NCOLS,
142 $ DCMPLX( -ONE ), A, LDA, X, LDX, DCMPLX( ONE ), C,
143 $ LDB )
144 NORMRS = ZLANGE( 'Max', NROWS, NRHS, C, LDB, RWORK )
145 IF( NORMRS.GT.SMLNUM ) THEN
146 ISCL = 1
147 CALL ZLASCL( 'General', 0, 0, NORMRS, ONE, NROWS, NRHS, C, LDB,
148 $ INFO )
149 END IF
150 *
151 * compute R'*A
152 *
153 CALL ZGEMM( 'Conjugate transpose', TRANS, NRHS, NCOLS, NROWS,
154 $ DCMPLX( ONE ), C, LDB, A, LDA, DCMPLX( ZERO ), WORK,
155 $ NRHS )
156 *
157 * compute and properly scale error
158 *
159 ERR = ZLANGE( 'One-norm', NRHS, NCOLS, WORK, NRHS, RWORK )
160 IF( NORMA.NE.ZERO )
161 $ ERR = ERR / NORMA
162 *
163 IF( ISCL.EQ.1 )
164 $ ERR = ERR*NORMRS
165 *
166 IF( IRESID.EQ.1 ) THEN
167 NORMB = ZLANGE( 'One-norm', NROWS, NRHS, B, LDB, RWORK )
168 IF( NORMB.NE.ZERO )
169 $ ERR = ERR / NORMB
170 ELSE
171 NORMX = ZLANGE( 'One-norm', NCOLS, NRHS, X, LDX, RWORK )
172 IF( NORMX.NE.ZERO )
173 $ ERR = ERR / NORMX
174 END IF
175 *
176 ZQRT17 = ERR / ( DLAMCH( 'Epsilon' )*DBLE( MAX( M, N, NRHS ) ) )
177 RETURN
178 *
179 * End of ZQRT17
180 *
181 END
2 $ LDA, X, LDX, B, LDB, C, WORK, LWORK )
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 CHARACTER TRANS
10 INTEGER IRESID, LDA, LDB, LDX, LWORK, M, N, NRHS
11 * ..
12 * .. Array Arguments ..
13 COMPLEX*16 A( LDA, * ), B( LDB, * ), C( LDB, * ),
14 $ WORK( LWORK ), X( LDX, * )
15 * ..
16 *
17 * Purpose
18 * =======
19 *
20 * ZQRT17 computes the ratio
21 *
22 * || R'*op(A) ||/(||A||*alpha*max(M,N,NRHS)*eps)
23 *
24 * where R = op(A)*X - B, op(A) is A or A', and
25 *
26 * alpha = ||B|| if IRESID = 1 (zero-residual problem)
27 * alpha = ||R|| if IRESID = 2 (otherwise).
28 *
29 * Arguments
30 * =========
31 *
32 * TRANS (input) CHARACTER*1
33 * Specifies whether or not the transpose of A is used.
34 * = 'N': No transpose, op(A) = A.
35 * = 'C': Conjugate transpose, op(A) = A'.
36 *
37 * IRESID (input) INTEGER
38 * IRESID = 1 indicates zero-residual problem.
39 * IRESID = 2 indicates non-zero residual.
40 *
41 * M (input) INTEGER
42 * The number of rows of the matrix A.
43 * If TRANS = 'N', the number of rows of the matrix B.
44 * If TRANS = 'C', the number of rows of the matrix X.
45 *
46 * N (input) INTEGER
47 * The number of columns of the matrix A.
48 * If TRANS = 'N', the number of rows of the matrix X.
49 * If TRANS = 'C', the number of rows of the matrix B.
50 *
51 * NRHS (input) INTEGER
52 * The number of columns of the matrices X and B.
53 *
54 * A (input) COMPLEX*16 array, dimension (LDA,N)
55 * The m-by-n matrix A.
56 *
57 * LDA (input) INTEGER
58 * The leading dimension of the array A. LDA >= M.
59 *
60 * X (input) COMPLEX*16 array, dimension (LDX,NRHS)
61 * If TRANS = 'N', the n-by-nrhs matrix X.
62 * If TRANS = 'C', the m-by-nrhs matrix X.
63 *
64 * LDX (input) INTEGER
65 * The leading dimension of the array X.
66 * If TRANS = 'N', LDX >= N.
67 * If TRANS = 'C', LDX >= M.
68 *
69 * B (input) COMPLEX*16 array, dimension (LDB,NRHS)
70 * If TRANS = 'N', the m-by-nrhs matrix B.
71 * If TRANS = 'C', the n-by-nrhs matrix B.
72 *
73 * LDB (input) INTEGER
74 * The leading dimension of the array B.
75 * If TRANS = 'N', LDB >= M.
76 * If TRANS = 'C', LDB >= N.
77 *
78 * C (workspace) COMPLEX*16 array, dimension (LDB,NRHS)
79 *
80 * WORK (workspace) COMPLEX*16 array, dimension (LWORK)
81 *
82 * LWORK (input) INTEGER
83 * The length of the array WORK. LWORK >= NRHS*(M+N).
84 *
85 * =====================================================================
86 *
87 * .. Parameters ..
88 DOUBLE PRECISION ZERO, ONE
89 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
90 * ..
91 * .. Local Scalars ..
92 INTEGER INFO, ISCL, NCOLS, NROWS
93 DOUBLE PRECISION BIGNUM, ERR, NORMA, NORMB, NORMRS, NORMX,
94 $ SMLNUM
95 * ..
96 * .. Local Arrays ..
97 DOUBLE PRECISION RWORK( 1 )
98 * ..
99 * .. External Functions ..
100 LOGICAL LSAME
101 DOUBLE PRECISION DLAMCH, ZLANGE
102 EXTERNAL LSAME, DLAMCH, ZLANGE
103 * ..
104 * .. External Subroutines ..
105 EXTERNAL XERBLA, ZGEMM, ZLACPY, ZLASCL
106 * ..
107 * .. Intrinsic Functions ..
108 INTRINSIC DBLE, DCMPLX, MAX
109 * ..
110 * .. Executable Statements ..
111 *
112 ZQRT17 = ZERO
113 *
114 IF( LSAME( TRANS, 'N' ) ) THEN
115 NROWS = M
116 NCOLS = N
117 ELSE IF( LSAME( TRANS, 'C' ) ) THEN
118 NROWS = N
119 NCOLS = M
120 ELSE
121 CALL XERBLA( 'ZQRT17', 1 )
122 RETURN
123 END IF
124 *
125 IF( LWORK.LT.NCOLS*NRHS ) THEN
126 CALL XERBLA( 'ZQRT17', 13 )
127 RETURN
128 END IF
129 *
130 IF( M.LE.0 .OR. N.LE.0 .OR. NRHS.LE.0 )
131 $ RETURN
132 *
133 NORMA = ZLANGE( 'One-norm', M, N, A, LDA, RWORK )
134 SMLNUM = DLAMCH( 'Safe minimum' ) / DLAMCH( 'Precision' )
135 BIGNUM = ONE / SMLNUM
136 ISCL = 0
137 *
138 * compute residual and scale it
139 *
140 CALL ZLACPY( 'All', NROWS, NRHS, B, LDB, C, LDB )
141 CALL ZGEMM( TRANS, 'No transpose', NROWS, NRHS, NCOLS,
142 $ DCMPLX( -ONE ), A, LDA, X, LDX, DCMPLX( ONE ), C,
143 $ LDB )
144 NORMRS = ZLANGE( 'Max', NROWS, NRHS, C, LDB, RWORK )
145 IF( NORMRS.GT.SMLNUM ) THEN
146 ISCL = 1
147 CALL ZLASCL( 'General', 0, 0, NORMRS, ONE, NROWS, NRHS, C, LDB,
148 $ INFO )
149 END IF
150 *
151 * compute R'*A
152 *
153 CALL ZGEMM( 'Conjugate transpose', TRANS, NRHS, NCOLS, NROWS,
154 $ DCMPLX( ONE ), C, LDB, A, LDA, DCMPLX( ZERO ), WORK,
155 $ NRHS )
156 *
157 * compute and properly scale error
158 *
159 ERR = ZLANGE( 'One-norm', NRHS, NCOLS, WORK, NRHS, RWORK )
160 IF( NORMA.NE.ZERO )
161 $ ERR = ERR / NORMA
162 *
163 IF( ISCL.EQ.1 )
164 $ ERR = ERR*NORMRS
165 *
166 IF( IRESID.EQ.1 ) THEN
167 NORMB = ZLANGE( 'One-norm', NROWS, NRHS, B, LDB, RWORK )
168 IF( NORMB.NE.ZERO )
169 $ ERR = ERR / NORMB
170 ELSE
171 NORMX = ZLANGE( 'One-norm', NCOLS, NRHS, X, LDX, RWORK )
172 IF( NORMX.NE.ZERO )
173 $ ERR = ERR / NORMX
174 END IF
175 *
176 ZQRT17 = ERR / ( DLAMCH( 'Epsilon' )*DBLE( MAX( M, N, NRHS ) ) )
177 RETURN
178 *
179 * End of ZQRT17
180 *
181 END