1 SUBROUTINE DGET36( 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
9 DOUBLE PRECISION RMAX
10 * ..
11 * .. Array Arguments ..
12 INTEGER NINFO( 3 )
13 * ..
14 *
15 * Purpose
16 * =======
17 *
18 * DGET36 tests DTREXC, a routine for moving blocks (either 1 by 1 or
19 * 2 by 2) on the diagonal of a matrix in real Schur form. Thus, DLAEXC
20 * computes an orthogonal matrix Q such that
21 *
22 * Q' * T1 * Q = T2
23 *
24 * and where one of the diagonal blocks of T1 (the one at row IFST) has
25 * been moved to position ILST.
26 *
27 * The test code verifies that the residual Q'*T1*Q-T2 is small, that T2
28 * is in Schur form, and that the final position of the IFST block is
29 * ILST (within +-1).
30 *
31 * The test matrices are read from a file with logical unit number NIN.
32 *
33 * Arguments
34 * ==========
35 *
36 * RMAX (output) DOUBLE PRECISION
37 * Value of the largest test ratio.
38 *
39 * LMAX (output) INTEGER
40 * Example number where largest test ratio achieved.
41 *
42 * NINFO (output) INTEGER array, dimension (3)
43 * NINFO(J) is the number of examples where INFO=J.
44 *
45 * KNT (output) INTEGER
46 * Total number of examples tested.
47 *
48 * NIN (input) INTEGER
49 * Input logical unit number.
50 *
51 * =====================================================================
52 *
53 * .. Parameters ..
54 DOUBLE PRECISION ZERO, ONE
55 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
56 INTEGER LDT, LWORK
57 PARAMETER ( LDT = 10, LWORK = 2*LDT*LDT )
58 * ..
59 * .. Local Scalars ..
60 INTEGER I, IFST, IFST1, IFST2, IFSTSV, ILST, ILST1,
61 $ ILST2, ILSTSV, INFO1, INFO2, J, LOC, N
62 DOUBLE PRECISION EPS, RES
63 * ..
64 * .. Local Arrays ..
65 DOUBLE PRECISION Q( LDT, LDT ), RESULT( 2 ), T1( LDT, LDT ),
66 $ T2( LDT, LDT ), TMP( LDT, LDT ), WORK( LWORK )
67 * ..
68 * .. External Functions ..
69 DOUBLE PRECISION DLAMCH
70 EXTERNAL DLAMCH
71 * ..
72 * .. External Subroutines ..
73 EXTERNAL DHST01, DLACPY, DLASET, DTREXC
74 * ..
75 * .. Intrinsic Functions ..
76 INTRINSIC ABS, SIGN
77 * ..
78 * .. Executable Statements ..
79 *
80 EPS = DLAMCH( 'P' )
81 RMAX = ZERO
82 LMAX = 0
83 KNT = 0
84 NINFO( 1 ) = 0
85 NINFO( 2 ) = 0
86 NINFO( 3 ) = 0
87 *
88 * Read input data until N=0
89 *
90 10 CONTINUE
91 READ( NIN, FMT = * )N, IFST, ILST
92 IF( N.EQ.0 )
93 $ RETURN
94 KNT = KNT + 1
95 DO 20 I = 1, N
96 READ( NIN, FMT = * )( TMP( I, J ), J = 1, N )
97 20 CONTINUE
98 CALL DLACPY( 'F', N, N, TMP, LDT, T1, LDT )
99 CALL DLACPY( 'F', N, N, TMP, LDT, T2, LDT )
100 IFSTSV = IFST
101 ILSTSV = ILST
102 IFST1 = IFST
103 ILST1 = ILST
104 IFST2 = IFST
105 ILST2 = ILST
106 RES = ZERO
107 *
108 * Test without accumulating Q
109 *
110 CALL DLASET( 'Full', N, N, ZERO, ONE, Q, LDT )
111 CALL DTREXC( 'N', N, T1, LDT, Q, LDT, IFST1, ILST1, WORK, INFO1 )
112 DO 40 I = 1, N
113 DO 30 J = 1, N
114 IF( I.EQ.J .AND. Q( I, J ).NE.ONE )
115 $ RES = RES + ONE / EPS
116 IF( I.NE.J .AND. Q( I, J ).NE.ZERO )
117 $ RES = RES + ONE / EPS
118 30 CONTINUE
119 40 CONTINUE
120 *
121 * Test with accumulating Q
122 *
123 CALL DLASET( 'Full', N, N, ZERO, ONE, Q, LDT )
124 CALL DTREXC( 'V', N, T2, LDT, Q, LDT, IFST2, ILST2, WORK, INFO2 )
125 *
126 * Compare T1 with T2
127 *
128 DO 60 I = 1, N
129 DO 50 J = 1, N
130 IF( T1( I, J ).NE.T2( I, J ) )
131 $ RES = RES + ONE / EPS
132 50 CONTINUE
133 60 CONTINUE
134 IF( IFST1.NE.IFST2 )
135 $ RES = RES + ONE / EPS
136 IF( ILST1.NE.ILST2 )
137 $ RES = RES + ONE / EPS
138 IF( INFO1.NE.INFO2 )
139 $ RES = RES + ONE / EPS
140 *
141 * Test for successful reordering of T2
142 *
143 IF( INFO2.NE.0 ) THEN
144 NINFO( INFO2 ) = NINFO( INFO2 ) + 1
145 ELSE
146 IF( ABS( IFST2-IFSTSV ).GT.1 )
147 $ RES = RES + ONE / EPS
148 IF( ABS( ILST2-ILSTSV ).GT.1 )
149 $ RES = RES + ONE / EPS
150 END IF
151 *
152 * Test for small residual, and orthogonality of Q
153 *
154 CALL DHST01( N, 1, N, TMP, LDT, T2, LDT, Q, LDT, WORK, LWORK,
155 $ RESULT )
156 RES = RES + RESULT( 1 ) + RESULT( 2 )
157 *
158 * Test for T2 being in Schur form
159 *
160 LOC = 1
161 70 CONTINUE
162 IF( T2( LOC+1, LOC ).NE.ZERO ) THEN
163 *
164 * 2 by 2 block
165 *
166 IF( T2( LOC, LOC+1 ).EQ.ZERO .OR. T2( LOC, LOC ).NE.
167 $ T2( LOC+1, LOC+1 ) .OR. SIGN( ONE, T2( LOC, LOC+1 ) ).EQ.
168 $ SIGN( ONE, T2( LOC+1, LOC ) ) )RES = RES + ONE / EPS
169 DO 80 I = LOC + 2, N
170 IF( T2( I, LOC ).NE.ZERO )
171 $ RES = RES + ONE / RES
172 IF( T2( I, LOC+1 ).NE.ZERO )
173 $ RES = RES + ONE / RES
174 80 CONTINUE
175 LOC = LOC + 2
176 ELSE
177 *
178 * 1 by 1 block
179 *
180 DO 90 I = LOC + 1, N
181 IF( T2( I, LOC ).NE.ZERO )
182 $ RES = RES + ONE / RES
183 90 CONTINUE
184 LOC = LOC + 1
185 END IF
186 IF( LOC.LT.N )
187 $ GO TO 70
188 IF( RES.GT.RMAX ) THEN
189 RMAX = RES
190 LMAX = KNT
191 END IF
192 GO TO 10
193 *
194 * End of DGET36
195 *
196 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
9 DOUBLE PRECISION RMAX
10 * ..
11 * .. Array Arguments ..
12 INTEGER NINFO( 3 )
13 * ..
14 *
15 * Purpose
16 * =======
17 *
18 * DGET36 tests DTREXC, a routine for moving blocks (either 1 by 1 or
19 * 2 by 2) on the diagonal of a matrix in real Schur form. Thus, DLAEXC
20 * computes an orthogonal matrix Q such that
21 *
22 * Q' * T1 * Q = T2
23 *
24 * and where one of the diagonal blocks of T1 (the one at row IFST) has
25 * been moved to position ILST.
26 *
27 * The test code verifies that the residual Q'*T1*Q-T2 is small, that T2
28 * is in Schur form, and that the final position of the IFST block is
29 * ILST (within +-1).
30 *
31 * The test matrices are read from a file with logical unit number NIN.
32 *
33 * Arguments
34 * ==========
35 *
36 * RMAX (output) DOUBLE PRECISION
37 * Value of the largest test ratio.
38 *
39 * LMAX (output) INTEGER
40 * Example number where largest test ratio achieved.
41 *
42 * NINFO (output) INTEGER array, dimension (3)
43 * NINFO(J) is the number of examples where INFO=J.
44 *
45 * KNT (output) INTEGER
46 * Total number of examples tested.
47 *
48 * NIN (input) INTEGER
49 * Input logical unit number.
50 *
51 * =====================================================================
52 *
53 * .. Parameters ..
54 DOUBLE PRECISION ZERO, ONE
55 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
56 INTEGER LDT, LWORK
57 PARAMETER ( LDT = 10, LWORK = 2*LDT*LDT )
58 * ..
59 * .. Local Scalars ..
60 INTEGER I, IFST, IFST1, IFST2, IFSTSV, ILST, ILST1,
61 $ ILST2, ILSTSV, INFO1, INFO2, J, LOC, N
62 DOUBLE PRECISION EPS, RES
63 * ..
64 * .. Local Arrays ..
65 DOUBLE PRECISION Q( LDT, LDT ), RESULT( 2 ), T1( LDT, LDT ),
66 $ T2( LDT, LDT ), TMP( LDT, LDT ), WORK( LWORK )
67 * ..
68 * .. External Functions ..
69 DOUBLE PRECISION DLAMCH
70 EXTERNAL DLAMCH
71 * ..
72 * .. External Subroutines ..
73 EXTERNAL DHST01, DLACPY, DLASET, DTREXC
74 * ..
75 * .. Intrinsic Functions ..
76 INTRINSIC ABS, SIGN
77 * ..
78 * .. Executable Statements ..
79 *
80 EPS = DLAMCH( 'P' )
81 RMAX = ZERO
82 LMAX = 0
83 KNT = 0
84 NINFO( 1 ) = 0
85 NINFO( 2 ) = 0
86 NINFO( 3 ) = 0
87 *
88 * Read input data until N=0
89 *
90 10 CONTINUE
91 READ( NIN, FMT = * )N, IFST, ILST
92 IF( N.EQ.0 )
93 $ RETURN
94 KNT = KNT + 1
95 DO 20 I = 1, N
96 READ( NIN, FMT = * )( TMP( I, J ), J = 1, N )
97 20 CONTINUE
98 CALL DLACPY( 'F', N, N, TMP, LDT, T1, LDT )
99 CALL DLACPY( 'F', N, N, TMP, LDT, T2, LDT )
100 IFSTSV = IFST
101 ILSTSV = ILST
102 IFST1 = IFST
103 ILST1 = ILST
104 IFST2 = IFST
105 ILST2 = ILST
106 RES = ZERO
107 *
108 * Test without accumulating Q
109 *
110 CALL DLASET( 'Full', N, N, ZERO, ONE, Q, LDT )
111 CALL DTREXC( 'N', N, T1, LDT, Q, LDT, IFST1, ILST1, WORK, INFO1 )
112 DO 40 I = 1, N
113 DO 30 J = 1, N
114 IF( I.EQ.J .AND. Q( I, J ).NE.ONE )
115 $ RES = RES + ONE / EPS
116 IF( I.NE.J .AND. Q( I, J ).NE.ZERO )
117 $ RES = RES + ONE / EPS
118 30 CONTINUE
119 40 CONTINUE
120 *
121 * Test with accumulating Q
122 *
123 CALL DLASET( 'Full', N, N, ZERO, ONE, Q, LDT )
124 CALL DTREXC( 'V', N, T2, LDT, Q, LDT, IFST2, ILST2, WORK, INFO2 )
125 *
126 * Compare T1 with T2
127 *
128 DO 60 I = 1, N
129 DO 50 J = 1, N
130 IF( T1( I, J ).NE.T2( I, J ) )
131 $ RES = RES + ONE / EPS
132 50 CONTINUE
133 60 CONTINUE
134 IF( IFST1.NE.IFST2 )
135 $ RES = RES + ONE / EPS
136 IF( ILST1.NE.ILST2 )
137 $ RES = RES + ONE / EPS
138 IF( INFO1.NE.INFO2 )
139 $ RES = RES + ONE / EPS
140 *
141 * Test for successful reordering of T2
142 *
143 IF( INFO2.NE.0 ) THEN
144 NINFO( INFO2 ) = NINFO( INFO2 ) + 1
145 ELSE
146 IF( ABS( IFST2-IFSTSV ).GT.1 )
147 $ RES = RES + ONE / EPS
148 IF( ABS( ILST2-ILSTSV ).GT.1 )
149 $ RES = RES + ONE / EPS
150 END IF
151 *
152 * Test for small residual, and orthogonality of Q
153 *
154 CALL DHST01( N, 1, N, TMP, LDT, T2, LDT, Q, LDT, WORK, LWORK,
155 $ RESULT )
156 RES = RES + RESULT( 1 ) + RESULT( 2 )
157 *
158 * Test for T2 being in Schur form
159 *
160 LOC = 1
161 70 CONTINUE
162 IF( T2( LOC+1, LOC ).NE.ZERO ) THEN
163 *
164 * 2 by 2 block
165 *
166 IF( T2( LOC, LOC+1 ).EQ.ZERO .OR. T2( LOC, LOC ).NE.
167 $ T2( LOC+1, LOC+1 ) .OR. SIGN( ONE, T2( LOC, LOC+1 ) ).EQ.
168 $ SIGN( ONE, T2( LOC+1, LOC ) ) )RES = RES + ONE / EPS
169 DO 80 I = LOC + 2, N
170 IF( T2( I, LOC ).NE.ZERO )
171 $ RES = RES + ONE / RES
172 IF( T2( I, LOC+1 ).NE.ZERO )
173 $ RES = RES + ONE / RES
174 80 CONTINUE
175 LOC = LOC + 2
176 ELSE
177 *
178 * 1 by 1 block
179 *
180 DO 90 I = LOC + 1, N
181 IF( T2( I, LOC ).NE.ZERO )
182 $ RES = RES + ONE / RES
183 90 CONTINUE
184 LOC = LOC + 1
185 END IF
186 IF( LOC.LT.N )
187 $ GO TO 70
188 IF( RES.GT.RMAX ) THEN
189 RMAX = RES
190 LMAX = KNT
191 END IF
192 GO TO 10
193 *
194 * End of DGET36
195 *
196 END