1*> \brief \b DLAQZ2
2*
3*  =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6*            http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download DLAQZ2 + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlaqz2.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlaqz2.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlaqz2.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18*  Definition:
19*  ===========
20*
21*      SUBROUTINE DLAQZ2( ILQ, ILZ, K, ISTARTM, ISTOPM, IHI, A, LDA, B,
22*     $    LDB, NQ, QSTART, Q, LDQ, NZ, ZSTART, Z, LDZ )
23*      IMPLICIT NONE
24*
25*      Arguments
26*      LOGICAL, INTENT( IN ) :: ILQ, ILZ
27*      INTEGER, INTENT( IN ) :: K, LDA, LDB, LDQ, LDZ, ISTARTM, ISTOPM,
28*     $    NQ, NZ, QSTART, ZSTART, IHI
29*      DOUBLE PRECISION :: A( LDA, * ), B( LDB, * ), Q( LDQ, * ), Z( LDZ,
30*     $    * )
31*      ..
32*
33*
34*> \par Purpose:
35*  =============
36*>
37*> \verbatim
38*>
39*>      DLAQZ2 chases a 2x2 shift bulge in a matrix pencil down a single position
40*> \endverbatim
41*
42*
43*  Arguments:
44*  ==========
45*
46*>
47*> \param[in] ILQ
48*> \verbatim
49*>          ILQ is LOGICAL
50*>              Determines whether or not to update the matrix Q
51*> \endverbatim
52*>
53*> \param[in] ILZ
54*> \verbatim
55*>          ILZ is LOGICAL
56*>              Determines whether or not to update the matrix Z
57*> \endverbatim
58*>
59*> \param[in] K
60*> \verbatim
61*>          K is INTEGER
62*>              Index indicating the position of the bulge.
63*>              On entry, the bulge is located in
64*>              (A(k+1:k+2,k:k+1),B(k+1:k+2,k:k+1)).
65*>              On exit, the bulge is located in
66*>              (A(k+2:k+3,k+1:k+2),B(k+2:k+3,k+1:k+2)).
67*> \endverbatim
68*>
69*> \param[in] ISTARTM
70*> \verbatim
71*>          ISTARTM is INTEGER
72*> \endverbatim
73*>
74*> \param[in] ISTOPM
75*> \verbatim
76*>          ISTOPM is INTEGER
77*>              Updates to (A,B) are restricted to
78*>              (istartm:k+3,k:istopm). It is assumed
79*>              without checking that istartm <= k+1 and
80*>              k+2 <= istopm
81*> \endverbatim
82*>
83*> \param[in] IHI
84*> \verbatim
85*>          IHI is INTEGER
86*> \endverbatim
87*>
88*> \param[inout] A
89*> \verbatim
90*>          A is DOUBLE PRECISION array, dimension (LDA,N)
91*> \endverbatim
92*>
93*> \param[in] LDA
94*> \verbatim
95*>          LDA is INTEGER
96*>              The leading dimension of A as declared in
97*>              the calling procedure.
98*> \endverbatim
99*
100*> \param[inout] B
101*> \verbatim
102*>          B is DOUBLE PRECISION array, dimension (LDB,N)
103*> \endverbatim
104*>
105*> \param[in] LDB
106*> \verbatim
107*>          LDB is INTEGER
108*>              The leading dimension of B as declared in
109*>              the calling procedure.
110*> \endverbatim
111*>
112*> \param[in] NQ
113*> \verbatim
114*>          NQ is INTEGER
115*>              The order of the matrix Q
116*> \endverbatim
117*>
118*> \param[in] QSTART
119*> \verbatim
120*>          QSTART is INTEGER
121*>              Start index of the matrix Q. Rotations are applied
122*>              To columns k+2-qStart:k+4-qStart of Q.
123*> \endverbatim
124*
125*> \param[inout] Q
126*> \verbatim
127*>          Q is DOUBLE PRECISION array, dimension (LDQ,NQ)
128*> \endverbatim
129*>
130*> \param[in] LDQ
131*> \verbatim
132*>          LDQ is INTEGER
133*>              The leading dimension of Q as declared in
134*>              the calling procedure.
135*> \endverbatim
136*>
137*> \param[in] NZ
138*> \verbatim
139*>          NZ is INTEGER
140*>              The order of the matrix Z
141*> \endverbatim
142*>
143*> \param[in] ZSTART
144*> \verbatim
145*>          ZSTART is INTEGER
146*>              Start index of the matrix Z. Rotations are applied
147*>              To columns k+1-qStart:k+3-qStart of Z.
148*> \endverbatim
149*
150*> \param[inout] Z
151*> \verbatim
152*>          Z is DOUBLE PRECISION array, dimension (LDZ,NZ)
153*> \endverbatim
154*>
155*> \param[in] LDZ
156*> \verbatim
157*>          LDZ is INTEGER
158*>              The leading dimension of Q as declared in
159*>              the calling procedure.
160*> \endverbatim
161*
162*  Authors:
163*  ========
164*
165*> \author Thijs Steel, KU Leuven
166*
167*> \date May 2020
168*
169*> \ingroup doubleGEcomputational
170*>
171*  =====================================================================
172      SUBROUTINE DLAQZ2( ILQ, ILZ, K, ISTARTM, ISTOPM, IHI, A, LDA, B,
173     $                   LDB, NQ, QSTART, Q, LDQ, NZ, ZSTART, Z, LDZ )
174      IMPLICIT NONE
175*
176*     Arguments
177      LOGICAL, INTENT( IN ) :: ILQ, ILZ
178      INTEGER, INTENT( IN ) :: K, LDA, LDB, LDQ, LDZ, ISTARTM, ISTOPM,
179     $         NQ, NZ, QSTART, ZSTART, IHI
180      DOUBLE PRECISION :: A( LDA, * ), B( LDB, * ), Q( LDQ, * ), Z( LDZ,
181     $                    * )
182*
183*     Parameters
184      DOUBLE PRECISION :: ZERO, ONE, HALF
185      PARAMETER( ZERO = 0.0D0, ONE = 1.0D0, HALF = 0.5D0 )
186*
187*     Local variables
188      DOUBLE PRECISION :: H( 2, 3 ), C1, S1, C2, S2, TEMP
189*
190*     External functions
191      EXTERNAL :: DLARTG, DROT
192*
193      IF( K+2 .EQ. IHI ) THEN
194*        Shift is located on the edge of the matrix, remove it
195         H = B( IHI-1:IHI, IHI-2:IHI )
196*        Make H upper triangular
197         CALL DLARTG( H( 1, 1 ), H( 2, 1 ), C1, S1, TEMP )
198         H( 2, 1 ) = ZERO
199         H( 1, 1 ) = TEMP
200         CALL DROT( 2, H( 1, 2 ), 2, H( 2, 2 ), 2, C1, S1 )
201*
202         CALL DLARTG( H( 2, 3 ), H( 2, 2 ), C1, S1, TEMP )
203         CALL DROT( 1, H( 1, 3 ), 1, H( 1, 2 ), 1, C1, S1 )
204         CALL DLARTG( H( 1, 2 ), H( 1, 1 ), C2, S2, TEMP )
205*
206         CALL DROT( IHI-ISTARTM+1, B( ISTARTM, IHI ), 1, B( ISTARTM,
207     $              IHI-1 ), 1, C1, S1 )
208         CALL DROT( IHI-ISTARTM+1, B( ISTARTM, IHI-1 ), 1, B( ISTARTM,
209     $              IHI-2 ), 1, C2, S2 )
210         B( IHI-1, IHI-2 ) = ZERO
211         B( IHI, IHI-2 ) = ZERO
212         CALL DROT( IHI-ISTARTM+1, A( ISTARTM, IHI ), 1, A( ISTARTM,
213     $              IHI-1 ), 1, C1, S1 )
214         CALL DROT( IHI-ISTARTM+1, A( ISTARTM, IHI-1 ), 1, A( ISTARTM,
215     $              IHI-2 ), 1, C2, S2 )
216         IF ( ILZ ) THEN
217            CALL DROT( NZ, Z( 1, IHI-ZSTART+1 ), 1, Z( 1, IHI-1-ZSTART+
218     $                 1 ), 1, C1, S1 )
219            CALL DROT( NZ, Z( 1, IHI-1-ZSTART+1 ), 1, Z( 1,
220     $                 IHI-2-ZSTART+1 ), 1, C2, S2 )
221         END IF
222*
223         CALL DLARTG( A( IHI-1, IHI-2 ), A( IHI, IHI-2 ), C1, S1,
224     $                TEMP )
225         A( IHI-1, IHI-2 ) = TEMP
226         A( IHI, IHI-2 ) = ZERO
227         CALL DROT( ISTOPM-IHI+2, A( IHI-1, IHI-1 ), LDA, A( IHI,
228     $              IHI-1 ), LDA, C1, S1 )
229         CALL DROT( ISTOPM-IHI+2, B( IHI-1, IHI-1 ), LDB, B( IHI,
230     $              IHI-1 ), LDB, C1, S1 )
231         IF ( ILQ ) THEN
232            CALL DROT( NQ, Q( 1, IHI-1-QSTART+1 ), 1, Q( 1, IHI-QSTART+
233     $                 1 ), 1, C1, S1 )
234         END IF
235*
236         CALL DLARTG( B( IHI, IHI ), B( IHI, IHI-1 ), C1, S1, TEMP )
237         B( IHI, IHI ) = TEMP
238         B( IHI, IHI-1 ) = ZERO
239         CALL DROT( IHI-ISTARTM, B( ISTARTM, IHI ), 1, B( ISTARTM,
240     $              IHI-1 ), 1, C1, S1 )
241         CALL DROT( IHI-ISTARTM+1, A( ISTARTM, IHI ), 1, A( ISTARTM,
242     $              IHI-1 ), 1, C1, S1 )
243         IF ( ILZ ) THEN
244            CALL DROT( NZ, Z( 1, IHI-ZSTART+1 ), 1, Z( 1, IHI-1-ZSTART+
245     $                 1 ), 1, C1, S1 )
246         END IF
247*
248      ELSE
249*
250*        Normal operation, move bulge down
251*
252         H = B( K+1:K+2, K:K+2 )
253*
254*        Make H upper triangular
255*
256         CALL DLARTG( H( 1, 1 ), H( 2, 1 ), C1, S1, TEMP )
257         H( 2, 1 ) = ZERO
258         H( 1, 1 ) = TEMP
259         CALL DROT( 2, H( 1, 2 ), 2, H( 2, 2 ), 2, C1, S1 )
260*
261*        Calculate Z1 and Z2
262*
263         CALL DLARTG( H( 2, 3 ), H( 2, 2 ), C1, S1, TEMP )
264         CALL DROT( 1, H( 1, 3 ), 1, H( 1, 2 ), 1, C1, S1 )
265         CALL DLARTG( H( 1, 2 ), H( 1, 1 ), C2, S2, TEMP )
266*
267*        Apply transformations from the right
268*
269         CALL DROT( K+3-ISTARTM+1, A( ISTARTM, K+2 ), 1, A( ISTARTM,
270     $              K+1 ), 1, C1, S1 )
271         CALL DROT( K+3-ISTARTM+1, A( ISTARTM, K+1 ), 1, A( ISTARTM,
272     $              K ), 1, C2, S2 )
273         CALL DROT( K+2-ISTARTM+1, B( ISTARTM, K+2 ), 1, B( ISTARTM,
274     $              K+1 ), 1, C1, S1 )
275         CALL DROT( K+2-ISTARTM+1, B( ISTARTM, K+1 ), 1, B( ISTARTM,
276     $              K ), 1, C2, S2 )
277         IF ( ILZ ) THEN
278            CALL DROT( NZ, Z( 1, K+2-ZSTART+1 ), 1, Z( 1, K+1-ZSTART+
279     $                 1 ), 1, C1, S1 )
280            CALL DROT( NZ, Z( 1, K+1-ZSTART+1 ), 1, Z( 1, K-ZSTART+1 ),
281     $                 1, C2, S2 )
282         END IF
283         B( K+1, K ) = ZERO
284         B( K+2, K ) = ZERO
285*
286*        Calculate Q1 and Q2
287*
288         CALL DLARTG( A( K+2, K ), A( K+3, K ), C1, S1, TEMP )
289         A( K+2, K ) = TEMP
290         A( K+3, K ) = ZERO
291         CALL DLARTG( A( K+1, K ), A( K+2, K ), C2, S2, TEMP )
292         A( K+1, K ) = TEMP
293         A( K+2, K ) = ZERO
294*
295*        Apply transformations from the left
296*
297         CALL DROT( ISTOPM-K, A( K+2, K+1 ), LDA, A( K+3, K+1 ), LDA,
298     $              C1, S1 )
299         CALL DROT( ISTOPM-K, A( K+1, K+1 ), LDA, A( K+2, K+1 ), LDA,
300     $              C2, S2 )
301*
302         CALL DROT( ISTOPM-K, B( K+2, K+1 ), LDB, B( K+3, K+1 ), LDB,
303     $              C1, S1 )
304         CALL DROT( ISTOPM-K, B( K+1, K+1 ), LDB, B( K+2, K+1 ), LDB,
305     $              C2, S2 )
306         IF ( ILQ ) THEN
307            CALL DROT( NQ, Q( 1, K+2-QSTART+1 ), 1, Q( 1, K+3-QSTART+
308     $                 1 ), 1, C1, S1 )
309            CALL DROT( NQ, Q( 1, K+1-QSTART+1 ), 1, Q( 1, K+2-QSTART+
310     $                 1 ), 1, C2, S2 )
311         END IF
312*
313      END IF
314*
315*     End of DLAQZ2
316*
317      END SUBROUTINE