1 SUBROUTINE ZGET35( RMAX, LMAX, NINFO, KNT, NIN )
2 *
3 * -- LAPACK test routine (version 3.1) --
4 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
5 * November 2006
6 *
7 * .. Scalar Arguments ..
8 INTEGER KNT, LMAX, NIN, NINFO
9 DOUBLE PRECISION RMAX
10 * ..
11 *
12 * Purpose
13 * =======
14 *
15 * ZGET35 tests ZTRSYL, a routine for solving the Sylvester matrix
16 * equation
17 *
18 * op(A)*X + ISGN*X*op(B) = scale*C,
19 *
20 * A and B are assumed to be in Schur canonical form, op() represents an
21 * optional transpose, and ISGN can be -1 or +1. Scale is an output
22 * less than or equal to 1, chosen to avoid overflow in X.
23 *
24 * The test code verifies that the following residual is order 1:
25 *
26 * norm(op(A)*X + ISGN*X*op(B) - scale*C) /
27 * (EPS*max(norm(A),norm(B))*norm(X))
28 *
29 * Arguments
30 * ==========
31 *
32 * RMAX (output) DOUBLE PRECISION
33 * Value of the largest test ratio.
34 *
35 * LMAX (output) INTEGER
36 * Example number where largest test ratio achieved.
37 *
38 * NINFO (output) INTEGER
39 * Number of examples where INFO is nonzero.
40 *
41 * KNT (output) INTEGER
42 * Total number of examples tested.
43 *
44 * NIN (input) INTEGER
45 * Input logical unit number.
46 *
47 * =====================================================================
48 *
49 * .. Parameters ..
50 INTEGER LDT
51 PARAMETER ( LDT = 10 )
52 DOUBLE PRECISION ZERO, ONE, TWO
53 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0, TWO = 2.0D0 )
54 DOUBLE PRECISION LARGE
55 PARAMETER ( LARGE = 1.0D6 )
56 COMPLEX*16 CONE
57 PARAMETER ( CONE = 1.0D0 )
58 * ..
59 * .. Local Scalars ..
60 CHARACTER TRANA, TRANB
61 INTEGER I, IMLA, IMLAD, IMLB, IMLC, INFO, ISGN, ITRANA,
62 $ ITRANB, J, M, N
63 DOUBLE PRECISION BIGNUM, EPS, RES, RES1, SCALE, SMLNUM, TNRM,
64 $ XNRM
65 COMPLEX*16 RMUL
66 * ..
67 * .. Local Arrays ..
68 DOUBLE PRECISION DUM( 1 ), VM1( 3 ), VM2( 3 )
69 COMPLEX*16 A( LDT, LDT ), ATMP( LDT, LDT ), B( LDT, LDT ),
70 $ BTMP( LDT, LDT ), C( LDT, LDT ),
71 $ CSAV( LDT, LDT ), CTMP( LDT, LDT )
72 * ..
73 * .. External Functions ..
74 DOUBLE PRECISION DLAMCH, ZLANGE
75 EXTERNAL DLAMCH, ZLANGE
76 * ..
77 * .. External Subroutines ..
78 EXTERNAL DLABAD, ZGEMM, ZTRSYL
79 * ..
80 * .. Intrinsic Functions ..
81 INTRINSIC ABS, DBLE, MAX, SQRT
82 * ..
83 * .. Executable Statements ..
84 *
85 * Get machine parameters
86 *
87 EPS = DLAMCH( 'P' )
88 SMLNUM = DLAMCH( 'S' ) / EPS
89 BIGNUM = ONE / SMLNUM
90 CALL DLABAD( SMLNUM, BIGNUM )
91 *
92 * Set up test case parameters
93 *
94 VM1( 1 ) = SQRT( SMLNUM )
95 VM1( 2 ) = ONE
96 VM1( 3 ) = LARGE
97 VM2( 1 ) = ONE
98 VM2( 2 ) = ONE + TWO*EPS
99 VM2( 3 ) = TWO
100 *
101 KNT = 0
102 NINFO = 0
103 LMAX = 0
104 RMAX = ZERO
105 *
106 * Begin test loop
107 *
108 10 CONTINUE
109 READ( NIN, FMT = * )M, N
110 IF( N.EQ.0 )
111 $ RETURN
112 DO 20 I = 1, M
113 READ( NIN, FMT = * )( ATMP( I, J ), J = 1, M )
114 20 CONTINUE
115 DO 30 I = 1, N
116 READ( NIN, FMT = * )( BTMP( I, J ), J = 1, N )
117 30 CONTINUE
118 DO 40 I = 1, M
119 READ( NIN, FMT = * )( CTMP( I, J ), J = 1, N )
120 40 CONTINUE
121 DO 170 IMLA = 1, 3
122 DO 160 IMLAD = 1, 3
123 DO 150 IMLB = 1, 3
124 DO 140 IMLC = 1, 3
125 DO 130 ITRANA = 1, 2
126 DO 120 ITRANB = 1, 2
127 DO 110 ISGN = -1, 1, 2
128 IF( ITRANA.EQ.1 )
129 $ TRANA = 'N'
130 IF( ITRANA.EQ.2 )
131 $ TRANA = 'C'
132 IF( ITRANB.EQ.1 )
133 $ TRANB = 'N'
134 IF( ITRANB.EQ.2 )
135 $ TRANB = 'C'
136 TNRM = ZERO
137 DO 60 I = 1, M
138 DO 50 J = 1, M
139 A( I, J ) = ATMP( I, J )*VM1( IMLA )
140 TNRM = MAX( TNRM, ABS( A( I, J ) ) )
141 50 CONTINUE
142 A( I, I ) = A( I, I )*VM2( IMLAD )
143 TNRM = MAX( TNRM, ABS( A( I, I ) ) )
144 60 CONTINUE
145 DO 80 I = 1, N
146 DO 70 J = 1, N
147 B( I, J ) = BTMP( I, J )*VM1( IMLB )
148 TNRM = MAX( TNRM, ABS( B( I, J ) ) )
149 70 CONTINUE
150 80 CONTINUE
151 IF( TNRM.EQ.ZERO )
152 $ TNRM = ONE
153 DO 100 I = 1, M
154 DO 90 J = 1, N
155 C( I, J ) = CTMP( I, J )*VM1( IMLC )
156 CSAV( I, J ) = C( I, J )
157 90 CONTINUE
158 100 CONTINUE
159 KNT = KNT + 1
160 CALL ZTRSYL( TRANA, TRANB, ISGN, M, N, A,
161 $ LDT, B, LDT, C, LDT, SCALE,
162 $ INFO )
163 IF( INFO.NE.0 )
164 $ NINFO = NINFO + 1
165 XNRM = ZLANGE( 'M', M, N, C, LDT, DUM )
166 RMUL = CONE
167 IF( XNRM.GT.ONE .AND. TNRM.GT.ONE ) THEN
168 IF( XNRM.GT.BIGNUM / TNRM ) THEN
169 RMUL = MAX( XNRM, TNRM )
170 RMUL = CONE / RMUL
171 END IF
172 END IF
173 CALL ZGEMM( TRANA, 'N', M, N, M, RMUL, A,
174 $ LDT, C, LDT, -SCALE*RMUL, CSAV,
175 $ LDT )
176 CALL ZGEMM( 'N', TRANB, M, N, N,
177 $ DBLE( ISGN )*RMUL, C, LDT, B,
178 $ LDT, CONE, CSAV, LDT )
179 RES1 = ZLANGE( 'M', M, N, CSAV, LDT, DUM )
180 RES = RES1 / MAX( SMLNUM, SMLNUM*XNRM,
181 $ ( ( ABS( RMUL )*TNRM )*EPS )*XNRM )
182 IF( RES.GT.RMAX ) THEN
183 LMAX = KNT
184 RMAX = RES
185 END IF
186 110 CONTINUE
187 120 CONTINUE
188 130 CONTINUE
189 140 CONTINUE
190 150 CONTINUE
191 160 CONTINUE
192 170 CONTINUE
193 GO TO 10
194 *
195 * End of ZGET35
196 *
197 END
2 *
3 * -- LAPACK test routine (version 3.1) --
4 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
5 * November 2006
6 *
7 * .. Scalar Arguments ..
8 INTEGER KNT, LMAX, NIN, NINFO
9 DOUBLE PRECISION RMAX
10 * ..
11 *
12 * Purpose
13 * =======
14 *
15 * ZGET35 tests ZTRSYL, a routine for solving the Sylvester matrix
16 * equation
17 *
18 * op(A)*X + ISGN*X*op(B) = scale*C,
19 *
20 * A and B are assumed to be in Schur canonical form, op() represents an
21 * optional transpose, and ISGN can be -1 or +1. Scale is an output
22 * less than or equal to 1, chosen to avoid overflow in X.
23 *
24 * The test code verifies that the following residual is order 1:
25 *
26 * norm(op(A)*X + ISGN*X*op(B) - scale*C) /
27 * (EPS*max(norm(A),norm(B))*norm(X))
28 *
29 * Arguments
30 * ==========
31 *
32 * RMAX (output) DOUBLE PRECISION
33 * Value of the largest test ratio.
34 *
35 * LMAX (output) INTEGER
36 * Example number where largest test ratio achieved.
37 *
38 * NINFO (output) INTEGER
39 * Number of examples where INFO is nonzero.
40 *
41 * KNT (output) INTEGER
42 * Total number of examples tested.
43 *
44 * NIN (input) INTEGER
45 * Input logical unit number.
46 *
47 * =====================================================================
48 *
49 * .. Parameters ..
50 INTEGER LDT
51 PARAMETER ( LDT = 10 )
52 DOUBLE PRECISION ZERO, ONE, TWO
53 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0, TWO = 2.0D0 )
54 DOUBLE PRECISION LARGE
55 PARAMETER ( LARGE = 1.0D6 )
56 COMPLEX*16 CONE
57 PARAMETER ( CONE = 1.0D0 )
58 * ..
59 * .. Local Scalars ..
60 CHARACTER TRANA, TRANB
61 INTEGER I, IMLA, IMLAD, IMLB, IMLC, INFO, ISGN, ITRANA,
62 $ ITRANB, J, M, N
63 DOUBLE PRECISION BIGNUM, EPS, RES, RES1, SCALE, SMLNUM, TNRM,
64 $ XNRM
65 COMPLEX*16 RMUL
66 * ..
67 * .. Local Arrays ..
68 DOUBLE PRECISION DUM( 1 ), VM1( 3 ), VM2( 3 )
69 COMPLEX*16 A( LDT, LDT ), ATMP( LDT, LDT ), B( LDT, LDT ),
70 $ BTMP( LDT, LDT ), C( LDT, LDT ),
71 $ CSAV( LDT, LDT ), CTMP( LDT, LDT )
72 * ..
73 * .. External Functions ..
74 DOUBLE PRECISION DLAMCH, ZLANGE
75 EXTERNAL DLAMCH, ZLANGE
76 * ..
77 * .. External Subroutines ..
78 EXTERNAL DLABAD, ZGEMM, ZTRSYL
79 * ..
80 * .. Intrinsic Functions ..
81 INTRINSIC ABS, DBLE, MAX, SQRT
82 * ..
83 * .. Executable Statements ..
84 *
85 * Get machine parameters
86 *
87 EPS = DLAMCH( 'P' )
88 SMLNUM = DLAMCH( 'S' ) / EPS
89 BIGNUM = ONE / SMLNUM
90 CALL DLABAD( SMLNUM, BIGNUM )
91 *
92 * Set up test case parameters
93 *
94 VM1( 1 ) = SQRT( SMLNUM )
95 VM1( 2 ) = ONE
96 VM1( 3 ) = LARGE
97 VM2( 1 ) = ONE
98 VM2( 2 ) = ONE + TWO*EPS
99 VM2( 3 ) = TWO
100 *
101 KNT = 0
102 NINFO = 0
103 LMAX = 0
104 RMAX = ZERO
105 *
106 * Begin test loop
107 *
108 10 CONTINUE
109 READ( NIN, FMT = * )M, N
110 IF( N.EQ.0 )
111 $ RETURN
112 DO 20 I = 1, M
113 READ( NIN, FMT = * )( ATMP( I, J ), J = 1, M )
114 20 CONTINUE
115 DO 30 I = 1, N
116 READ( NIN, FMT = * )( BTMP( I, J ), J = 1, N )
117 30 CONTINUE
118 DO 40 I = 1, M
119 READ( NIN, FMT = * )( CTMP( I, J ), J = 1, N )
120 40 CONTINUE
121 DO 170 IMLA = 1, 3
122 DO 160 IMLAD = 1, 3
123 DO 150 IMLB = 1, 3
124 DO 140 IMLC = 1, 3
125 DO 130 ITRANA = 1, 2
126 DO 120 ITRANB = 1, 2
127 DO 110 ISGN = -1, 1, 2
128 IF( ITRANA.EQ.1 )
129 $ TRANA = 'N'
130 IF( ITRANA.EQ.2 )
131 $ TRANA = 'C'
132 IF( ITRANB.EQ.1 )
133 $ TRANB = 'N'
134 IF( ITRANB.EQ.2 )
135 $ TRANB = 'C'
136 TNRM = ZERO
137 DO 60 I = 1, M
138 DO 50 J = 1, M
139 A( I, J ) = ATMP( I, J )*VM1( IMLA )
140 TNRM = MAX( TNRM, ABS( A( I, J ) ) )
141 50 CONTINUE
142 A( I, I ) = A( I, I )*VM2( IMLAD )
143 TNRM = MAX( TNRM, ABS( A( I, I ) ) )
144 60 CONTINUE
145 DO 80 I = 1, N
146 DO 70 J = 1, N
147 B( I, J ) = BTMP( I, J )*VM1( IMLB )
148 TNRM = MAX( TNRM, ABS( B( I, J ) ) )
149 70 CONTINUE
150 80 CONTINUE
151 IF( TNRM.EQ.ZERO )
152 $ TNRM = ONE
153 DO 100 I = 1, M
154 DO 90 J = 1, N
155 C( I, J ) = CTMP( I, J )*VM1( IMLC )
156 CSAV( I, J ) = C( I, J )
157 90 CONTINUE
158 100 CONTINUE
159 KNT = KNT + 1
160 CALL ZTRSYL( TRANA, TRANB, ISGN, M, N, A,
161 $ LDT, B, LDT, C, LDT, SCALE,
162 $ INFO )
163 IF( INFO.NE.0 )
164 $ NINFO = NINFO + 1
165 XNRM = ZLANGE( 'M', M, N, C, LDT, DUM )
166 RMUL = CONE
167 IF( XNRM.GT.ONE .AND. TNRM.GT.ONE ) THEN
168 IF( XNRM.GT.BIGNUM / TNRM ) THEN
169 RMUL = MAX( XNRM, TNRM )
170 RMUL = CONE / RMUL
171 END IF
172 END IF
173 CALL ZGEMM( TRANA, 'N', M, N, M, RMUL, A,
174 $ LDT, C, LDT, -SCALE*RMUL, CSAV,
175 $ LDT )
176 CALL ZGEMM( 'N', TRANB, M, N, N,
177 $ DBLE( ISGN )*RMUL, C, LDT, B,
178 $ LDT, CONE, CSAV, LDT )
179 RES1 = ZLANGE( 'M', M, N, CSAV, LDT, DUM )
180 RES = RES1 / MAX( SMLNUM, SMLNUM*XNRM,
181 $ ( ( ABS( RMUL )*TNRM )*EPS )*XNRM )
182 IF( RES.GT.RMAX ) THEN
183 LMAX = KNT
184 RMAX = RES
185 END IF
186 110 CONTINUE
187 120 CONTINUE
188 130 CONTINUE
189 140 CONTINUE
190 150 CONTINUE
191 160 CONTINUE
192 170 CONTINUE
193 GO TO 10
194 *
195 * End of ZGET35
196 *
197 END