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