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