LAPACK  3.10.1 LAPACK: Linear Algebra PACKage
dget40.f
Go to the documentation of this file.
1 *> \brief \b DGET40
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 * Definition:
9 * ===========
10 *
11 * SUBROUTINE DGET40( RMAX, LMAX, NINFO, KNT, NIN )
12 *
13 * .. Scalar Arguments ..
14 * INTEGER KNT, LMAX, NIN
15 * DOUBLE PRECISION RMAX
16 * ..
17 * .. Array Arguments ..
18 * INTEGER NINFO( 3 )
19 *
20 *
21 *> \par Purpose:
22 * =============
23 *>
24 *> \verbatim
25 *>
26 *> DGET40 tests DTGEXC, a routine for swapping adjacent blocks (either
27 *> 1 by 1 or 2 by 2) on the diagonal of a pencil in real generalized Schur form.
28 *> Thus, DTGEXC computes an orthogonal matrices Q and Z such that
29 *>
30 *> Q' * ( [ A B ], [ D E ] ) * Z = ( [ C1 B1 ], [ F1 E1 ] )
31 *> ( [ 0 C ] [ F ] ) ( [ 0 A1 ] [ D1] )
32 *>
33 *> where (C1,F1) is similar to (C,F) and (A1,D1) is similar to (A,D).
34 *> Both (A,D) and (C,F) are assumed to be in standard form
35 *> and (A1,D1) and (C1,F1) are returned with the
36 *> same properties.
37 *> \endverbatim
38 *
39 * Arguments:
40 * ==========
41 *
42 *> \param[out] RMAX
43 *> \verbatim
44 *> RMAX is DOUBLE PRECISION
45 *> Value of the largest test ratio.
46 *> \endverbatim
47 *>
48 *> \param[out] LMAX
49 *> \verbatim
50 *> LMAX is INTEGER
51 *> Example number where largest test ratio achieved.
52 *> \endverbatim
53 *>
54 *> \param[out] NINFO
55 *> \verbatim
56 *> NINFO is INTEGER(3)
57 *> Number of examples where INFO is nonzero.
58 *> \endverbatim
59 *>
60 *> \param[out] KNT
61 *> \verbatim
62 *> KNT is INTEGER
63 *> Total number of examples tested.
64 *> \endverbatim
65 *>
66 *> \param[out] NIN
67 *> \verbatim
68 *> NINFO is INTEGER
69 *> \endverbatim
70 *
71 * Authors:
72 * ========
73 *
74 *> \author Univ. of Tennessee
75 *> \author Univ. of California Berkeley
76 *> \author Univ. of Colorado Denver
77 *> \author NAG Ltd.
78 *
79 *> \ingroup double_eig
80 *
81 * =====================================================================
82  SUBROUTINE dget40( RMAX, LMAX, NINFO, KNT, NIN )
83 *
84 * -- LAPACK test routine --
85 * -- LAPACK is a software package provided by Univ. of Tennessee, --
86 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
87 *
88 * .. Scalar Arguments ..
89  INTEGER KNT, LMAX, NIN
90  DOUBLE PRECISION RMAX
91 * ..
92 * .. Array Arguments ..
93  INTEGER NINFO( 3 )
94 * ..
95 *
96 * =====================================================================
97 *
98 * .. Parameters ..
99  DOUBLE PRECISION ZERO, ONE
100  parameter( zero = 0.0d0, one = 1.0d0 )
101  INTEGER LDT, LWORK
102  parameter( ldt = 10, lwork = 100 + 4*ldt + 16 )
103 * ..
104 * .. Local Scalars ..
105  INTEGER I, IFST, IFST1, IFST2, IFSTSV, ILST, ILST1,
106  \$ ILST2, ILSTSV, INFO1, INFO2, J, LOC, N
107  DOUBLE PRECISION EPS, RES
108 * ..
109 * .. Local Arrays ..
110  DOUBLE PRECISION Q( LDT, LDT ), Z( LDT, LDT ), RESULT( 4 ),
111  \$ T( LDT, LDT ), T1( LDT, LDT ), T2( LDT, LDT ),
112  \$ S( LDT, LDT ), S1( LDT, LDT ), S2( LDT, LDT ),
113  \$ TMP( LDT, LDT ), WORK( LWORK )
114 * ..
115 * .. External Functions ..
116  DOUBLE PRECISION DLAMCH
117  EXTERNAL dlamch
118 * ..
119 * .. External Subroutines ..
120  EXTERNAL dhst01, dlacpy, dlaset, dtgexc
121 * ..
122 * .. Intrinsic Functions ..
123  INTRINSIC abs, sign
124 * ..
125 * .. Executable Statements ..
126 *
127  eps = dlamch( 'P' )
128  rmax = zero
129  lmax = 0
130  knt = 0
131  ninfo( 1 ) = 0
132  ninfo( 2 ) = 0
133  ninfo( 3 ) = 0
134 *
135 * Read input data until N=0
136 *
137  10 CONTINUE
138  READ( nin, fmt = * )n, ifst, ilst
139  IF( n.EQ.0 )
140  \$ RETURN
141  knt = knt + 1
142  DO 20 i = 1, n
143  READ( nin, fmt = * )( tmp( i, j ), j = 1, n )
144  20 CONTINUE
145  CALL dlacpy( 'F', n, n, tmp, ldt, t, ldt )
146  CALL dlacpy( 'F', n, n, tmp, ldt, t1, ldt )
147  CALL dlacpy( 'F', n, n, tmp, ldt, t2, ldt )
148  DO 25 i = 1, n
149  READ( nin, fmt = * )( tmp( i, j ), j = 1, n )
150  25 CONTINUE
151  CALL dlacpy( 'F', n, n, tmp, ldt, s, ldt )
152  CALL dlacpy( 'F', n, n, tmp, ldt, s1, ldt )
153  CALL dlacpy( 'F', n, n, tmp, ldt, s2, ldt )
154  ifstsv = ifst
155  ilstsv = ilst
156  ifst1 = ifst
157  ilst1 = ilst
158  ifst2 = ifst
159  ilst2 = ilst
160  res = zero
161 *
162 * Test without accumulating Q and Z
163 *
164  CALL dlaset( 'Full', n, n, zero, one, q, ldt )
165  CALL dlaset( 'Full', n, n, zero, one, z, ldt )
166  CALL dtgexc( .false., .false., n, t1, ldt, s1, ldt, q, ldt,
167  \$ z, ldt, ifst1, ilst1, work, lwork, info1 )
168  DO 40 i = 1, n
169  DO 30 j = 1, n
170  IF( i.EQ.j .AND. q( i, j ).NE.one )
171  \$ res = res + one / eps
172  IF( i.NE.j .AND. q( i, j ).NE.zero )
173  \$ res = res + one / eps
174  IF( i.EQ.j .AND. z( i, j ).NE.one )
175  \$ res = res + one / eps
176  IF( i.NE.j .AND. z( i, j ).NE.zero )
177  \$ res = res + one / eps
178  30 CONTINUE
179  40 CONTINUE
180 *
181 * Test with accumulating Q
182 *
183  CALL dlaset( 'Full', n, n, zero, one, q, ldt )
184  CALL dlaset( 'Full', n, n, zero, one, z, ldt )
185  CALL dtgexc( .true., .true., n, t2, ldt, s2, ldt, q, ldt,
186  \$ z, ldt, ifst2, ilst2, work, lwork, info2 )
187 *
188 * Compare T1 with T2 and S1 with S2
189 *
190  DO 60 i = 1, n
191  DO 50 j = 1, n
192  IF( t1( i, j ).NE.t2( i, j ) )
193  \$ res = res + one / eps
194  IF( s1( i, j ).NE.s2( i, j ) )
195  \$ res = res + one / eps
196  50 CONTINUE
197  60 CONTINUE
198  IF( ifst1.NE.ifst2 )
199  \$ res = res + one / eps
200  IF( ilst1.NE.ilst2 )
201  \$ res = res + one / eps
202  IF( info1.NE.info2 )
203  \$ res = res + one / eps
204 *
205 * Test orthogonality of Q and Z and backward error on T2 and S2
206 *
207  CALL dget51( 1, n, t, ldt, t2, ldt, q, ldt, z, ldt, work,
208  \$ result( 1 ) )
209  CALL dget51( 1, n, s, ldt, s2, ldt, q, ldt, z, ldt, work,
210  \$ result( 2 ) )
211  CALL dget51( 3, n, t, ldt, t2, ldt, q, ldt, q, ldt, work,
212  \$ result( 3 ) )
213  CALL dget51( 3, n, t, ldt, t2, ldt, z, ldt, z, ldt, work,
214  \$ result( 4 ) )
215  res = res + result( 1 ) + result( 2 ) + result( 3 ) + result( 4 )
216 *
217 * Read next matrix pair
218 *
219  GO TO 10
220 *
221 * End of DGET40
222 *
223  END
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
Definition: dlacpy.f:103
subroutine dlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: dlaset.f:110
subroutine dhst01(N, ILO, IHI, A, LDA, H, LDH, Q, LDQ, WORK, LWORK, RESULT)
DHST01
Definition: dhst01.f:134
subroutine dget40(RMAX, LMAX, NINFO, KNT, NIN)
DGET40
Definition: dget40.f:83
subroutine dget51(ITYPE, N, A, LDA, B, LDB, U, LDU, V, LDV, WORK, RESULT)
DGET51
Definition: dget51.f:149
subroutine dtgexc(WANTQ, WANTZ, N, A, LDA, B, LDB, Q, LDQ, Z, LDZ, IFST, ILST, WORK, LWORK, INFO)
DTGEXC
Definition: dtgexc.f:220