LAPACK 3.11.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dggev3.f
Go to the documentation of this file.
1*> \brief <b> DGGEV3 computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices (blocked algorithm)</b>
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download DGGEV3 + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dggev3.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dggev3.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dggev3.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE DGGEV3( JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR,
22* $ ALPHAI, BETA, VL, LDVL, VR, LDVR, WORK, LWORK,
23* $ INFO )
24*
25* .. Scalar Arguments ..
26* CHARACTER JOBVL, JOBVR
27* INTEGER INFO, LDA, LDB, LDVL, LDVR, LWORK, N
28* ..
29* .. Array Arguments ..
30* DOUBLE PRECISION A( LDA, * ), ALPHAI( * ), ALPHAR( * ),
31* $ B( LDB, * ), BETA( * ), VL( LDVL, * ),
32* $ VR( LDVR, * ), WORK( * )
33* ..
34*
35*
36*> \par Purpose:
37* =============
38*>
39*> \verbatim
40*>
41*> DGGEV3 computes for a pair of N-by-N real nonsymmetric matrices (A,B)
42*> the generalized eigenvalues, and optionally, the left and/or right
43*> generalized eigenvectors.
44*>
45*> A generalized eigenvalue for a pair of matrices (A,B) is a scalar
46*> lambda or a ratio alpha/beta = lambda, such that A - lambda*B is
47*> singular. It is usually represented as the pair (alpha,beta), as
48*> there is a reasonable interpretation for beta=0, and even for both
49*> being zero.
50*>
51*> The right eigenvector v(j) corresponding to the eigenvalue lambda(j)
52*> of (A,B) satisfies
53*>
54*> A * v(j) = lambda(j) * B * v(j).
55*>
56*> The left eigenvector u(j) corresponding to the eigenvalue lambda(j)
57*> of (A,B) satisfies
58*>
59*> u(j)**H * A = lambda(j) * u(j)**H * B .
60*>
61*> where u(j)**H is the conjugate-transpose of u(j).
62*>
63*> \endverbatim
64*
65* Arguments:
66* ==========
67*
68*> \param[in] JOBVL
69*> \verbatim
70*> JOBVL is CHARACTER*1
71*> = 'N': do not compute the left generalized eigenvectors;
72*> = 'V': compute the left generalized eigenvectors.
73*> \endverbatim
74*>
75*> \param[in] JOBVR
76*> \verbatim
77*> JOBVR is CHARACTER*1
78*> = 'N': do not compute the right generalized eigenvectors;
79*> = 'V': compute the right generalized eigenvectors.
80*> \endverbatim
81*>
82*> \param[in] N
83*> \verbatim
84*> N is INTEGER
85*> The order of the matrices A, B, VL, and VR. N >= 0.
86*> \endverbatim
87*>
88*> \param[in,out] A
89*> \verbatim
90*> A is DOUBLE PRECISION array, dimension (LDA, N)
91*> On entry, the matrix A in the pair (A,B).
92*> On exit, A has been overwritten.
93*> \endverbatim
94*>
95*> \param[in] LDA
96*> \verbatim
97*> LDA is INTEGER
98*> The leading dimension of A. LDA >= max(1,N).
99*> \endverbatim
100*>
101*> \param[in,out] B
102*> \verbatim
103*> B is DOUBLE PRECISION array, dimension (LDB, N)
104*> On entry, the matrix B in the pair (A,B).
105*> On exit, B has been overwritten.
106*> \endverbatim
107*>
108*> \param[in] LDB
109*> \verbatim
110*> LDB is INTEGER
111*> The leading dimension of B. LDB >= max(1,N).
112*> \endverbatim
113*>
114*> \param[out] ALPHAR
115*> \verbatim
116*> ALPHAR is DOUBLE PRECISION array, dimension (N)
117*> \endverbatim
118*>
119*> \param[out] ALPHAI
120*> \verbatim
121*> ALPHAI is DOUBLE PRECISION array, dimension (N)
122*> \endverbatim
123*>
124*> \param[out] BETA
125*> \verbatim
126*> BETA is DOUBLE PRECISION array, dimension (N)
127*> On exit, (ALPHAR(j) + ALPHAI(j)*i)/BETA(j), j=1,...,N, will
128*> be the generalized eigenvalues. If ALPHAI(j) is zero, then
129*> the j-th eigenvalue is real; if positive, then the j-th and
130*> (j+1)-st eigenvalues are a complex conjugate pair, with
131*> ALPHAI(j+1) negative.
132*>
133*> Note: the quotients ALPHAR(j)/BETA(j) and ALPHAI(j)/BETA(j)
134*> may easily over- or underflow, and BETA(j) may even be zero.
135*> Thus, the user should avoid naively computing the ratio
136*> alpha/beta. However, ALPHAR and ALPHAI will be always less
137*> than and usually comparable with norm(A) in magnitude, and
138*> BETA always less than and usually comparable with norm(B).
139*> \endverbatim
140*>
141*> \param[out] VL
142*> \verbatim
143*> VL is DOUBLE PRECISION array, dimension (LDVL,N)
144*> If JOBVL = 'V', the left eigenvectors u(j) are stored one
145*> after another in the columns of VL, in the same order as
146*> their eigenvalues. If the j-th eigenvalue is real, then
147*> u(j) = VL(:,j), the j-th column of VL. If the j-th and
148*> (j+1)-th eigenvalues form a complex conjugate pair, then
149*> u(j) = VL(:,j)+i*VL(:,j+1) and u(j+1) = VL(:,j)-i*VL(:,j+1).
150*> Each eigenvector is scaled so the largest component has
151*> abs(real part)+abs(imag. part)=1.
152*> Not referenced if JOBVL = 'N'.
153*> \endverbatim
154*>
155*> \param[in] LDVL
156*> \verbatim
157*> LDVL is INTEGER
158*> The leading dimension of the matrix VL. LDVL >= 1, and
159*> if JOBVL = 'V', LDVL >= N.
160*> \endverbatim
161*>
162*> \param[out] VR
163*> \verbatim
164*> VR is DOUBLE PRECISION array, dimension (LDVR,N)
165*> If JOBVR = 'V', the right eigenvectors v(j) are stored one
166*> after another in the columns of VR, in the same order as
167*> their eigenvalues. If the j-th eigenvalue is real, then
168*> v(j) = VR(:,j), the j-th column of VR. If the j-th and
169*> (j+1)-th eigenvalues form a complex conjugate pair, then
170*> v(j) = VR(:,j)+i*VR(:,j+1) and v(j+1) = VR(:,j)-i*VR(:,j+1).
171*> Each eigenvector is scaled so the largest component has
172*> abs(real part)+abs(imag. part)=1.
173*> Not referenced if JOBVR = 'N'.
174*> \endverbatim
175*>
176*> \param[in] LDVR
177*> \verbatim
178*> LDVR is INTEGER
179*> The leading dimension of the matrix VR. LDVR >= 1, and
180*> if JOBVR = 'V', LDVR >= N.
181*> \endverbatim
182*>
183*> \param[out] WORK
184*> \verbatim
185*> WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
186*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
187*> \endverbatim
188*>
189*> \param[in] LWORK
190*> \verbatim
191*> LWORK is INTEGER
192*>
193*> If LWORK = -1, then a workspace query is assumed; the routine
194*> only calculates the optimal size of the WORK array, returns
195*> this value as the first entry of the WORK array, and no error
196*> message related to LWORK is issued by XERBLA.
197*> \endverbatim
198*>
199*> \param[out] INFO
200*> \verbatim
201*> INFO is INTEGER
202*> = 0: successful exit
203*> < 0: if INFO = -i, the i-th argument had an illegal value.
204*> = 1,...,N:
205*> The QZ iteration failed. No eigenvectors have been
206*> calculated, but ALPHAR(j), ALPHAI(j), and BETA(j)
207*> should be correct for j=INFO+1,...,N.
208*> > N: =N+1: other than QZ iteration failed in DLAQZ0.
209*> =N+2: error return from DTGEVC.
210*> \endverbatim
211*
212* Authors:
213* ========
214*
215*> \author Univ. of Tennessee
216*> \author Univ. of California Berkeley
217*> \author Univ. of Colorado Denver
218*> \author NAG Ltd.
219*
220*> \ingroup doubleGEeigen
221*
222* =====================================================================
223 SUBROUTINE dggev3( JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR,
224 $ ALPHAI, BETA, VL, LDVL, VR, LDVR, WORK, LWORK,
225 $ INFO )
226*
227* -- LAPACK driver routine --
228* -- LAPACK is a software package provided by Univ. of Tennessee, --
229* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
230*
231* .. Scalar Arguments ..
232 CHARACTER JOBVL, JOBVR
233 INTEGER INFO, LDA, LDB, LDVL, LDVR, LWORK, N
234* ..
235* .. Array Arguments ..
236 DOUBLE PRECISION A( LDA, * ), ALPHAI( * ), ALPHAR( * ),
237 $ B( LDB, * ), BETA( * ), VL( LDVL, * ),
238 $ vr( ldvr, * ), work( * )
239* ..
240*
241* =====================================================================
242*
243* .. Parameters ..
244 DOUBLE PRECISION ZERO, ONE
245 PARAMETER ( ZERO = 0.0d+0, one = 1.0d+0 )
246* ..
247* .. Local Scalars ..
248 LOGICAL ILASCL, ILBSCL, ILV, ILVL, ILVR, LQUERY
249 CHARACTER CHTEMP
250 INTEGER ICOLS, IERR, IHI, IJOBVL, IJOBVR, ILEFT, ILO,
251 $ in, iright, irows, itau, iwrk, jc, jr, lwkopt
252 DOUBLE PRECISION ANRM, ANRMTO, BIGNUM, BNRM, BNRMTO, EPS,
253 $ SMLNUM, TEMP
254* ..
255* .. Local Arrays ..
256 LOGICAL LDUMMA( 1 )
257* ..
258* .. External Subroutines ..
259 EXTERNAL dgeqrf, dggbak, dggbal, dgghd3, dlaqz0, dlabad,
261 $ xerbla
262* ..
263* .. External Functions ..
264 LOGICAL LSAME
265 DOUBLE PRECISION DLAMCH, DLANGE
266 EXTERNAL lsame, dlamch, dlange
267* ..
268* .. Intrinsic Functions ..
269 INTRINSIC abs, max, sqrt
270* ..
271* .. Executable Statements ..
272*
273* Decode the input arguments
274*
275 IF( lsame( jobvl, 'N' ) ) THEN
276 ijobvl = 1
277 ilvl = .false.
278 ELSE IF( lsame( jobvl, 'V' ) ) THEN
279 ijobvl = 2
280 ilvl = .true.
281 ELSE
282 ijobvl = -1
283 ilvl = .false.
284 END IF
285*
286 IF( lsame( jobvr, 'N' ) ) THEN
287 ijobvr = 1
288 ilvr = .false.
289 ELSE IF( lsame( jobvr, 'V' ) ) THEN
290 ijobvr = 2
291 ilvr = .true.
292 ELSE
293 ijobvr = -1
294 ilvr = .false.
295 END IF
296 ilv = ilvl .OR. ilvr
297*
298* Test the input arguments
299*
300 info = 0
301 lquery = ( lwork.EQ.-1 )
302 IF( ijobvl.LE.0 ) THEN
303 info = -1
304 ELSE IF( ijobvr.LE.0 ) THEN
305 info = -2
306 ELSE IF( n.LT.0 ) THEN
307 info = -3
308 ELSE IF( lda.LT.max( 1, n ) ) THEN
309 info = -5
310 ELSE IF( ldb.LT.max( 1, n ) ) THEN
311 info = -7
312 ELSE IF( ldvl.LT.1 .OR. ( ilvl .AND. ldvl.LT.n ) ) THEN
313 info = -12
314 ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) ) THEN
315 info = -14
316 ELSE IF( lwork.LT.max( 1, 8*n ) .AND. .NOT.lquery ) THEN
317 info = -16
318 END IF
319*
320* Compute workspace
321*
322 IF( info.EQ.0 ) THEN
323 CALL dgeqrf( n, n, b, ldb, work, work, -1, ierr )
324 lwkopt = max(1, 8*n, 3*n+int( work( 1 ) ) )
325 CALL dormqr( 'L', 'T', n, n, n, b, ldb, work, a, lda, work, -1,
326 $ ierr )
327 lwkopt = max( lwkopt, 3*n+int( work( 1 ) ) )
328 IF( ilvl ) THEN
329 CALL dorgqr( n, n, n, vl, ldvl, work, work, -1, ierr )
330 lwkopt = max( lwkopt, 3*n+int( work( 1 ) ) )
331 END IF
332 IF( ilv ) THEN
333 CALL dgghd3( jobvl, jobvr, n, 1, n, a, lda, b, ldb, vl,
334 $ ldvl, vr, ldvr, work, -1, ierr )
335 lwkopt = max( lwkopt, 3*n+int( work( 1 ) ) )
336 CALL dlaqz0( 'S', jobvl, jobvr, n, 1, n, a, lda, b, ldb,
337 $ alphar, alphai, beta, vl, ldvl, vr, ldvr,
338 $ work, -1, 0, ierr )
339 lwkopt = max( lwkopt, 2*n+int( work( 1 ) ) )
340 ELSE
341 CALL dgghd3( 'N', 'N', n, 1, n, a, lda, b, ldb, vl, ldvl,
342 $ vr, ldvr, work, -1, ierr )
343 lwkopt = max( lwkopt, 3*n+int( work( 1 ) ) )
344 CALL dlaqz0( 'E', jobvl, jobvr, n, 1, n, a, lda, b, ldb,
345 $ alphar, alphai, beta, vl, ldvl, vr, ldvr,
346 $ work, -1, 0, ierr )
347 lwkopt = max( lwkopt, 2*n+int( work( 1 ) ) )
348 END IF
349
350 work( 1 ) = lwkopt
351 END IF
352*
353 IF( info.NE.0 ) THEN
354 CALL xerbla( 'DGGEV3 ', -info )
355 RETURN
356 ELSE IF( lquery ) THEN
357 RETURN
358 END IF
359*
360* Quick return if possible
361*
362 IF( n.EQ.0 )
363 $ RETURN
364*
365* Get machine constants
366*
367 eps = dlamch( 'P' )
368 smlnum = dlamch( 'S' )
369 bignum = one / smlnum
370 CALL dlabad( smlnum, bignum )
371 smlnum = sqrt( smlnum ) / eps
372 bignum = one / smlnum
373*
374* Scale A if max element outside range [SMLNUM,BIGNUM]
375*
376 anrm = dlange( 'M', n, n, a, lda, work )
377 ilascl = .false.
378 IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
379 anrmto = smlnum
380 ilascl = .true.
381 ELSE IF( anrm.GT.bignum ) THEN
382 anrmto = bignum
383 ilascl = .true.
384 END IF
385 IF( ilascl )
386 $ CALL dlascl( 'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
387*
388* Scale B if max element outside range [SMLNUM,BIGNUM]
389*
390 bnrm = dlange( 'M', n, n, b, ldb, work )
391 ilbscl = .false.
392 IF( bnrm.GT.zero .AND. bnrm.LT.smlnum ) THEN
393 bnrmto = smlnum
394 ilbscl = .true.
395 ELSE IF( bnrm.GT.bignum ) THEN
396 bnrmto = bignum
397 ilbscl = .true.
398 END IF
399 IF( ilbscl )
400 $ CALL dlascl( 'G', 0, 0, bnrm, bnrmto, n, n, b, ldb, ierr )
401*
402* Permute the matrices A, B to isolate eigenvalues if possible
403*
404 ileft = 1
405 iright = n + 1
406 iwrk = iright + n
407 CALL dggbal( 'P', n, a, lda, b, ldb, ilo, ihi, work( ileft ),
408 $ work( iright ), work( iwrk ), ierr )
409*
410* Reduce B to triangular form (QR decomposition of B)
411*
412 irows = ihi + 1 - ilo
413 IF( ilv ) THEN
414 icols = n + 1 - ilo
415 ELSE
416 icols = irows
417 END IF
418 itau = iwrk
419 iwrk = itau + irows
420 CALL dgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
421 $ work( iwrk ), lwork+1-iwrk, ierr )
422*
423* Apply the orthogonal transformation to matrix A
424*
425 CALL dormqr( 'L', 'T', irows, icols, irows, b( ilo, ilo ), ldb,
426 $ work( itau ), a( ilo, ilo ), lda, work( iwrk ),
427 $ lwork+1-iwrk, ierr )
428*
429* Initialize VL
430*
431 IF( ilvl ) THEN
432 CALL dlaset( 'Full', n, n, zero, one, vl, ldvl )
433 IF( irows.GT.1 ) THEN
434 CALL dlacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
435 $ vl( ilo+1, ilo ), ldvl )
436 END IF
437 CALL dorgqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
438 $ work( itau ), work( iwrk ), lwork+1-iwrk, ierr )
439 END IF
440*
441* Initialize VR
442*
443 IF( ilvr )
444 $ CALL dlaset( 'Full', n, n, zero, one, vr, ldvr )
445*
446* Reduce to generalized Hessenberg form
447*
448 IF( ilv ) THEN
449*
450* Eigenvectors requested -- work on whole matrix.
451*
452 CALL dgghd3( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
453 $ ldvl, vr, ldvr, work( iwrk ), lwork+1-iwrk, ierr )
454 ELSE
455 CALL dgghd3( 'N', 'N', irows, 1, irows, a( ilo, ilo ), lda,
456 $ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr,
457 $ work( iwrk ), lwork+1-iwrk, ierr )
458 END IF
459*
460* Perform QZ algorithm (Compute eigenvalues, and optionally, the
461* Schur forms and Schur vectors)
462*
463 iwrk = itau
464 IF( ilv ) THEN
465 chtemp = 'S'
466 ELSE
467 chtemp = 'E'
468 END IF
469 CALL dlaqz0( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
470 $ alphar, alphai, beta, vl, ldvl, vr, ldvr,
471 $ work( iwrk ), lwork+1-iwrk, 0, ierr )
472 IF( ierr.NE.0 ) THEN
473 IF( ierr.GT.0 .AND. ierr.LE.n ) THEN
474 info = ierr
475 ELSE IF( ierr.GT.n .AND. ierr.LE.2*n ) THEN
476 info = ierr - n
477 ELSE
478 info = n + 1
479 END IF
480 GO TO 110
481 END IF
482*
483* Compute Eigenvectors
484*
485 IF( ilv ) THEN
486 IF( ilvl ) THEN
487 IF( ilvr ) THEN
488 chtemp = 'B'
489 ELSE
490 chtemp = 'L'
491 END IF
492 ELSE
493 chtemp = 'R'
494 END IF
495 CALL dtgevc( chtemp, 'B', ldumma, n, a, lda, b, ldb, vl, ldvl,
496 $ vr, ldvr, n, in, work( iwrk ), ierr )
497 IF( ierr.NE.0 ) THEN
498 info = n + 2
499 GO TO 110
500 END IF
501*
502* Undo balancing on VL and VR and normalization
503*
504 IF( ilvl ) THEN
505 CALL dggbak( 'P', 'L', n, ilo, ihi, work( ileft ),
506 $ work( iright ), n, vl, ldvl, ierr )
507 DO 50 jc = 1, n
508 IF( alphai( jc ).LT.zero )
509 $ GO TO 50
510 temp = zero
511 IF( alphai( jc ).EQ.zero ) THEN
512 DO 10 jr = 1, n
513 temp = max( temp, abs( vl( jr, jc ) ) )
514 10 CONTINUE
515 ELSE
516 DO 20 jr = 1, n
517 temp = max( temp, abs( vl( jr, jc ) )+
518 $ abs( vl( jr, jc+1 ) ) )
519 20 CONTINUE
520 END IF
521 IF( temp.LT.smlnum )
522 $ GO TO 50
523 temp = one / temp
524 IF( alphai( jc ).EQ.zero ) THEN
525 DO 30 jr = 1, n
526 vl( jr, jc ) = vl( jr, jc )*temp
527 30 CONTINUE
528 ELSE
529 DO 40 jr = 1, n
530 vl( jr, jc ) = vl( jr, jc )*temp
531 vl( jr, jc+1 ) = vl( jr, jc+1 )*temp
532 40 CONTINUE
533 END IF
534 50 CONTINUE
535 END IF
536 IF( ilvr ) THEN
537 CALL dggbak( 'P', 'R', n, ilo, ihi, work( ileft ),
538 $ work( iright ), n, vr, ldvr, ierr )
539 DO 100 jc = 1, n
540 IF( alphai( jc ).LT.zero )
541 $ GO TO 100
542 temp = zero
543 IF( alphai( jc ).EQ.zero ) THEN
544 DO 60 jr = 1, n
545 temp = max( temp, abs( vr( jr, jc ) ) )
546 60 CONTINUE
547 ELSE
548 DO 70 jr = 1, n
549 temp = max( temp, abs( vr( jr, jc ) )+
550 $ abs( vr( jr, jc+1 ) ) )
551 70 CONTINUE
552 END IF
553 IF( temp.LT.smlnum )
554 $ GO TO 100
555 temp = one / temp
556 IF( alphai( jc ).EQ.zero ) THEN
557 DO 80 jr = 1, n
558 vr( jr, jc ) = vr( jr, jc )*temp
559 80 CONTINUE
560 ELSE
561 DO 90 jr = 1, n
562 vr( jr, jc ) = vr( jr, jc )*temp
563 vr( jr, jc+1 ) = vr( jr, jc+1 )*temp
564 90 CONTINUE
565 END IF
566 100 CONTINUE
567 END IF
568*
569* End of eigenvector calculation
570*
571 END IF
572*
573* Undo scaling if necessary
574*
575 110 CONTINUE
576*
577 IF( ilascl ) THEN
578 CALL dlascl( 'G', 0, 0, anrmto, anrm, n, 1, alphar, n, ierr )
579 CALL dlascl( 'G', 0, 0, anrmto, anrm, n, 1, alphai, n, ierr )
580 END IF
581*
582 IF( ilbscl ) THEN
583 CALL dlascl( 'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
584 END IF
585*
586 work( 1 ) = lwkopt
587 RETURN
588*
589* End of DGGEV3
590*
591 END
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:74
subroutine dlascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
DLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: dlascl.f:143
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 xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine dggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
DGGBAK
Definition: dggbak.f:147
subroutine dggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
DGGBAL
Definition: dggbal.f:177
subroutine dgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
DGEQRF
Definition: dgeqrf.f:146
subroutine dtgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, INFO)
DTGEVC
Definition: dtgevc.f:295
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 dggev3(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR, ALPHAI, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, INFO)
DGGEV3 computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices (...
Definition: dggev3.f:226
subroutine dorgqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
DORGQR
Definition: dorgqr.f:128
subroutine dormqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
DORMQR
Definition: dormqr.f:167
subroutine dgghd3(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, WORK, LWORK, INFO)
DGGHD3
Definition: dgghd3.f:230