1*> \brief \b ZLATM6
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 ZLATM6( TYPE, N, A, LDA, B, X, LDX, Y, LDY, ALPHA,
12*                          BETA, WX, WY, S, DIF )
13*
14*       .. Scalar Arguments ..
15*       INTEGER            LDA, LDX, LDY, N, TYPE
16*       COMPLEX*16         ALPHA, BETA, WX, WY
17*       ..
18*       .. Array Arguments ..
19*       DOUBLE PRECISION   DIF( * ), S( * )
20*       COMPLEX*16         A( LDA, * ), B( LDA, * ), X( LDX, * ),
21*      $                   Y( LDY, * )
22*       ..
23*
24*
25*> \par Purpose:
26*  =============
27*>
28*> \verbatim
29*>
30*> ZLATM6 generates test matrices for the generalized eigenvalue
31*> problem, their corresponding right and left eigenvector matrices,
32*> and also reciprocal condition numbers for all eigenvalues and
33*> the reciprocal condition numbers of eigenvectors corresponding to
34*> the 1th and 5th eigenvalues.
35*>
36*> Test Matrices
37*> =============
38*>
39*> Two kinds of test matrix pairs
40*>          (A, B) = inverse(YH) * (Da, Db) * inverse(X)
41*> are used in the tests:
42*>
43*> Type 1:
44*>    Da = 1+a   0    0    0    0    Db = 1   0   0   0   0
45*>          0   2+a   0    0    0         0   1   0   0   0
46*>          0    0   3+a   0    0         0   0   1   0   0
47*>          0    0    0   4+a   0         0   0   0   1   0
48*>          0    0    0    0   5+a ,      0   0   0   0   1
49*> and Type 2:
50*>    Da = 1+i   0    0       0       0    Db = 1   0   0   0   0
51*>          0   1-i   0       0       0         0   1   0   0   0
52*>          0    0    1       0       0         0   0   1   0   0
53*>          0    0    0 (1+a)+(1+b)i  0         0   0   0   1   0
54*>          0    0    0       0 (1+a)-(1+b)i,   0   0   0   0   1 .
55*>
56*> In both cases the same inverse(YH) and inverse(X) are used to compute
57*> (A, B), giving the exact eigenvectors to (A,B) as (YH, X):
58*>
59*> YH:  =  1    0   -y    y   -y    X =  1   0  -x  -x   x
60*>         0    1   -y    y   -y         0   1   x  -x  -x
61*>         0    0    1    0    0         0   0   1   0   0
62*>         0    0    0    1    0         0   0   0   1   0
63*>         0    0    0    0    1,        0   0   0   0   1 , where
64*>
65*> a, b, x and y will have all values independently of each other.
66*> \endverbatim
67*
68*  Arguments:
69*  ==========
70*
71*> \param[in] TYPE
72*> \verbatim
73*>          TYPE is INTEGER
74*>          Specifies the problem type (see futher details).
75*> \endverbatim
76*>
77*> \param[in] N
78*> \verbatim
79*>          N is INTEGER
80*>          Size of the matrices A and B.
81*> \endverbatim
82*>
83*> \param[out] A
84*> \verbatim
85*>          A is COMPLEX*16 array, dimension (LDA, N).
86*>          On exit A N-by-N is initialized according to TYPE.
87*> \endverbatim
88*>
89*> \param[in] LDA
90*> \verbatim
91*>          LDA is INTEGER
92*>          The leading dimension of A and of B.
93*> \endverbatim
94*>
95*> \param[out] B
96*> \verbatim
97*>          B is COMPLEX*16 array, dimension (LDA, N).
98*>          On exit B N-by-N is initialized according to TYPE.
99*> \endverbatim
100*>
101*> \param[out] X
102*> \verbatim
103*>          X is COMPLEX*16 array, dimension (LDX, N).
104*>          On exit X is the N-by-N matrix of right eigenvectors.
105*> \endverbatim
106*>
107*> \param[in] LDX
108*> \verbatim
109*>          LDX is INTEGER
110*>          The leading dimension of X.
111*> \endverbatim
112*>
113*> \param[out] Y
114*> \verbatim
115*>          Y is COMPLEX*16 array, dimension (LDY, N).
116*>          On exit Y is the N-by-N matrix of left eigenvectors.
117*> \endverbatim
118*>
119*> \param[in] LDY
120*> \verbatim
121*>          LDY is INTEGER
122*>          The leading dimension of Y.
123*> \endverbatim
124*>
125*> \param[in] ALPHA
126*> \verbatim
127*>          ALPHA is COMPLEX*16
128*> \endverbatim
129*>
130*> \param[in] BETA
131*> \verbatim
132*>          BETA is COMPLEX*16
133*> \verbatim
134*>          Weighting constants for matrix A.
135*> \endverbatim
136*>
137*> \param[in] WX
138*> \verbatim
139*>          WX is COMPLEX*16
140*>          Constant for right eigenvector matrix.
141*> \endverbatim
142*>
143*> \param[in] WY
144*> \verbatim
145*>          WY is COMPLEX*16
146*>          Constant for left eigenvector matrix.
147*> \endverbatim
148*>
149*> \param[out] S
150*> \verbatim
151*>          S is DOUBLE PRECISION array, dimension (N)
152*>          S(i) is the reciprocal condition number for eigenvalue i.
153*> \endverbatim
154*>
155*> \param[out] DIF
156*> \verbatim
157*>          DIF is DOUBLE PRECISION array, dimension (N)
158*>          DIF(i) is the reciprocal condition number for eigenvector i.
159*> \endverbatim
160*
161*  Authors:
162*  ========
163*
164*> \author Univ. of Tennessee
165*> \author Univ. of California Berkeley
166*> \author Univ. of Colorado Denver
167*> \author NAG Ltd.
168*
169*> \date November 2011
170*
171*> \ingroup complex16_matgen
172*
173*  =====================================================================
174      SUBROUTINE ZLATM6( TYPE, N, A, LDA, B, X, LDX, Y, LDY, ALPHA,
175     $                   BETA, WX, WY, S, DIF )
176*
177*  -- LAPACK computational routine (version 3.4.0) --
178*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
179*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
180*     November 2011
181*
182*     .. Scalar Arguments ..
183      INTEGER            LDA, LDX, LDY, N, TYPE
184      COMPLEX*16         ALPHA, BETA, WX, WY
185*     ..
186*     .. Array Arguments ..
187      DOUBLE PRECISION   DIF( * ), S( * )
188      COMPLEX*16         A( LDA, * ), B( LDA, * ), X( LDX, * ),
189     $                   Y( LDY, * )
190*     ..
191*
192*  =====================================================================
193*
194*     .. Parameters ..
195      DOUBLE PRECISION   RONE, TWO, THREE
196      PARAMETER          ( RONE = 1.0D+0, TWO = 2.0D+0, THREE = 3.0D+0 )
197      COMPLEX*16         ZERO, ONE
198      PARAMETER          ( ZERO = ( 0.0D+0, 0.0D+0 ),
199     $                   ONE = ( 1.0D+0, 0.0D+0 ) )
200*     ..
201*     .. Local Scalars ..
202      INTEGER            I, INFO, J
203*     ..
204*     .. Local Arrays ..
205      DOUBLE PRECISION   RWORK( 50 )
206      COMPLEX*16         WORK( 26 ), Z( 8, 8 )
207*     ..
208*     .. Intrinsic Functions ..
209      INTRINSIC          CDABS, DBLE, DCMPLX, DCONJG, SQRT
210*     ..
211*     .. External Subroutines ..
212      EXTERNAL           ZGESVD, ZLACPY, ZLAKF2
213*     ..
214*     .. Executable Statements ..
215*
216*     Generate test problem ...
217*     (Da, Db) ...
218*
219      DO 20 I = 1, N
220         DO 10 J = 1, N
221*
222            IF( I.EQ.J ) THEN
223               A( I, I ) = DCMPLX( I ) + ALPHA
224               B( I, I ) = ONE
225            ELSE
226               A( I, J ) = ZERO
227               B( I, J ) = ZERO
228            END IF
229*
230   10    CONTINUE
231   20 CONTINUE
232      IF( TYPE.EQ.2 ) THEN
233         A( 1, 1 ) = DCMPLX( RONE, RONE )
234         A( 2, 2 ) = DCONJG( A( 1, 1 ) )
235         A( 3, 3 ) = ONE
236         A( 4, 4 ) = DCMPLX( DBLE( ONE+ALPHA ), DBLE( ONE+BETA ) )
237         A( 5, 5 ) = DCONJG( A( 4, 4 ) )
238      END IF
239*
240*     Form X and Y
241*
242      CALL ZLACPY( 'F', N, N, B, LDA, Y, LDY )
243      Y( 3, 1 ) = -DCONJG( WY )
244      Y( 4, 1 ) = DCONJG( WY )
245      Y( 5, 1 ) = -DCONJG( WY )
246      Y( 3, 2 ) = -DCONJG( WY )
247      Y( 4, 2 ) = DCONJG( WY )
248      Y( 5, 2 ) = -DCONJG( WY )
249*
250      CALL ZLACPY( 'F', N, N, B, LDA, X, LDX )
251      X( 1, 3 ) = -WX
252      X( 1, 4 ) = -WX
253      X( 1, 5 ) = WX
254      X( 2, 3 ) = WX
255      X( 2, 4 ) = -WX
256      X( 2, 5 ) = -WX
257*
258*     Form (A, B)
259*
260      B( 1, 3 ) = WX + WY
261      B( 2, 3 ) = -WX + WY
262      B( 1, 4 ) = WX - WY
263      B( 2, 4 ) = WX - WY
264      B( 1, 5 ) = -WX + WY
265      B( 2, 5 ) = WX + WY
266      A( 1, 3 ) = WX*A( 1, 1 ) + WY*A( 3, 3 )
267      A( 2, 3 ) = -WX*A( 2, 2 ) + WY*A( 3, 3 )
268      A( 1, 4 ) = WX*A( 1, 1 ) - WY*A( 4, 4 )
269      A( 2, 4 ) = WX*A( 2, 2 ) - WY*A( 4, 4 )
270      A( 1, 5 ) = -WX*A( 1, 1 ) + WY*A( 5, 5 )
271      A( 2, 5 ) = WX*A( 2, 2 ) + WY*A( 5, 5 )
272*
273*     Compute condition numbers
274*
275      S( 1 ) = RONE / SQRT( ( RONE+THREE*CDABS( WY )*CDABS( WY ) ) /
276     $         ( RONE+CDABS( A( 1, 1 ) )*CDABS( A( 1, 1 ) ) ) )
277      S( 2 ) = RONE / SQRT( ( RONE+THREE*CDABS( WY )*CDABS( WY ) ) /
278     $         ( RONE+CDABS( A( 2, 2 ) )*CDABS( A( 2, 2 ) ) ) )
279      S( 3 ) = RONE / SQRT( ( RONE+TWO*CDABS( WX )*CDABS( WX ) ) /
280     $         ( RONE+CDABS( A( 3, 3 ) )*CDABS( A( 3, 3 ) ) ) )
281      S( 4 ) = RONE / SQRT( ( RONE+TWO*CDABS( WX )*CDABS( WX ) ) /
282     $         ( RONE+CDABS( A( 4, 4 ) )*CDABS( A( 4, 4 ) ) ) )
283      S( 5 ) = RONE / SQRT( ( RONE+TWO*CDABS( WX )*CDABS( WX ) ) /
284     $         ( RONE+CDABS( A( 5, 5 ) )*CDABS( A( 5, 5 ) ) ) )
285*
286      CALL ZLAKF2( 1, 4, A, LDA, A( 2, 2 ), B, B( 2, 2 ), Z, 8 )
287      CALL ZGESVD( 'N', 'N', 8, 8, Z, 8, RWORK, WORK, 1, WORK( 2 ), 1,
288     $             WORK( 3 ), 24, RWORK( 9 ), INFO )
289      DIF( 1 ) = RWORK( 8 )
290*
291      CALL ZLAKF2( 4, 1, A, LDA, A( 5, 5 ), B, B( 5, 5 ), Z, 8 )
292      CALL ZGESVD( 'N', 'N', 8, 8, Z, 8, RWORK, WORK, 1, WORK( 2 ), 1,
293     $             WORK( 3 ), 24, RWORK( 9 ), INFO )
294      DIF( 5 ) = RWORK( 8 )
295*
296      RETURN
297*
298*     End of ZLATM6
299*
300      END
301