LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dlaqz3.f
Go to the documentation of this file.
1*> \brief \b DLAQZ3
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download DLAQZ3 + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlaqz3.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlaqz3.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlaqz3.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE DLAQZ3( ILSCHUR, ILQ, ILZ, N, ILO, IHI, NW, A, LDA, B,
22* $ LDB, Q, LDQ, Z, LDZ, NS, ND, ALPHAR, ALPHAI, BETA, QC, LDQC,
23* $ ZC, LDZC, WORK, LWORK, REC, INFO )
24* IMPLICIT NONE
25*
26* Arguments
27* LOGICAL, INTENT( IN ) :: ILSCHUR, ILQ, ILZ
28* INTEGER, INTENT( IN ) :: N, ILO, IHI, NW, LDA, LDB, LDQ, LDZ,
29* $ LDQC, LDZC, LWORK, REC
30*
31* DOUBLE PRECISION, INTENT( INOUT ) :: A( LDA, * ), B( LDB, * ),
32* $ Q( LDQ, * ), Z( LDZ, * ), ALPHAR( * ), ALPHAI( * ), BETA( * )
33* INTEGER, INTENT( OUT ) :: NS, ND, INFO
34* DOUBLE PRECISION :: QC( LDQC, * ), ZC( LDZC, * ), WORK( * )
35* ..
36*
37*
38*> \par Purpose:
39* =============
40*>
41*> \verbatim
42*>
43*> DLAQZ3 performs AED
44*> \endverbatim
45*
46* Arguments:
47* ==========
48*
49*> \param[in] ILSCHUR
50*> \verbatim
51*> ILSCHUR is LOGICAL
52*> Determines whether or not to update the full Schur form
53*> \endverbatim
54*>
55*> \param[in] ILQ
56*> \verbatim
57*> ILQ is LOGICAL
58*> Determines whether or not to update the matrix Q
59*> \endverbatim
60*>
61*> \param[in] ILZ
62*> \verbatim
63*> ILZ is LOGICAL
64*> Determines whether or not to update the matrix Z
65*> \endverbatim
66*>
67*> \param[in] N
68*> \verbatim
69*> N is INTEGER
70*> The order of the matrices A, B, Q, and Z. N >= 0.
71*> \endverbatim
72*>
73*> \param[in] ILO
74*> \verbatim
75*> ILO is INTEGER
76*> \endverbatim
77*>
78*> \param[in] IHI
79*> \verbatim
80*> IHI is INTEGER
81*> ILO and IHI mark the rows and columns of (A,B) which
82*> are to be normalized
83*> \endverbatim
84*>
85*> \param[in] NW
86*> \verbatim
87*> NW is INTEGER
88*> The desired size of the deflation window.
89*> \endverbatim
90*>
91*> \param[in,out] A
92*> \verbatim
93*> A is DOUBLE PRECISION array, dimension (LDA, N)
94*> \endverbatim
95*>
96*> \param[in] LDA
97*> \verbatim
98*> LDA is INTEGER
99*> The leading dimension of the array A. LDA >= max( 1, N ).
100*> \endverbatim
101*>
102*> \param[in,out] B
103*> \verbatim
104*> B is DOUBLE PRECISION array, dimension (LDB, N)
105*> \endverbatim
106*>
107*> \param[in] LDB
108*> \verbatim
109*> LDB is INTEGER
110*> The leading dimension of the array B. LDB >= max( 1, N ).
111*> \endverbatim
112*>
113*> \param[in,out] Q
114*> \verbatim
115*> Q is DOUBLE PRECISION array, dimension (LDQ, N)
116*> \endverbatim
117*>
118*> \param[in] LDQ
119*> \verbatim
120*> LDQ is INTEGER
121*> \endverbatim
122*>
123*> \param[in,out] Z
124*> \verbatim
125*> Z is DOUBLE PRECISION array, dimension (LDZ, N)
126*> \endverbatim
127*>
128*> \param[in] LDZ
129*> \verbatim
130*> LDZ is INTEGER
131*> \endverbatim
132*>
133*> \param[out] NS
134*> \verbatim
135*> NS is INTEGER
136*> The number of unconverged eigenvalues available to
137*> use as shifts.
138*> \endverbatim
139*>
140*> \param[out] ND
141*> \verbatim
142*> ND is INTEGER
143*> The number of converged eigenvalues found.
144*> \endverbatim
145*>
146*> \param[out] ALPHAR
147*> \verbatim
148*> ALPHAR is DOUBLE PRECISION array, dimension (N)
149*> The real parts of each scalar alpha defining an eigenvalue
150*> of GNEP.
151*> \endverbatim
152*>
153*> \param[out] ALPHAI
154*> \verbatim
155*> ALPHAI is DOUBLE PRECISION array, dimension (N)
156*> The imaginary parts of each scalar alpha defining an
157*> eigenvalue of GNEP.
158*> If ALPHAI(j) is zero, then the j-th eigenvalue is real; if
159*> positive, then the j-th and (j+1)-st eigenvalues are a
160*> complex conjugate pair, with ALPHAI(j+1) = -ALPHAI(j).
161*> \endverbatim
162*>
163*> \param[out] BETA
164*> \verbatim
165*> BETA is DOUBLE PRECISION array, dimension (N)
166*> The scalars beta that define the eigenvalues of GNEP.
167*> Together, the quantities alpha = (ALPHAR(j),ALPHAI(j)) and
168*> beta = BETA(j) represent the j-th eigenvalue of the matrix
169*> pair (A,B), in one of the forms lambda = alpha/beta or
170*> mu = beta/alpha. Since either lambda or mu may overflow,
171*> they should not, in general, be computed.
172*> \endverbatim
173*>
174*> \param[in,out] QC
175*> \verbatim
176*> QC is DOUBLE PRECISION array, dimension (LDQC, NW)
177*> \endverbatim
178*>
179*> \param[in] LDQC
180*> \verbatim
181*> LDQC is INTEGER
182*> \endverbatim
183*>
184*> \param[in,out] ZC
185*> \verbatim
186*> ZC is DOUBLE PRECISION array, dimension (LDZC, NW)
187*> \endverbatim
188*>
189*> \param[in] LDZC
190*> \verbatim
191*> LDZ is INTEGER
192*> \endverbatim
193*>
194*> \param[out] WORK
195*> \verbatim
196*> WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
197*> On exit, if INFO >= 0, WORK(1) returns the optimal LWORK.
198*> \endverbatim
199*>
200*> \param[in] LWORK
201*> \verbatim
202*> LWORK is INTEGER
203*> The dimension of the array WORK. LWORK >= max(1,N).
204*>
205*> If LWORK = -1, then a workspace query is assumed; the routine
206*> only calculates the optimal size of the WORK array, returns
207*> this value as the first entry of the WORK array, and no error
208*> message related to LWORK is issued by XERBLA.
209*> \endverbatim
210*>
211*> \param[in] REC
212*> \verbatim
213*> REC is INTEGER
214*> REC indicates the current recursion level. Should be set
215*> to 0 on first call.
216*> \endverbatim
217*>
218*> \param[out] INFO
219*> \verbatim
220*> INFO is INTEGER
221*> = 0: successful exit
222*> < 0: if INFO = -i, the i-th argument had an illegal value
223*> \endverbatim
224*
225* Authors:
226* ========
227*
228*> \author Thijs Steel, KU Leuven
229*
230*> \date May 2020
231*
232*> \ingroup laqz3
233*>
234* =====================================================================
235 RECURSIVE SUBROUTINE dlaqz3( ILSCHUR, ILQ, ILZ, N, ILO, IHI, NW,
236 $ A, LDA, B, LDB, Q, LDQ, Z, LDZ, NS,
237 $ ND, ALPHAR, ALPHAI, BETA, QC, LDQC,
238 $ ZC, LDZC, WORK, LWORK, REC, INFO )
239 IMPLICIT NONE
240
241* Arguments
242 LOGICAL, INTENT( IN ) :: ilschur, ilq, ilz
243 INTEGER, INTENT( IN ) :: n, ilo, ihi, nw, lda, ldb, ldq, ldz,
244 $ ldqc, ldzc, lwork, rec
245
246 DOUBLE PRECISION, INTENT( INOUT ) :: a( lda, * ), b( ldb, * ),
247 $ q( ldq, * ), z( ldz, * ), alphar( * ),
248 $ alphai( * ), beta( * )
249 INTEGER, INTENT( OUT ) :: ns, nd, info
250 DOUBLE PRECISION :: qc( ldqc, * ), zc( ldzc, * ), work( * )
251
252* Parameters
253 DOUBLE PRECISION :: zero, one, half
254 PARAMETER( zero = 0.0d0, one = 1.0d0, half = 0.5d0 )
255
256* Local Scalars
257 LOGICAL :: bulge
258 INTEGER :: jw, kwtop, kwbot, istopm, istartm, k, k2, dtgexc_info,
259 $ ifst, ilst, lworkreq, qz_small_info
260 DOUBLE PRECISION :: s, smlnum, ulp, safmin, safmax, c1, s1, temp
261
262* External Functions
263 EXTERNAL :: xerbla, dtgexc, dlaqz0, dlacpy, dlaset,
265 DOUBLE PRECISION, EXTERNAL :: dlamch
266
267 info = 0
268
269* Set up deflation window
270 jw = min( nw, ihi-ilo+1 )
271 kwtop = ihi-jw+1
272 IF ( kwtop .EQ. ilo ) THEN
273 s = zero
274 ELSE
275 s = a( kwtop, kwtop-1 )
276 END IF
277
278* Determine required workspace
279 ifst = 1
280 ilst = jw
281 CALL dtgexc( .true., .true., jw, a, lda, b, ldb, qc, ldqc, zc,
282 $ ldzc, ifst, ilst, work, -1, dtgexc_info )
283 lworkreq = int( work( 1 ) )
284 CALL dlaqz0( 'S', 'V', 'V', jw, 1, jw, a( kwtop, kwtop ), lda,
285 $ b( kwtop, kwtop ), ldb, alphar, alphai, beta, qc,
286 $ ldqc, zc, ldzc, work, -1, rec+1, qz_small_info )
287 lworkreq = max( lworkreq, int( work( 1 ) )+2*jw**2 )
288 lworkreq = max( lworkreq, n*nw, 2*nw**2+n )
289 IF ( lwork .EQ.-1 ) THEN
290* workspace query, quick return
291 work( 1 ) = lworkreq
292 RETURN
293 ELSE IF ( lwork .LT. lworkreq ) THEN
294 info = -26
295 END IF
296
297 IF( info.NE.0 ) THEN
298 CALL xerbla( 'DLAQZ3', -info )
299 RETURN
300 END IF
301
302* Get machine constants
303 safmin = dlamch( 'SAFE MINIMUM' )
304 safmax = one/safmin
305 ulp = dlamch( 'PRECISION' )
306 smlnum = safmin*( dble( n )/ulp )
307
308 IF ( ihi .EQ. kwtop ) THEN
309* 1 by 1 deflation window, just try a regular deflation
310 alphar( kwtop ) = a( kwtop, kwtop )
311 alphai( kwtop ) = zero
312 beta( kwtop ) = b( kwtop, kwtop )
313 ns = 1
314 nd = 0
315 IF ( abs( s ) .LE. max( smlnum, ulp*abs( a( kwtop,
316 $ kwtop ) ) ) ) THEN
317 ns = 0
318 nd = 1
319 IF ( kwtop .GT. ilo ) THEN
320 a( kwtop, kwtop-1 ) = zero
321 END IF
322 END IF
323 END IF
324
325
326* Store window in case of convergence failure
327 CALL dlacpy( 'ALL', jw, jw, a( kwtop, kwtop ), lda, work, jw )
328 CALL dlacpy( 'ALL', jw, jw, b( kwtop, kwtop ), ldb, work( jw**2+
329 $ 1 ), jw )
330
331* Transform window to real schur form
332 CALL dlaset( 'FULL', jw, jw, zero, one, qc, ldqc )
333 CALL dlaset( 'FULL', jw, jw, zero, one, zc, ldzc )
334 CALL dlaqz0( 'S', 'V', 'V', jw, 1, jw, a( kwtop, kwtop ), lda,
335 $ b( kwtop, kwtop ), ldb, alphar, alphai, beta, qc,
336 $ ldqc, zc, ldzc, work( 2*jw**2+1 ), lwork-2*jw**2,
337 $ rec+1, qz_small_info )
338
339 IF( qz_small_info .NE. 0 ) THEN
340* Convergence failure, restore the window and exit
341 nd = 0
342 ns = jw-qz_small_info
343 CALL dlacpy( 'ALL', jw, jw, work, jw, a( kwtop, kwtop ), lda )
344 CALL dlacpy( 'ALL', jw, jw, work( jw**2+1 ), jw, b( kwtop,
345 $ kwtop ), ldb )
346 RETURN
347 END IF
348
349* Deflation detection loop
350 IF ( kwtop .EQ. ilo .OR. s .EQ. zero ) THEN
351 kwbot = kwtop-1
352 ELSE
353 kwbot = ihi
354 k = 1
355 k2 = 1
356 DO WHILE ( k .LE. jw )
357 bulge = .false.
358 IF ( kwbot-kwtop+1 .GE. 2 ) THEN
359 bulge = a( kwbot, kwbot-1 ) .NE. zero
360 END IF
361 IF ( bulge ) THEN
362
363* Try to deflate complex conjugate eigenvalue pair
364 temp = abs( a( kwbot, kwbot ) )+sqrt( abs( a( kwbot,
365 $ kwbot-1 ) ) )*sqrt( abs( a( kwbot-1, kwbot ) ) )
366 IF( temp .EQ. zero )THEN
367 temp = abs( s )
368 END IF
369 IF ( max( abs( s*qc( 1, kwbot-kwtop ) ), abs( s*qc( 1,
370 $ kwbot-kwtop+1 ) ) ) .LE. max( smlnum,
371 $ ulp*temp ) ) THEN
372* Deflatable
373 kwbot = kwbot-2
374 ELSE
375* Not deflatable, move out of the way
376 ifst = kwbot-kwtop+1
377 ilst = k2
378 CALL dtgexc( .true., .true., jw, a( kwtop, kwtop ),
379 $ lda, b( kwtop, kwtop ), ldb, qc, ldqc,
380 $ zc, ldzc, ifst, ilst, work, lwork,
381 $ dtgexc_info )
382 k2 = k2+2
383 END IF
384 k = k+2
385 ELSE
386
387* Try to deflate real eigenvalue
388 temp = abs( a( kwbot, kwbot ) )
389 IF( temp .EQ. zero ) THEN
390 temp = abs( s )
391 END IF
392 IF ( ( abs( s*qc( 1, kwbot-kwtop+1 ) ) ) .LE. max( ulp*
393 $ temp, smlnum ) ) THEN
394* Deflatable
395 kwbot = kwbot-1
396 ELSE
397* Not deflatable, move out of the way
398 ifst = kwbot-kwtop+1
399 ilst = k2
400 CALL dtgexc( .true., .true., jw, a( kwtop, kwtop ),
401 $ lda, b( kwtop, kwtop ), ldb, qc, ldqc,
402 $ zc, ldzc, ifst, ilst, work, lwork,
403 $ dtgexc_info )
404 k2 = k2+1
405 END IF
406
407 k = k+1
408
409 END IF
410 END DO
411 END IF
412
413* Store eigenvalues
414 nd = ihi-kwbot
415 ns = jw-nd
416 k = kwtop
417 DO WHILE ( k .LE. ihi )
418 bulge = .false.
419 IF ( k .LT. ihi ) THEN
420 IF ( a( k+1, k ) .NE. zero ) THEN
421 bulge = .true.
422 END IF
423 END IF
424 IF ( bulge ) THEN
425* 2x2 eigenvalue block
426 CALL dlag2( a( k, k ), lda, b( k, k ), ldb, safmin,
427 $ beta( k ), beta( k+1 ), alphar( k ),
428 $ alphar( k+1 ), alphai( k ) )
429 alphai( k+1 ) = -alphai( k )
430 k = k+2
431 ELSE
432* 1x1 eigenvalue block
433 alphar( k ) = a( k, k )
434 alphai( k ) = zero
435 beta( k ) = b( k, k )
436 k = k+1
437 END IF
438 END DO
439
440 IF ( kwtop .NE. ilo .AND. s .NE. zero ) THEN
441* Reflect spike back, this will create optimally packed bulges
442 a( kwtop:kwbot, kwtop-1 ) = a( kwtop, kwtop-1 )*qc( 1,
443 $ 1:jw-nd )
444 DO k = kwbot-1, kwtop, -1
445 CALL dlartg( a( k, kwtop-1 ), a( k+1, kwtop-1 ), c1, s1,
446 $ temp )
447 a( k, kwtop-1 ) = temp
448 a( k+1, kwtop-1 ) = zero
449 k2 = max( kwtop, k-1 )
450 CALL drot( ihi-k2+1, a( k, k2 ), lda, a( k+1, k2 ), lda, c1,
451 $ s1 )
452 CALL drot( ihi-( k-1 )+1, b( k, k-1 ), ldb, b( k+1, k-1 ),
453 $ ldb, c1, s1 )
454 CALL drot( jw, qc( 1, k-kwtop+1 ), 1, qc( 1, k+1-kwtop+1 ),
455 $ 1, c1, s1 )
456 END DO
457
458* Chase bulges down
459 istartm = kwtop
460 istopm = ihi
461 k = kwbot-1
462 DO WHILE ( k .GE. kwtop )
463 IF ( ( k .GE. kwtop+1 ) .AND. a( k+1, k-1 ) .NE. zero ) THEN
464
465* Move double pole block down and remove it
466 DO k2 = k-1, kwbot-2
467 CALL dlaqz2( .true., .true., k2, kwtop, kwtop+jw-1,
468 $ kwbot, a, lda, b, ldb, jw, kwtop, qc,
469 $ ldqc, jw, kwtop, zc, ldzc )
470 END DO
471
472 k = k-2
473 ELSE
474
475* k points to single shift
476 DO k2 = k, kwbot-2
477
478* Move shift down
479 CALL dlartg( b( k2+1, k2+1 ), b( k2+1, k2 ), c1, s1,
480 $ temp )
481 b( k2+1, k2+1 ) = temp
482 b( k2+1, k2 ) = zero
483 CALL drot( k2+2-istartm+1, a( istartm, k2+1 ), 1,
484 $ a( istartm, k2 ), 1, c1, s1 )
485 CALL drot( k2-istartm+1, b( istartm, k2+1 ), 1,
486 $ b( istartm, k2 ), 1, c1, s1 )
487 CALL drot( jw, zc( 1, k2+1-kwtop+1 ), 1, zc( 1,
488 $ k2-kwtop+1 ), 1, c1, s1 )
489
490 CALL dlartg( a( k2+1, k2 ), a( k2+2, k2 ), c1, s1,
491 $ temp )
492 a( k2+1, k2 ) = temp
493 a( k2+2, k2 ) = zero
494 CALL drot( istopm-k2, a( k2+1, k2+1 ), lda, a( k2+2,
495 $ k2+1 ), lda, c1, s1 )
496 CALL drot( istopm-k2, b( k2+1, k2+1 ), ldb, b( k2+2,
497 $ k2+1 ), ldb, c1, s1 )
498 CALL drot( jw, qc( 1, k2+1-kwtop+1 ), 1, qc( 1,
499 $ k2+2-kwtop+1 ), 1, c1, s1 )
500
501 END DO
502
503* Remove the shift
504 CALL dlartg( b( kwbot, kwbot ), b( kwbot, kwbot-1 ), c1,
505 $ s1, temp )
506 b( kwbot, kwbot ) = temp
507 b( kwbot, kwbot-1 ) = zero
508 CALL drot( kwbot-istartm, b( istartm, kwbot ), 1,
509 $ b( istartm, kwbot-1 ), 1, c1, s1 )
510 CALL drot( kwbot-istartm+1, a( istartm, kwbot ), 1,
511 $ a( istartm, kwbot-1 ), 1, c1, s1 )
512 CALL drot( jw, zc( 1, kwbot-kwtop+1 ), 1, zc( 1,
513 $ kwbot-1-kwtop+1 ), 1, c1, s1 )
514
515 k = k-1
516 END IF
517 END DO
518
519 END IF
520
521* Apply Qc and Zc to rest of the matrix
522 IF ( ilschur ) THEN
523 istartm = 1
524 istopm = n
525 ELSE
526 istartm = ilo
527 istopm = ihi
528 END IF
529
530 IF ( istopm-ihi > 0 ) THEN
531 CALL dgemm( 'T', 'N', jw, istopm-ihi, jw, one, qc, ldqc,
532 $ a( kwtop, ihi+1 ), lda, zero, work, jw )
533 CALL dlacpy( 'ALL', jw, istopm-ihi, work, jw, a( kwtop,
534 $ ihi+1 ), lda )
535 CALL dgemm( 'T', 'N', jw, istopm-ihi, jw, one, qc, ldqc,
536 $ b( kwtop, ihi+1 ), ldb, zero, work, jw )
537 CALL dlacpy( 'ALL', jw, istopm-ihi, work, jw, b( kwtop,
538 $ ihi+1 ), ldb )
539 END IF
540 IF ( ilq ) THEN
541 CALL dgemm( 'N', 'N', n, jw, jw, one, q( 1, kwtop ), ldq, qc,
542 $ ldqc, zero, work, n )
543 CALL dlacpy( 'ALL', n, jw, work, n, q( 1, kwtop ), ldq )
544 END IF
545
546 IF ( kwtop-1-istartm+1 > 0 ) THEN
547 CALL dgemm( 'N', 'N', kwtop-istartm, jw, jw, one, a( istartm,
548 $ kwtop ), lda, zc, ldzc, zero, work,
549 $ kwtop-istartm )
550 CALL dlacpy( 'ALL', kwtop-istartm, jw, work, kwtop-istartm,
551 $ a( istartm, kwtop ), lda )
552 CALL dgemm( 'N', 'N', kwtop-istartm, jw, jw, one, b( istartm,
553 $ kwtop ), ldb, zc, ldzc, zero, work,
554 $ kwtop-istartm )
555 CALL dlacpy( 'ALL', kwtop-istartm, jw, work, kwtop-istartm,
556 $ b( istartm, kwtop ), ldb )
557 END IF
558 IF ( ilz ) THEN
559 CALL dgemm( 'N', 'N', n, jw, jw, one, z( 1, kwtop ), ldz, zc,
560 $ ldzc, zero, work, n )
561 CALL dlacpy( 'ALL', n, jw, work, n, z( 1, kwtop ), ldz )
562 END IF
563
564 END SUBROUTINE
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
DGEMM
Definition dgemm.f:188
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 dlag2(a, lda, b, ldb, safmin, scale1, scale2, wr1, wr2, wi)
DLAG2 computes the eigenvalues of a 2-by-2 generalized eigenvalue problem, with scaling as necessary ...
Definition dlag2.f:156
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
recursive subroutine dlaqz0(wants, wantq, wantz, n, ilo, ihi, a, lda, b, ldb, alphar, alphai, beta, q, ldq, z, ldz, work, lwork, rec, info)
DLAQZ0
Definition dlaqz0.f:306
subroutine dlaqz2(ilq, ilz, k, istartm, istopm, ihi, a, lda, b, ldb, nq, qstart, q, ldq, nz, zstart, z, ldz)
DLAQZ2
Definition dlaqz2.f:174
recursive subroutine dlaqz3(ilschur, ilq, ilz, n, ilo, ihi, nw, a, lda, b, ldb, q, ldq, z, ldz, ns, nd, alphar, alphai, beta, qc, ldqc, zc, ldzc, work, lwork, rec, info)
DLAQZ3
Definition dlaqz3.f:239
subroutine dlartg(f, g, c, s, r)
DLARTG generates a plane rotation with real cosine and real sine.
Definition dlartg.f90:111
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 drot(n, dx, incx, dy, incy, c, s)
DROT
Definition drot.f:92
subroutine dtgexc(wantq, wantz, n, a, lda, b, ldb, q, ldq, z, ldz, ifst, ilst, work, lwork, info)
DTGEXC
Definition dtgexc.f:220