SCALAPACK 2.2.2
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
pslahqr.f
Go to the documentation of this file.
1 SUBROUTINE pslahqr( WANTT, WANTZ, N, ILO, IHI, A, DESCA, WR, WI,
2 $ ILOZ, IHIZ, Z, DESCZ, WORK, LWORK, IWORK,
3 $ ILWORK, INFO )
4*
5* -- ScaLAPACK routine (version 2.0.2) --
6* Univ. of Tennessee, Univ. of California Berkeley, Univ. of Colorado Denver
7* May 1 2012
8*
9* .. Scalar Arguments ..
10 LOGICAL WANTT, WANTZ
11 INTEGER IHI, IHIZ, ILO, ILOZ, ILWORK, INFO, LWORK, N
12* ..
13* .. Array Arguments ..
14 INTEGER DESCA( * ), DESCZ( * ), IWORK( * )
15 REAL A( * ), WI( * ), WORK( * ), WR( * ), Z( * )
16* ..
17*
18* Purpose
19* =======
20*
21* PSLAHQR is an auxiliary routine used to find the Schur decomposition
22* and or eigenvalues of a matrix already in Hessenberg form from
23* cols ILO to IHI.
24*
25* Notes
26* =====
27*
28* Each global data object is described by an associated description
29* vector. This vector stores the information required to establish
30* the mapping between an object element and its corresponding process
31* and memory location.
32*
33* Let A be a generic term for any 2D block cyclicly distributed array.
34* Such a global array has an associated description vector DESCA.
35* In the following comments, the character _ should be read as
36* "of the global array".
37*
38* NOTATION STORED IN EXPLANATION
39* --------------- -------------- --------------------------------------
40* DTYPE_A(global) DESCA( DTYPE_ )The descriptor type. In this case,
41* DTYPE_A = 1.
42* CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
43* the BLACS process grid A is distribu-
44* ted over. The context itself is glo-
45* bal, but the handle (the integer
46* value) may vary.
47* M_A (global) DESCA( M_ ) The number of rows in the global
48* array A.
49* N_A (global) DESCA( N_ ) The number of columns in the global
50* array A.
51* MB_A (global) DESCA( MB_ ) The blocking factor used to distribute
52* the rows of the array.
53* NB_A (global) DESCA( NB_ ) The blocking factor used to distribute
54* the columns of the array.
55* RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
56* row of the array A is distributed.
57* CSRC_A (global) DESCA( CSRC_ ) The process column over which the
58* first column of the array A is
59* distributed.
60* LLD_A (local) DESCA( LLD_ ) The leading dimension of the local
61* array. LLD_A >= MAX(1,LOCr(M_A)).
62*
63* Let K be the number of rows or columns of a distributed matrix,
64* and assume that its process grid has dimension p x q.
65* LOCr( K ) denotes the number of elements of K that a process
66* would receive if K were distributed over the p processes of its
67* process column.
68* Similarly, LOCc( K ) denotes the number of elements of K that a
69* process would receive if K were distributed over the q processes of
70* its process row.
71* The values of LOCr() and LOCc() may be determined via a call to the
72* ScaLAPACK tool function, NUMROC:
73* LOCr( M ) = NUMROC( M, MB_A, MYROW, RSRC_A, NPROW ),
74* LOCc( N ) = NUMROC( N, NB_A, MYCOL, CSRC_A, NPCOL ).
75* An upper bound for these quantities may be computed by:
76* LOCr( M ) <= ceil( ceil(M/MB_A)/NPROW )*MB_A
77* LOCc( N ) <= ceil( ceil(N/NB_A)/NPCOL )*NB_A
78*
79* Arguments
80* =========
81*
82* WANTT (global input) LOGICAL
83* = .TRUE. : the full Schur form T is required;
84* = .FALSE.: only eigenvalues are required.
85*
86* WANTZ (global input) LOGICAL
87* = .TRUE. : the matrix of Schur vectors Z is required;
88* = .FALSE.: Schur vectors are not required.
89*
90* N (global input) INTEGER
91* The order of the Hessenberg matrix A (and Z if WANTZ).
92* N >= 0.
93*
94* ILO (global input) INTEGER
95* IHI (global input) INTEGER
96* It is assumed that A is already upper quasi-triangular in
97* rows and columns IHI+1:N, and that A(ILO,ILO-1) = 0 (unless
98* ILO = 1). PSLAHQR works primarily with the Hessenberg
99* submatrix in rows and columns ILO to IHI, but applies
100* transformations to all of H if WANTT is .TRUE..
101* 1 <= ILO <= max(1,IHI); IHI <= N.
102*
103* A (global input/output) REAL array, dimension
104* (DESCA(LLD_),*)
105* On entry, the upper Hessenberg matrix A.
106* On exit, if WANTT is .TRUE., A is upper quasi-triangular in
107* rows and columns ILO:IHI, with any 2-by-2 or larger diagonal
108* blocks not yet in standard form. If WANTT is .FALSE., the
109* contents of A are unspecified on exit.
110*
111* DESCA (global and local input) INTEGER array of dimension DLEN_.
112* The array descriptor for the distributed matrix A.
113*
114* WR (global replicated output) REAL array,
115* dimension (N)
116* WI (global replicated output) REAL array,
117* dimension (N)
118* The real and imaginary parts, respectively, of the computed
119* eigenvalues ILO to IHI are stored in the corresponding
120* elements of WR and WI. If two eigenvalues are computed as a
121* complex conjugate pair, they are stored in consecutive
122* elements of WR and WI, say the i-th and (i+1)th, with
123* WI(i) > 0 and WI(i+1) < 0. If WANTT is .TRUE., the
124* eigenvalues are stored in the same order as on the diagonal
125* of the Schur form returned in A. A may be returned with
126* larger diagonal blocks until the next release.
127*
128* ILOZ (global input) INTEGER
129* IHIZ (global input) INTEGER
130* Specify the rows of Z to which transformations must be
131* applied if WANTZ is .TRUE..
132* 1 <= ILOZ <= ILO; IHI <= IHIZ <= N.
133*
134* Z (global input/output) REAL array.
135* If WANTZ is .TRUE., on entry Z must contain the current
136* matrix Z of transformations accumulated by PDHSEQR, and on
137* exit Z has been updated; transformations are applied only to
138* the submatrix Z(ILOZ:IHIZ,ILO:IHI).
139* If WANTZ is .FALSE., Z is not referenced.
140*
141* DESCZ (global and local input) INTEGER array of dimension DLEN_.
142* The array descriptor for the distributed matrix Z.
143*
144* WORK (local output) REAL array of size LWORK
145*
146* LWORK (local input) INTEGER
147* WORK(LWORK) is a local array and LWORK is assumed big enough
148* so that LWORK >= 3*N +
149* MAX( 2*MAX(DESCZ(LLD_),DESCA(LLD_)) + 2*LOCc(N),
150* 7*Ceil(N/HBL)/LCM(NPROW,NPCOL)) )
151*
152* IWORK (global and local input) INTEGER array of size ILWORK
153*
154* ILWORK (local input) INTEGER
155* This holds the some of the IBLK integer arrays. This is held
156* as a place holder for the next release.
157*
158* INFO (global output) INTEGER
159* < 0: parameter number -INFO incorrect or inconsistent
160* = 0: successful exit
161* > 0: PSLAHQR failed to compute all the eigenvalues ILO to IHI
162* in a total of 30*(IHI-ILO+1) iterations; if INFO = i,
163* elements i+1:ihi of WR and WI contain those eigenvalues
164* which have been successfully computed.
165*
166* Logic:
167* This algorithm is very similar to _LAHQR. Unlike _LAHQR,
168* instead of sending one double shift through the largest
169* unreduced submatrix, this algorithm sends multiple double shifts
170* and spaces them apart so that there can be parallelism across
171* several processor row/columns. Another critical difference is
172* that this algorithm aggregrates multiple transforms together in
173* order to apply them in a block fashion.
174*
175* Important Local Variables:
176* IBLK = The maximum number of bulges that can be computed.
177* Currently fixed. Future releases this won't be fixed.
178* HBL = The square block size (HBL=DESCA(MB_)=DESCA(NB_))
179* ROTN = The number of transforms to block together
180* NBULGE = The number of bulges that will be attempted on the
181* current submatrix.
182* IBULGE = The current number of bulges started.
183* K1(*),K2(*) = The current bulge loops from K1(*) to K2(*).
184*
185* Subroutines:
186* This routine calls:
187* PSLACONSB -> To determine where to start each iteration
188* PSLAWIL -> Given the shift, get the transformation
189* SLASORTE -> Pair up eigenvalues so that reals are paired.
190* PSLACP3 -> Parallel array to local replicated array copy &
191* back.
192* SLAREF -> Row/column reflector applier. Core routine
193* here.
194* PSLASMSUB -> Finds negligible subdiagonal elements.
195*
196* Current Notes and/or Restrictions:
197* 1.) This code requires the distributed block size to be square
198* and at least six (6); unlike simpler codes like LU, this
199* algorithm is extremely sensitive to block size. Unwise
200* choices of too small a block size can lead to bad
201* performance.
202* 2.) This code requires A and Z to be distributed identically
203* and have identical contxts.
204* 3.) This release currently does not have a routine for
205* resolving the Schur blocks into regular 2x2 form after
206* this code is completed. Because of this, a significant
207* performance impact is required while the deflation is done
208* by sometimes a single column of processors.
209* 4.) This code does not currently block the initial transforms
210* so that none of the rows or columns for any bulge are
211* completed until all are started. To offset pipeline
212* start-up it is recommended that at least 2*LCM(NPROW,NPCOL)
213* bulges are used (if possible)
214* 5.) The maximum number of bulges currently supported is fixed at
215* 32. In future versions this will be limited only by the
216* incoming WORK array.
217* 6.) The matrix A must be in upper Hessenberg form. If elements
218* below the subdiagonal are nonzero, the resulting transforms
219* may be nonsimilar. This is also true with the LAPACK
220* routine.
221* 7.) For this release, it is assumed RSRC_=CSRC_=0
222* 8.) Currently, all the eigenvalues are distributed to all the
223* nodes. Future releases will probably distribute the
224* eigenvalues by the column partitioning.
225* 9.) The internals of this routine are subject to change.
226*
227* Implemented by: G. Henry, November 17, 1996
228*
229* =====================================================================
230*
231* .. Parameters ..
232 INTEGER BLOCK_CYCLIC_2D, CSRC_, CTXT_, DLEN_, DTYPE_,
233 $ LLD_, MB_, M_, NB_, N_, RSRC_
234 parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
235 $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
236 $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
237 REAL ZERO, ONE, HALF
238 PARAMETER ( ZERO = 0.0, one = 1.0, half = 0.5 )
239 REAL CONST
240 parameter( const = 1.50 )
241 INTEGER IBLK
242 parameter( iblk = 32 )
243* ..
244* .. Local Scalars ..
245 INTEGER CONTXT, DOWN, HBL, I, I1, I2, IAFIRST, IBULGE,
246 $ ICBUF, ICOL, ICOL1, ICOL2, IDIA, IERR, II,
247 $ irbuf, irow, irow1, irow2, ispec, istart,
248 $ istartcol, istartrow, istop, isub, isup,
249 $ itermax, itmp1, itmp2, itn, its, j, jafirst,
250 $ jblk, jj, k, ki, l, lcmrc, lda, ldz, left,
251 $ lihih, lihiz, liloh, liloz, locali1, locali2,
252 $ localk, localm, m, modkm1, mycol, myrow,
253 $ nbulge, nh, node, npcol, nprow, nr, num, nz,
254 $ right, rotn, up, vecsidx
255 REAL AVE, DISC, H00, H10, H11, H12, H21, H22, H33,
256 $ H43H34, H44, OVFL, S, SMLNUM, SUM, T1, T1COPY,
257 $ t2, t3, ulp, unfl, v1save, v2, v2save, v3,
258 $ v3save, cs, sn
259* ..
260* .. Local Arrays ..
261 INTEGER ICURCOL( IBLK ), ICURROW( IBLK ), K1( IBLK ),
262 $ K2( IBLK ), KCOL( IBLK ), KP2COL( IBLK ),
263 $ kp2row( iblk ), krow( iblk ), localk2( iblk )
264 REAL S1( 2*IBLK, 2*IBLK ), SMALLA( 6, 6, IBLK ),
265 $ VCOPY( 3 )
266* ..
267* .. External Functions ..
268 INTEGER ILCM, NUMROC
269 REAL PSLAMCH
270 EXTERNAL ilcm, numroc, pslamch
271* ..
272* .. External Subroutines ..
273 EXTERNAL blacs_gridinfo, scopy, sgebr2d, sgebs2d,
274 $ sgerv2d, sgesd2d, sgsum2d, slahqr, slaref,
275 $ slarfg, slasorte, igamn2d, infog1l, infog2l,
277 $ pslawil, pxerbla, slanv2
278* ..
279* .. Intrinsic Functions ..
280 INTRINSIC abs, max, min, mod, sign, sqrt
281* ..
282* .. Executable Statements ..
283*
284 info = 0
285*
286 itermax = 30*( ihi-ilo+1 )
287* ITERMAX = 0
288 IF( n.EQ.0 )
289 $ RETURN
290*
291* NODE (IAFIRST,JAFIRST) OWNS A(1,1)
292*
293 hbl = desca( mb_ )
294 contxt = desca( ctxt_ )
295 lda = desca( lld_ )
296 iafirst = desca( rsrc_ )
297 jafirst = desca( csrc_ )
298 ldz = descz( lld_ )
299 CALL blacs_gridinfo( contxt, nprow, npcol, myrow, mycol )
300 node = myrow*npcol + mycol
301 num = nprow*npcol
302 left = mod( mycol+npcol-1, npcol )
303 right = mod( mycol+1, npcol )
304 up = mod( myrow+nprow-1, nprow )
305 down = mod( myrow+1, nprow )
306 lcmrc = ilcm( nprow, npcol )
307*
308* Determine the number of columns we have so we can check workspace
309*
310 localk = numroc( n, hbl, mycol, jafirst, npcol )
311 jj = n / hbl
312 IF( jj*hbl.LT.n )
313 $ jj = jj + 1
314 jj = 7*jj / lcmrc
315 IF( lwork.LT.3*n+max( 2*max( lda, ldz )+2*localk, jj ) ) THEN
316 info = -15
317 END IF
318 IF( descz( ctxt_ ).NE.desca( ctxt_ ) ) THEN
319 info = -( 1300+ctxt_ )
320 END IF
321 IF( desca( mb_ ).NE.desca( nb_ ) ) THEN
322 info = -( 700+nb_ )
323 END IF
324 IF( descz( mb_ ).NE.descz( nb_ ) ) THEN
325 info = -( 1300+nb_ )
326 END IF
327 IF( desca( mb_ ).NE.descz( mb_ ) ) THEN
328 info = -( 1300+mb_ )
329 END IF
330 IF( ( desca( rsrc_ ).NE.0 ) .OR. ( desca( csrc_ ).NE.0 ) ) THEN
331 info = -( 700+rsrc_ )
332 END IF
333 IF( ( descz( rsrc_ ).NE.0 ) .OR. ( descz( csrc_ ).NE.0 ) ) THEN
334 info = -( 1300+rsrc_ )
335 END IF
336 IF( ( ilo.GT.n ) .OR. ( ilo.LT.1 ) ) THEN
337 info = -4
338 END IF
339 IF( ( ihi.GT.n ) .OR. ( ihi.LT.1 ) ) THEN
340 info = -5
341 END IF
342 IF( hbl.LT.5 ) THEN
343 info = -( 700+mb_ )
344 END IF
345 CALL igamn2d( contxt, 'ALL', ' ', 1, 1, info, 1, itmp1, itmp2, -1,
346 $ -1, -1 )
347 IF( info.LT.0 ) THEN
348 CALL pxerbla( contxt, 'PSLAHQR', -info )
349 RETURN
350 END IF
351*
352* Set work array indices
353*
354 vecsidx = 0
355 idia = 3*n
356 isub = 3*n
357 isup = 3*n
358 irbuf = 3*n
359 icbuf = 3*n
360*
361* Find a value for ROTN
362*
363 rotn = hbl / 3
364 rotn = max( rotn, hbl-2 )
365 rotn = min( rotn, 1 )
366*
367 IF( ilo.EQ.ihi ) THEN
368 CALL infog2l( ilo, ilo, desca, nprow, npcol, myrow, mycol,
369 $ irow, icol, ii, jj )
370 IF( ( myrow.EQ.ii ) .AND. ( mycol.EQ.jj ) ) THEN
371 wr( ilo ) = a( ( icol-1 )*lda+irow )
372 ELSE
373 wr( ilo ) = zero
374 END IF
375 wi( ilo ) = zero
376 RETURN
377 END IF
378*
379 nh = ihi - ilo + 1
380 nz = ihiz - iloz + 1
381*
382 CALL infog1l( iloz, hbl, nprow, myrow, 0, liloz, lihiz )
383 lihiz = numroc( ihiz, hbl, myrow, 0, nprow )
384*
385* Set machine-dependent constants for the stopping criterion.
386* If NORM(H) <= SQRT(OVFL), overflow should not occur.
387*
388 unfl = pslamch( contxt, 'SAFE MINIMUM' )
389 ovfl = one / unfl
390 CALL pslabad( contxt, unfl, ovfl )
391 ulp = pslamch( contxt, 'PRECISION' )
392 smlnum = unfl*( nh / ulp )
393*
394* I1 and I2 are the indices of the first row and last column of H
395* to which transformations must be applied. If eigenvalues only are
396* being computed, I1 and I2 are set inside the main loop.
397*
398 IF( wantt ) THEN
399 i1 = 1
400 i2 = n
401 END IF
402*
403* ITN is the total number of QR iterations allowed.
404*
405 itn = itermax
406*
407* The main loop begins here. I is the loop index and decreases from
408* IHI to ILO in steps of our schur block size (<=2*IBLK). Each
409* iteration of the loop works with the active submatrix in rows
410* and columns L to I. Eigenvalues I+1 to IHI have already
411* converged. Either L = ILO or the global A(L,L-1) is negligible
412* so that the matrix splits.
413*
414 i = ihi
415 10 CONTINUE
416 l = ilo
417 IF( i.LT.ilo )
418 $ GO TO 450
419*
420* Perform QR iterations on rows and columns ILO to I until a
421* submatrix of order 1 or 2 splits off at the bottom because a
422* subdiagonal element has become negligible.
423*
424 DO 420 its = 0, itn
425*
426* Look for a single small subdiagonal element.
427*
428 CALL pslasmsub( a, desca, i, l, k, smlnum, work( irbuf+1 ),
429 $ lwork-irbuf )
430 l = k
431*
432 IF( l.GT.ilo ) THEN
433*
434* H(L,L-1) is negligible
435*
436 CALL infog2l( l, l-1, desca, nprow, npcol, myrow, mycol,
437 $ irow, icol, itmp1, itmp2 )
438 IF( ( myrow.EQ.itmp1 ) .AND. ( mycol.EQ.itmp2 ) ) THEN
439 a( ( icol-1 )*lda+irow ) = zero
440 END IF
441 work( isub+l-1 ) = zero
442 END IF
443*
444* Exit from loop if a submatrix of order 1 or 2 has split off.
445*
446 m = l - 10
447* IF ( L .GE. I - (2*IBLK-1) )
448* IF ( L .GE. I - MAX(2*IBLK-1,HBL) )
449 IF( l.GE.i-1 )
450 $ GO TO 430
451*
452* Now the active submatrix is in rows and columns L to I. If
453* eigenvalues only are being computed, only the active submatrix
454* need be transformed.
455*
456 IF( .NOT.wantt ) THEN
457 i1 = l
458 i2 = i
459 END IF
460*
461* Copy submatrix of size 2*JBLK and prepare to do generalized
462* Wilkinson shift or an exceptional shift
463*
464 jblk = min( iblk, ( ( i-l+1 ) / 2 )-1 )
465 IF( jblk.GT.lcmrc ) THEN
466*
467* Make sure it's divisible by LCM (we want even workloads!)
468*
469 jblk = jblk - mod( jblk, lcmrc )
470 END IF
471 jblk = min( jblk, 2*lcmrc )
472 jblk = max( jblk, 1 )
473*
474 CALL pslacp3( 2*jblk, i-2*jblk+1, a, desca, s1, 2*iblk, -1, -1,
475 $ 0 )
476 IF( its.EQ.20 .OR. its.EQ.40 ) THEN
477*
478* Exceptional shift.
479*
480 DO 20 ii = 2*jblk, 2, -1
481 s1( ii, ii ) = const*( abs( s1( ii, ii ) )+
482 $ abs( s1( ii, ii-1 ) ) )
483 s1( ii, ii-1 ) = zero
484 s1( ii-1, ii ) = zero
485 20 CONTINUE
486 s1( 1, 1 ) = const*abs( s1( 1, 1 ) )
487 ELSE
488 CALL slahqr( .false., .false., 2*jblk, 1, 2*jblk, s1,
489 $ 2*iblk, work( irbuf+1 ), work( icbuf+1 ), 1,
490 $ 2*jblk, z, ldz, ierr )
491*
492* Prepare to use Wilkinson's double shift
493*
494 h44 = s1( 2*jblk, 2*jblk )
495 h33 = s1( 2*jblk-1, 2*jblk-1 )
496 h43h34 = s1( 2*jblk-1, 2*jblk )*s1( 2*jblk, 2*jblk-1 )
497 IF( ( jblk.GT.1 ) .AND. ( its.GT.30 ) ) THEN
498 s = s1( 2*jblk-1, 2*jblk-2 )
499 disc = ( h33-h44 )*half
500 disc = disc*disc + h43h34
501 IF( disc.GT.zero ) THEN
502*
503* Real roots: Use Wilkinson's shift twice
504*
505 disc = sqrt( disc )
506 ave = half*( h33+h44 )
507 IF( abs( h33 )-abs( h44 ).GT.zero ) THEN
508 h33 = h33*h44 - h43h34
509 h44 = h33 / ( sign( disc, ave )+ave )
510 ELSE
511 h44 = sign( disc, ave ) + ave
512 END IF
513 h33 = h44
514 h43h34 = zero
515 END IF
516 END IF
517 END IF
518*
519* Look for two consecutive small subdiagonal elements:
520* PSLACONSB is the routine that does this.
521*
522c CALL PSLACONSB( A, DESCA, I, L, M, H44, H33, H43H34,
523c $ WORK( IRBUF+1 ), LWORK-IRBUF )
524*
525* Skip small submatrices
526*
527* IF ( M .GE. I - 5 )
528* $ GO TO 80
529*
530* In principle PSLACONSB needs to check all shifts to decide
531* whether two consecutive small subdiagonal entries are suitable
532* as the starting position of the bulge chasing phase. It can be
533* dangerous to check the first pair of shifts only. Moreover it
534* is quite rare to obtain an M which is much larger than L. This
535* process is a bit expensive compared with the benefit.
536* Therefore it is sensible to abandon this routine. Total amount
537* of communications is saved in average.
538*
539 m = l
540* Double-shift QR step
541*
542* NBULGE is the number of bulges that will be attempted
543*
544 istop = min( m+rotn-mod( m, rotn ), i-2 )
545 istop = min( istop, m+hbl-3-mod( m-1, hbl ) )
546 istop = min( istop, i2-2 )
547 istop = max( istop, m )
548 nbulge = ( i-1-istop ) / hbl
549*
550* Do not exceed maximum determined.
551*
552 nbulge = min( nbulge, jblk )
553 IF( nbulge.GT.lcmrc ) THEN
554*
555* Make sure it's divisible by LCM (we want even workloads!)
556*
557 nbulge = nbulge - mod( nbulge, lcmrc )
558 END IF
559 nbulge = max( nbulge, 1 )
560*
561 IF( ( its.NE.20 ) .AND. ( its.NE.40 ) .AND. ( nbulge.GT.1 ) )
562 $ THEN
563*
564* sort the eigenpairs so that they are in twos for double
565* shifts. only call if several need sorting
566*
567 CALL slasorte( s1( 2*( jblk-nbulge )+1,
568 $ 2*( jblk-nbulge )+1 ), 2*iblk, 2*nbulge,
569 $ work( irbuf+1 ), ierr )
570 END IF
571*
572* IBULGE is the number of bulges going so far
573*
574 ibulge = 1
575*
576* "A" row defs : main row transforms from LOCALK to LOCALI2
577*
578 CALL infog1l( m, hbl, npcol, mycol, 0, itmp1, localk )
579 localk = numroc( n, hbl, mycol, 0, npcol )
580 CALL infog1l( 1, hbl, npcol, mycol, 0, icol1, locali2 )
581 locali2 = numroc( i2, hbl, mycol, 0, npcol )
582*
583* "A" col defs : main col transforms from LOCALI1 to LOCALM
584*
585 CALL infog1l( i1, hbl, nprow, myrow, 0, locali1, icol1 )
586 icol1 = numroc( n, hbl, myrow, 0, nprow )
587 CALL infog1l( 1, hbl, nprow, myrow, 0, localm, icol1 )
588 icol1 = numroc( min( m+3, i ), hbl, myrow, 0, nprow )
589*
590* Which row & column will start the bulges
591*
592 istartrow = mod( ( m+1 ) / hbl, nprow ) + iafirst
593 istartcol = mod( ( m+1 ) / hbl, npcol ) + jafirst
594*
595 CALL infog1l( m, hbl, nprow, myrow, 0, ii, itmp2 )
596 itmp2 = numroc( n, hbl, myrow, 0, nprow )
597 CALL infog1l( m, hbl, npcol, mycol, 0, jj, itmp2 )
598 itmp2 = numroc( n, hbl, mycol, 0, npcol )
599 CALL infog1l( 1, hbl, nprow, myrow, 0, istop, kp2row( 1 ) )
600 kp2row( 1 ) = numroc( m+2, hbl, myrow, 0, nprow )
601 CALL infog1l( 1, hbl, npcol, mycol, 0, istop, kp2col( 1 ) )
602 kp2col( 1 ) = numroc( m+2, hbl, mycol, 0, npcol )
603*
604* Set all values for bulges. All bulges are stored in
605* intermediate steps as loops over KI. Their current "task"
606* over the global M to I-1 values is always K1(KI) to K2(KI).
607* However, because there are many bulges, K1(KI) & K2(KI) might
608* go past that range while later bulges (KI+1,KI+2,etc..) are
609* finishing up.
610*
611* Rules:
612* If MOD(K1(KI)-1,HBL) < HBL-2 then MOD(K2(KI)-1,HBL)<HBL-2
613* If MOD(K1(KI)-1,HBL) = HBL-2 then MOD(K2(KI)-1,HBL)=HBL-2
614* If MOD(K1(KI)-1,HBL) = HBL-1 then MOD(K2(KI)-1,HBL)=HBL-1
615* K2(KI)-K1(KI) <= ROTN
616*
617* We first hit a border when MOD(K1(KI)-1,HBL)=HBL-2 and we hit
618* it again when MOD(K1(KI)-1,HBL)=HBL-1.
619*
620 DO 30 ki = 1, nbulge
621 k1( ki ) = m
622 istop = min( m+rotn-mod( m, rotn ), i-2 )
623 istop = min( istop, m+hbl-3-mod( m-1, hbl ) )
624 istop = min( istop, i2-2 )
625 istop = max( istop, m )
626 k2( ki ) = istop
627 icurrow( ki ) = istartrow
628 icurcol( ki ) = istartcol
629 localk2( ki ) = itmp1
630 krow( ki ) = ii
631 kcol( ki ) = jj
632 IF( ki.GT.1 )
633 $ kp2row( ki ) = kp2row( 1 )
634 IF( ki.GT.1 )
635 $ kp2col( ki ) = kp2col( 1 )
636 30 CONTINUE
637*
638* Get first transform on node who owns M+2,M+2
639*
640 DO 31 itmp1 = 1, 3
641 vcopy(itmp1) = zero
642 31 CONTINUE
643 itmp1 = istartrow
644 itmp2 = istartcol
645 CALL pslawil( itmp1, itmp2, m, a, desca, h44, h33, h43h34,
646 $ vcopy )
647 v1save = vcopy( 1 )
648 v2save = vcopy( 2 )
649 v3save = vcopy( 3 )
650 IF( k2( ibulge ).LE.i-1 ) THEN
651 40 CONTINUE
652 IF( ( k1( ibulge ).GE.m+5 ) .AND. ( ibulge.LT.nbulge ) )
653 $ THEN
654 IF( ( mod( k2( ibulge )+2, hbl ).EQ.mod( k2( ibulge+1 )+
655 $ 2, hbl ) ) .AND. ( k1( 1 ).LE.i-1 ) ) THEN
656 h44 = s1( 2*jblk-2*ibulge, 2*jblk-2*ibulge )
657 h33 = s1( 2*jblk-2*ibulge-1, 2*jblk-2*ibulge-1 )
658 h43h34 = s1( 2*jblk-2*ibulge-1, 2*jblk-2*ibulge )*
659 $ s1( 2*jblk-2*ibulge, 2*jblk-2*ibulge-1 )
660 itmp1 = istartrow
661 itmp2 = istartcol
662 CALL pslawil( itmp1, itmp2, m, a, desca, h44, h33,
663 $ h43h34, vcopy )
664 v1save = vcopy( 1 )
665 v2save = vcopy( 2 )
666 v3save = vcopy( 3 )
667 ibulge = ibulge + 1
668 END IF
669 END IF
670*
671* When we hit a border, there are row and column transforms that
672* overlap over several processors and the code gets very
673* "congested." As a remedy, when we first hit a border, a 6x6
674* *local* matrix is generated on one node (called SMALLA) and
675* work is done on that. At the end of the border, the data is
676* passed back and everything stays a lot simpler.
677*
678 DO 80 ki = 1, ibulge
679*
680 istart = max( k1( ki ), m )
681 istop = min( k2( ki ), i-1 )
682 k = istart
683 modkm1 = mod( k-1, hbl )
684 IF( ( modkm1.GE.hbl-2 ) .AND. ( k.LE.i-1 ) ) THEN
685 DO 81 itmp1 = 1, 6
686 DO 82 itmp2 = 1, 6
687 smalla(itmp1, itmp2, ki) = zero
688 82 CONTINUE
689 81 CONTINUE
690 IF( ( modkm1.EQ.hbl-2 ) .AND. ( k.LT.i-1 ) ) THEN
691*
692* Copy 6 elements from global A(K-1:K+4,K-1:K+4)
693*
694 CALL infog2l( k+2, k+2, desca, nprow, npcol, myrow,
695 $ mycol, irow1, icol1, itmp1, itmp2 )
696 CALL pslacp3( min( 6, n-k+2 ), k-1, a, desca,
697 $ smalla( 1, 1, ki ), 6, itmp1, itmp2,
698 $ 0 )
699 END IF
700 IF( modkm1.EQ.hbl-1 ) THEN
701*
702* Copy 6 elements from global A(K-2:K+3,K-2:K+3)
703*
704 CALL infog2l( k+1, k+1, desca, nprow, npcol, myrow,
705 $ mycol, irow1, icol1, itmp1, itmp2 )
706 CALL pslacp3( min( 6, n-k+3 ), k-2, a, desca,
707 $ smalla( 1, 1, ki ), 6, itmp1, itmp2,
708 $ 0 )
709 END IF
710 END IF
711*
712* SLAHQR used to have a single row application and a single
713* column application to H. Here we do something a little
714* more clever. We break each transformation down into 3
715* parts:
716* 1.) The minimum amount of work it takes to determine
717* a group of ROTN transformations (this is on
718* the critical path.) (Loops 130-180)
719* 2.) The small work it takes so that each of the rows
720* and columns is at the same place. For example,
721* all ROTN row transforms are all complete
722* through some column TMP. (Loops within 190)
723* 3.) The majority of the row and column transforms
724* are then applied in a block fashion.
725* (Loops 290 on.)
726*
727* Each of these three parts are further subdivided into 3
728* parts:
729* A.) Work at the start of a border when
730* MOD(ISTART-1,HBL) = HBL-2
731* B.) Work at the end of a border when
732* MOD(ISTART-1,HBL) = HBL-1
733* C.) Work in the middle of the block when
734* MOD(ISTART-1,HBL) < HBL-2
735*
736 IF( ( myrow.EQ.icurrow( ki ) ) .AND.
737 $ ( mycol.EQ.icurcol( ki ) ) .AND.
738 $ ( modkm1.EQ.hbl-2 ) .AND.
739 $ ( istart.LT.min( i-1, istop+1 ) ) ) THEN
740 k = istart
741 nr = min( 3, i-k+1 )
742 IF( k.GT.m ) THEN
743 CALL scopy( nr, smalla( 2, 1, ki ), 1, vcopy, 1 )
744 ELSE
745 vcopy( 1 ) = v1save
746 vcopy( 2 ) = v2save
747 vcopy( 3 ) = v3save
748 END IF
749 CALL slarfg( nr, vcopy( 1 ), vcopy( 2 ), 1, t1copy )
750 IF( k.GT.m ) THEN
751 smalla( 2, 1, ki ) = vcopy( 1 )
752 smalla( 3, 1, ki ) = zero
753 IF( k.LT.i-1 )
754 $ smalla( 4, 1, ki ) = zero
755 ELSE IF( m.GT.l ) THEN
756 smalla( 2, 1, ki ) = -smalla( 2, 1, ki )
757 END IF
758 v2 = vcopy( 2 )
759 t2 = t1copy*v2
760 work( vecsidx+( k-1 )*3+1 ) = vcopy( 2 )
761 work( vecsidx+( k-1 )*3+2 ) = vcopy( 3 )
762 work( vecsidx+( k-1 )*3+3 ) = t1copy
763 END IF
764*
765 IF( ( mod( istop-1, hbl ).EQ.hbl-1 ) .AND.
766 $ ( myrow.EQ.icurrow( ki ) ) .AND.
767 $ ( mycol.EQ.icurcol( ki ) ) .AND.
768 $ ( istart.LE.min( i, istop ) ) ) THEN
769 k = istart
770 nr = min( 3, i-k+1 )
771 IF( k.GT.m ) THEN
772 CALL scopy( nr, smalla( 3, 2, ki ), 1, vcopy, 1 )
773 ELSE
774 vcopy( 1 ) = v1save
775 vcopy( 2 ) = v2save
776 vcopy( 3 ) = v3save
777 END IF
778 CALL slarfg( nr, vcopy( 1 ), vcopy( 2 ), 1, t1copy )
779 IF( k.GT.m ) THEN
780 smalla( 3, 2, ki ) = vcopy( 1 )
781 smalla( 4, 2, ki ) = zero
782 IF( k.LT.i-1 )
783 $ smalla( 5, 2, ki ) = zero
784*
785* Set a subdiagonal to zero now if it's possible
786*
787* H11 = SMALLA(1,1,KI)
788* H10 = SMALLA(2,1,KI)
789* H22 = SMALLA(2,2,KI)
790* IF ( ABS(H10) .LE. MAX(ULP*(ABS(H11)+ABS(H22)),
791* $ SMLNUM) ) THEN
792* SMALLA(2,1,KI) = ZERO
793* WORK(ISUB+K-2) = ZERO
794* END IF
795 ELSE IF( m.GT.l ) THEN
796 smalla( 3, 2, ki ) = -smalla( 3, 2, ki )
797 END IF
798 v2 = vcopy( 2 )
799 t2 = t1copy*v2
800 work( vecsidx+( k-1 )*3+1 ) = vcopy( 2 )
801 work( vecsidx+( k-1 )*3+2 ) = vcopy( 3 )
802 work( vecsidx+( k-1 )*3+3 ) = t1copy
803 END IF
804*
805 IF( ( modkm1.EQ.0 ) .AND. ( istart.LE.i-1 ) .AND.
806 $ ( myrow.EQ.icurrow( ki ) ) .AND.
807 $ ( right.EQ.icurcol( ki ) ) ) THEN
808*
809* (IROW1,ICOL1) is (I,J)-coordinates of H(ISTART,ISTART)
810*
811 irow1 = krow( ki )
812 icol1 = localk2( ki )
813 IF( istart.GT.m ) THEN
814 vcopy( 1 ) = smalla( 4, 3, ki )
815 vcopy( 2 ) = smalla( 5, 3, ki )
816 vcopy( 3 ) = smalla( 6, 3, ki )
817 nr = min( 3, i-istart+1 )
818 CALL slarfg( nr, vcopy( 1 ), vcopy( 2 ), 1,
819 $ t1copy )
820 a( ( icol1-2 )*lda+irow1 ) = vcopy( 1 )
821 a( ( icol1-2 )*lda+irow1+1 ) = zero
822 IF( istart.LT.i-1 ) THEN
823 a( ( icol1-2 )*lda+irow1+2 ) = zero
824 END IF
825 ELSE
826 IF( m.GT.l ) THEN
827 a( ( icol1-2 )*lda+irow1 ) = -a( ( icol1-2 )*
828 $ lda+irow1 )
829 END IF
830 END IF
831 END IF
832*
833 IF( ( myrow.EQ.icurrow( ki ) ) .AND.
834 $ ( mycol.EQ.icurcol( ki ) ) .AND.
835 $ ( ( ( modkm1.EQ.hbl-2 ) .AND. ( istart.EQ.i-
836 $ 1 ) ) .OR. ( ( modkm1.LT.hbl-2 ) .AND. ( istart.LE.i-
837 $ 1 ) ) ) ) THEN
838*
839* (IROW1,ICOL1) is (I,J)-coordinates of H(ISTART,ISTART)
840*
841 irow1 = krow( ki )
842 icol1 = localk2( ki )
843 DO 70 k = istart, istop
844*
845* Create and do these transforms
846*
847 nr = min( 3, i-k+1 )
848 IF( k.GT.m ) THEN
849 IF( mod( k-1, hbl ).EQ.0 ) THEN
850 vcopy( 1 ) = smalla( 4, 3, ki )
851 vcopy( 2 ) = smalla( 5, 3, ki )
852 vcopy( 3 ) = smalla( 6, 3, ki )
853 ELSE
854 vcopy( 1 ) = a( ( icol1-2 )*lda+irow1 )
855 vcopy( 2 ) = a( ( icol1-2 )*lda+irow1+1 )
856 IF( nr.EQ.3 ) THEN
857 vcopy( 3 ) = a( ( icol1-2 )*lda+irow1+2 )
858 END IF
859 END IF
860 ELSE
861 vcopy( 1 ) = v1save
862 vcopy( 2 ) = v2save
863 vcopy( 3 ) = v3save
864 END IF
865 CALL slarfg( nr, vcopy( 1 ), vcopy( 2 ), 1,
866 $ t1copy )
867 IF( k.GT.m ) THEN
868 IF( mod( k-1, hbl ).GT.0 ) THEN
869 a( ( icol1-2 )*lda+irow1 ) = vcopy( 1 )
870 a( ( icol1-2 )*lda+irow1+1 ) = zero
871 IF( k.LT.i-1 ) THEN
872 a( ( icol1-2 )*lda+irow1+2 ) = zero
873 END IF
874*
875* Set a subdiagonal to zero now if it's possible
876*
877* IF ( (IROW1.GT.2) .AND. (ICOL1.GT.2) .AND.
878* $ (MOD(K-1,HBL) .GT. 1) ) THEN
879* H11 = A((ICOL1-3)*LDA+IROW1-2)
880* H10 = A((ICOL1-3)*LDA+IROW1-1)
881* H22 = A((ICOL1-2)*LDA+IROW1-1)
882* IF ( ABS(H10).LE.MAX(ULP*(ABS(H11)+ABS(H22)),
883* $ SMLNUM) ) THEN
884* A((ICOL1-3)*LDA+IROW1-1) = ZERO
885* END IF
886* END IF
887 END IF
888 ELSE IF( m.GT.l ) THEN
889 IF( mod( k-1, hbl ).GT.0 ) THEN
890 a( ( icol1-2 )*lda+irow1 ) = -a( ( icol1-2 )*
891 $ lda+irow1 )
892 END IF
893 END IF
894 v2 = vcopy( 2 )
895 t2 = t1copy*v2
896 work( vecsidx+( k-1 )*3+1 ) = vcopy( 2 )
897 work( vecsidx+( k-1 )*3+2 ) = vcopy( 3 )
898 work( vecsidx+( k-1 )*3+3 ) = t1copy
899 t1 = t1copy
900 IF( k.LT.istop ) THEN
901*
902* Do some work so next step is ready...
903*
904 v3 = vcopy( 3 )
905 t3 = t1*v3
906 DO 50 j = icol1, min( k2( ki )+1, i-1 ) +
907 $ icol1 - k
908 sum = a( ( j-1 )*lda+irow1 ) +
909 $ v2*a( ( j-1 )*lda+irow1+1 ) +
910 $ v3*a( ( j-1 )*lda+irow1+2 )
911 a( ( j-1 )*lda+irow1 ) = a( ( j-1 )*lda+
912 $ irow1 ) - sum*t1
913 a( ( j-1 )*lda+irow1+1 ) = a( ( j-1 )*lda+
914 $ irow1+1 ) - sum*t2
915 a( ( j-1 )*lda+irow1+2 ) = a( ( j-1 )*lda+
916 $ irow1+2 ) - sum*t3
917 50 CONTINUE
918 itmp1 = localk2( ki )
919 DO 60 j = irow1 + 1, irow1 + 3
920 sum = a( ( icol1-1 )*lda+j ) +
921 $ v2*a( icol1*lda+j ) +
922 $ v3*a( ( icol1+1 )*lda+j )
923 a( ( icol1-1 )*lda+j ) = a( ( icol1-1 )*lda+
924 $ j ) - sum*t1
925 a( icol1*lda+j ) = a( icol1*lda+j ) - sum*t2
926 a( ( icol1+1 )*lda+j ) = a( ( icol1+1 )*lda+
927 $ j ) - sum*t3
928 60 CONTINUE
929 END IF
930 irow1 = irow1 + 1
931 icol1 = icol1 + 1
932 70 CONTINUE
933 END IF
934*
935 IF( modkm1.EQ.hbl-2 ) THEN
936 IF( ( down.EQ.icurrow( ki ) ) .AND.
937 $ ( right.EQ.icurcol( ki ) ) .AND. ( num.GT.1 ) )
938 $ THEN
939 CALL sgerv2d( contxt, 3, 1,
940 $ work( vecsidx+( istart-1 )*3+1 ), 3,
941 $ down, right )
942 END IF
943 IF( ( myrow.EQ.icurrow( ki ) ) .AND.
944 $ ( mycol.EQ.icurcol( ki ) ) .AND. ( num.GT.1 ) )
945 $ THEN
946 CALL sgesd2d( contxt, 3, 1,
947 $ work( vecsidx+( istart-1 )*3+1 ), 3,
948 $ up, left )
949 END IF
950 IF( ( down.EQ.icurrow( ki ) ) .AND.
951 $ ( npcol.GT.1 ) .AND. ( istart.LE.istop ) ) THEN
952 jj = mod( icurcol( ki )+npcol-1, npcol )
953 IF( mycol.NE.jj ) THEN
954 CALL sgebr2d( contxt, 'ROW', ' ',
955 $ 3*( istop-istart+1 ), 1,
956 $ work( vecsidx+( istart-1 )*3+1 ),
957 $ 3*( istop-istart+1 ), myrow, jj )
958 ELSE
959 CALL sgebs2d( contxt, 'ROW', ' ',
960 $ 3*( istop-istart+1 ), 1,
961 $ work( vecsidx+( istart-1 )*3+1 ),
962 $ 3*( istop-istart+1 ) )
963 END IF
964 END IF
965 END IF
966*
967* Broadcast Householder information from the block
968*
969 IF( ( myrow.EQ.icurrow( ki ) ) .AND. ( npcol.GT.1 ) .AND.
970 $ ( istart.LE.istop ) ) THEN
971 IF( mycol.NE.icurcol( ki ) ) THEN
972 CALL sgebr2d( contxt, 'ROW', ' ',
973 $ 3*( istop-istart+1 ), 1,
974 $ work( vecsidx+( istart-1 )*3+1 ),
975 $ 3*( istop-istart+1 ), myrow,
976 $ icurcol( ki ) )
977 ELSE
978 CALL sgebs2d( contxt, 'ROW', ' ',
979 $ 3*( istop-istart+1 ), 1,
980 $ work( vecsidx+( istart-1 )*3+1 ),
981 $ 3*( istop-istart+1 ) )
982 END IF
983 END IF
984 80 CONTINUE
985*
986* Now do column transforms and finish work
987*
988 DO 90 ki = 1, ibulge
989*
990 istart = max( k1( ki ), m )
991 istop = min( k2( ki ), i-1 )
992*
993 IF( mod( istart-1, hbl ).EQ.hbl-2 ) THEN
994 IF( ( right.EQ.icurcol( ki ) ) .AND.
995 $ ( nprow.GT.1 ) .AND. ( istart.LE.istop ) ) THEN
996 jj = mod( icurrow( ki )+nprow-1, nprow )
997 IF( myrow.NE.jj ) THEN
998 CALL sgebr2d( contxt, 'COL', ' ',
999 $ 3*( istop-istart+1 ), 1,
1000 $ work( vecsidx+( istart-1 )*3+1 ),
1001 $ 3*( istop-istart+1 ), jj, mycol )
1002 ELSE
1003 CALL sgebs2d( contxt, 'COL', ' ',
1004 $ 3*( istop-istart+1 ), 1,
1005 $ work( vecsidx+( istart-1 )*3+1 ),
1006 $ 3*( istop-istart+1 ) )
1007 END IF
1008 END IF
1009 END IF
1010*
1011 IF( ( mycol.EQ.icurcol( ki ) ) .AND. ( nprow.GT.1 ) .AND.
1012 $ ( istart.LE.istop ) ) THEN
1013 IF( myrow.NE.icurrow( ki ) ) THEN
1014 CALL sgebr2d( contxt, 'COL', ' ',
1015 $ 3*( istop-istart+1 ), 1,
1016 $ work( vecsidx+( istart-1 )*3+1 ),
1017 $ 3*( istop-istart+1 ), icurrow( ki ),
1018 $ mycol )
1019 ELSE
1020 CALL sgebs2d( contxt, 'COL', ' ',
1021 $ 3*( istop-istart+1 ), 1,
1022 $ work( vecsidx+( istart-1 )*3+1 ),
1023 $ 3*( istop-istart+1 ) )
1024 END IF
1025 END IF
1026 90 CONTINUE
1027*
1028* Now do make up work to have things in block fashion
1029*
1030 DO 150 ki = 1, ibulge
1031 istart = max( k1( ki ), m )
1032 istop = min( k2( ki ), i-1 )
1033*
1034 modkm1 = mod( istart-1, hbl )
1035 IF( ( myrow.EQ.icurrow( ki ) ) .AND.
1036 $ ( mycol.EQ.icurcol( ki ) ) .AND.
1037 $ ( modkm1.EQ.hbl-2 ) .AND. ( istart.LT.i-1 ) ) THEN
1038 k = istart
1039*
1040* Catch up on column & border work
1041*
1042 nr = min( 3, i-k+1 )
1043 v2 = work( vecsidx+( k-1 )*3+1 )
1044 v3 = work( vecsidx+( k-1 )*3+2 )
1045 t1 = work( vecsidx+( k-1 )*3+3 )
1046 IF( nr.EQ.3 ) THEN
1047*
1048* Do some work so next step is ready...
1049*
1050* V3 = VCOPY( 3 )
1051 t2 = t1*v2
1052 t3 = t1*v3
1053 itmp1 = min( 6, i2+2-k )
1054 itmp2 = max( i1-k+2, 1 )
1055 DO 100 j = 2, itmp1
1056 sum = smalla( 2, j, ki ) +
1057 $ v2*smalla( 3, j, ki ) +
1058 $ v3*smalla( 4, j, ki )
1059 smalla( 2, j, ki ) = smalla( 2, j, ki ) - sum*t1
1060 smalla( 3, j, ki ) = smalla( 3, j, ki ) - sum*t2
1061 smalla( 4, j, ki ) = smalla( 4, j, ki ) - sum*t3
1062 100 CONTINUE
1063 DO 110 j = itmp2, 5
1064 sum = smalla( j, 2, ki ) +
1065 $ v2*smalla( j, 3, ki ) +
1066 $ v3*smalla( j, 4, ki )
1067 smalla( j, 2, ki ) = smalla( j, 2, ki ) - sum*t1
1068 smalla( j, 3, ki ) = smalla( j, 3, ki ) - sum*t2
1069 smalla( j, 4, ki ) = smalla( j, 4, ki ) - sum*t3
1070 110 CONTINUE
1071 END IF
1072 END IF
1073*
1074 IF( ( mod( istart-1, hbl ).EQ.hbl-1 ) .AND.
1075 $ ( istart.LE.istop ) .AND.
1076 $ ( myrow.EQ.icurrow( ki ) ) .AND.
1077 $ ( mycol.EQ.icurcol( ki ) ) ) THEN
1078 k = istop
1079*
1080* Catch up on column & border work
1081*
1082 nr = min( 3, i-k+1 )
1083 v2 = work( vecsidx+( k-1 )*3+1 )
1084 v3 = work( vecsidx+( k-1 )*3+2 )
1085 t1 = work( vecsidx+( k-1 )*3+3 )
1086 IF( nr.EQ.3 ) THEN
1087*
1088* Do some work so next step is ready...
1089*
1090* V3 = VCOPY( 3 )
1091 t2 = t1*v2
1092 t3 = t1*v3
1093 itmp1 = min( 6, i2-k+3 )
1094 itmp2 = max( i1-k+3, 1 )
1095 DO 120 j = 3, itmp1
1096 sum = smalla( 3, j, ki ) +
1097 $ v2*smalla( 4, j, ki ) +
1098 $ v3*smalla( 5, j, ki )
1099 smalla( 3, j, ki ) = smalla( 3, j, ki ) - sum*t1
1100 smalla( 4, j, ki ) = smalla( 4, j, ki ) - sum*t2
1101 smalla( 5, j, ki ) = smalla( 5, j, ki ) - sum*t3
1102 120 CONTINUE
1103 DO 130 j = itmp2, 6
1104 sum = smalla( j, 3, ki ) +
1105 $ v2*smalla( j, 4, ki ) +
1106 $ v3*smalla( j, 5, ki )
1107 smalla( j, 3, ki ) = smalla( j, 3, ki ) - sum*t1
1108 smalla( j, 4, ki ) = smalla( j, 4, ki ) - sum*t2
1109 smalla( j, 5, ki ) = smalla( j, 5, ki ) - sum*t3
1110 130 CONTINUE
1111 END IF
1112 END IF
1113*
1114 modkm1 = mod( istart-1, hbl )
1115 IF( ( myrow.EQ.icurrow( ki ) ) .AND.
1116 $ ( mycol.EQ.icurcol( ki ) ) .AND.
1117 $ ( ( ( modkm1.EQ.hbl-2 ) .AND. ( istart.EQ.i-
1118 $ 1 ) ) .OR. ( ( modkm1.LT.hbl-2 ) .AND. ( istart.LE.i-
1119 $ 1 ) ) ) ) THEN
1120*
1121* (IROW1,ICOL1) is (I,J)-coordinates of H(ISTART,ISTART)
1122*
1123 irow1 = krow( ki )
1124 icol1 = localk2( ki )
1125 DO 140 k = istart, istop
1126*
1127* Catch up on column & border work
1128*
1129 nr = min( 3, i-k+1 )
1130 v2 = work( vecsidx+( k-1 )*3+1 )
1131 v3 = work( vecsidx+( k-1 )*3+2 )
1132 t1 = work( vecsidx+( k-1 )*3+3 )
1133 IF( k.LT.istop ) THEN
1134*
1135* Do some work so next step is ready...
1136*
1137 t2 = t1*v2
1138 t3 = t1*v3
1139 CALL slaref( 'Col', a, lda, .false., z, ldz,
1140 $ .false., icol1, icol1, istart,
1141 $ istop, min( istart+1, i )-k+irow1,
1142 $ irow1, liloz, lihiz,
1143 $ work( vecsidx+1 ), v2, v3, t1, t2,
1144 $ t3 )
1145 irow1 = irow1 + 1
1146 icol1 = icol1 + 1
1147 ELSE
1148 IF( ( nr.EQ.3 ) .AND. ( mod( k-1,
1149 $ hbl ).LT.hbl-2 ) ) THEN
1150 t2 = t1*v2
1151 t3 = t1*v3
1152 CALL slaref( 'Row', a, lda, .false., z, ldz,
1153 $ .false., irow1, irow1, istart,
1154 $ istop, icol1, min( min( k2( ki )
1155 $ +1, i-1 ), i2 )-k+icol1, liloz,
1156 $ lihiz, work( vecsidx+1 ), v2,
1157 $ v3, t1, t2, t3 )
1158 END IF
1159 END IF
1160 140 CONTINUE
1161 END IF
1162*
1163* Send SMALLA back again.
1164*
1165 k = istart
1166 modkm1 = mod( k-1, hbl )
1167 IF( ( modkm1.GE.hbl-2 ) .AND. ( k.LE.i-1 ) ) THEN
1168 IF( ( modkm1.EQ.hbl-2 ) .AND. ( k.LT.i-1 ) ) THEN
1169*
1170* Copy 6 elements from global A(K-1:K+4,K-1:K+4)
1171*
1172 CALL infog2l( k+2, k+2, desca, nprow, npcol, myrow,
1173 $ mycol, irow1, icol1, itmp1, itmp2 )
1174 CALL pslacp3( min( 6, n-k+2 ), k-1, a, desca,
1175 $ smalla( 1, 1, ki ), 6, itmp1, itmp2,
1176 $ 1 )
1177*
1178 END IF
1179 IF( modkm1.EQ.hbl-1 ) THEN
1180*
1181* Copy 6 elements from global A(K-2:K+3,K-2:K+3)
1182*
1183 CALL infog2l( k+1, k+1, desca, nprow, npcol, myrow,
1184 $ mycol, irow1, icol1, itmp1, itmp2 )
1185 CALL pslacp3( min( 6, n-k+3 ), k-2, a, desca,
1186 $ smalla( 1, 1, ki ), 6, itmp1, itmp2,
1187 $ 1 )
1188 END IF
1189 END IF
1190*
1191 150 CONTINUE
1192*
1193* Now start major set of block ROW reflections
1194*
1195 DO 160 ki = 1, ibulge
1196 IF( ( myrow.NE.icurrow( ki ) ) .AND.
1197 $ ( down.NE.icurrow( ki ) ) )GO TO 160
1198 istart = max( k1( ki ), m )
1199 istop = min( k2( ki ), i-1 )
1200*
1201 IF( ( istop.GT.istart ) .AND.
1202 $ ( mod( istart-1, hbl ).LT.hbl-2 ) .AND.
1203 $ ( icurrow( ki ).EQ.myrow ) ) THEN
1204 irow1 = min( k2( ki )+1, i-1 ) + 1
1205 CALL infog1l( irow1, hbl, npcol, mycol, 0, itmp1,
1206 $ itmp2 )
1207 itmp2 = numroc( i2, hbl, mycol, 0, npcol )
1208 ii = krow( ki )
1209 CALL slaref( 'Row', a, lda, wantz, z, ldz, .true., ii,
1210 $ ii, istart, istop, itmp1, itmp2, liloz,
1211 $ lihiz, work( vecsidx+1 ), v2, v3, t1, t2,
1212 $ t3 )
1213 END IF
1214 160 CONTINUE
1215*
1216 DO 180 ki = 1, ibulge
1217 IF( krow( ki ).GT.kp2row( ki ) )
1218 $ GO TO 180
1219 IF( ( myrow.NE.icurrow( ki ) ) .AND.
1220 $ ( down.NE.icurrow( ki ) ) )GO TO 180
1221 istart = max( k1( ki ), m )
1222 istop = min( k2( ki ), i-1 )
1223 IF( ( istart.EQ.istop ) .OR.
1224 $ ( mod( istart-1, hbl ).GE.hbl-2 ) .OR.
1225 $ ( icurrow( ki ).NE.myrow ) ) THEN
1226 DO 170 k = istart, istop
1227 v2 = work( vecsidx+( k-1 )*3+1 )
1228 v3 = work( vecsidx+( k-1 )*3+2 )
1229 t1 = work( vecsidx+( k-1 )*3+3 )
1230 nr = min( 3, i-k+1 )
1231 IF( ( nr.EQ.3 ) .AND. ( krow( ki ).LE.
1232 $ kp2row( ki ) ) ) THEN
1233 IF( ( k.LT.istop ) .AND.
1234 $ ( mod( k-1, hbl ).LT.hbl-2 ) ) THEN
1235 itmp1 = min( k2( ki )+1, i-1 ) + 1
1236 ELSE
1237 IF( mod( k-1, hbl ).LT.hbl-2 ) THEN
1238 itmp1 = min( k2( ki )+1, i-1 ) + 1
1239 END IF
1240 IF( mod( k-1, hbl ).EQ.hbl-2 ) THEN
1241 itmp1 = min( k+4, i2 ) + 1
1242 END IF
1243 IF( mod( k-1, hbl ).EQ.hbl-1 ) THEN
1244 itmp1 = min( k+3, i2 ) + 1
1245 END IF
1246 END IF
1247*
1248* Find local coor of rows K through K+2
1249*
1250 irow1 = krow( ki )
1251 irow2 = kp2row( ki )
1252 CALL infog1l( itmp1, hbl, npcol, mycol, 0,
1253 $ icol1, icol2 )
1254 icol2 = numroc( i2, hbl, mycol, 0, npcol )
1255 IF( ( mod( k-1, hbl ).LT.hbl-2 ) .OR.
1256 $ ( nprow.EQ.1 ) ) THEN
1257 t2 = t1*v2
1258 t3 = t1*v3
1259 CALL slaref( 'Row', a, lda, wantz, z, ldz,
1260 $ .false., irow1, irow1, istart,
1261 $ istop, icol1, icol2, liloz,
1262 $ lihiz, work( vecsidx+1 ), v2,
1263 $ v3, t1, t2, t3 )
1264 END IF
1265 IF( ( mod( k-1, hbl ).EQ.hbl-2 ) .AND.
1266 $ ( nprow.GT.1 ) ) THEN
1267 IF( irow1.EQ.irow2 ) THEN
1268 CALL sgesd2d( contxt, 1, icol2-icol1+1,
1269 $ a( ( icol1-1 )*lda+irow2 ),
1270 $ lda, up, mycol )
1271 END IF
1272 END IF
1273 IF( ( mod( k-1, hbl ).EQ.hbl-1 ) .AND.
1274 $ ( nprow.GT.1 ) ) THEN
1275 IF( irow1.EQ.irow2 ) THEN
1276 CALL sgesd2d( contxt, 1, icol2-icol1+1,
1277 $ a( ( icol1-1 )*lda+irow1 ),
1278 $ lda, down, mycol )
1279 END IF
1280 END IF
1281 END IF
1282 170 CONTINUE
1283 END IF
1284 180 CONTINUE
1285*
1286 DO 220 ki = 1, ibulge
1287 IF( krow( ki ).GT.kp2row( ki ) )
1288 $ GO TO 220
1289 IF( ( myrow.NE.icurrow( ki ) ) .AND.
1290 $ ( down.NE.icurrow( ki ) ) )GO TO 220
1291 istart = max( k1( ki ), m )
1292 istop = min( k2( ki ), i-1 )
1293 IF( ( istart.EQ.istop ) .OR.
1294 $ ( mod( istart-1, hbl ).GE.hbl-2 ) .OR.
1295 $ ( icurrow( ki ).NE.myrow ) ) THEN
1296 DO 210 k = istart, istop
1297 v2 = work( vecsidx+( k-1 )*3+1 )
1298 v3 = work( vecsidx+( k-1 )*3+2 )
1299 t1 = work( vecsidx+( k-1 )*3+3 )
1300 nr = min( 3, i-k+1 )
1301 IF( ( nr.EQ.3 ) .AND. ( krow( ki ).LE.
1302 $ kp2row( ki ) ) ) THEN
1303 IF( ( k.LT.istop ) .AND.
1304 $ ( mod( k-1, hbl ).LT.hbl-2 ) ) THEN
1305 itmp1 = min( k2( ki )+1, i-1 ) + 1
1306 ELSE
1307 IF( mod( k-1, hbl ).LT.hbl-2 ) THEN
1308 itmp1 = min( k2( ki )+1, i-1 ) + 1
1309 END IF
1310 IF( mod( k-1, hbl ).EQ.hbl-2 ) THEN
1311 itmp1 = min( k+4, i2 ) + 1
1312 END IF
1313 IF( mod( k-1, hbl ).EQ.hbl-1 ) THEN
1314 itmp1 = min( k+3, i2 ) + 1
1315 END IF
1316 END IF
1317*
1318 irow1 = krow( ki ) + k - istart
1319 irow2 = kp2row( ki ) + k - istart
1320 CALL infog1l( itmp1, hbl, npcol, mycol, 0,
1321 $ icol1, icol2 )
1322 icol2 = numroc( i2, hbl, mycol, 0, npcol )
1323 IF( ( mod( k-1, hbl ).EQ.hbl-2 ) .AND.
1324 $ ( nprow.GT.1 ) ) THEN
1325 IF( irow1.NE.irow2 ) THEN
1326 CALL sgerv2d( contxt, 1, icol2-icol1+1,
1327 $ work( irbuf+1 ), 1, down,
1328 $ mycol )
1329 t2 = t1*v2
1330 t3 = t1*v3
1331 DO 190 j = icol1, icol2
1332 sum = a( ( j-1 )*lda+irow1 ) +
1333 $ v2*a( ( j-1 )*lda+irow1+1 ) +
1334 $ v3*work( irbuf+j-icol1+1 )
1335 a( ( j-1 )*lda+irow1 ) = a( ( j-1 )*
1336 $ lda+irow1 ) - sum*t1
1337 a( ( j-1 )*lda+irow1+1 ) = a( ( j-1 )*
1338 $ lda+irow1+1 ) - sum*t2
1339 work( irbuf+j-icol1+1 ) = work( irbuf+
1340 $ j-icol1+1 ) - sum*t3
1341 190 CONTINUE
1342 CALL sgesd2d( contxt, 1, icol2-icol1+1,
1343 $ work( irbuf+1 ), 1, down,
1344 $ mycol )
1345 END IF
1346 END IF
1347 IF( ( mod( k-1, hbl ).EQ.hbl-1 ) .AND.
1348 $ ( nprow.GT.1 ) ) THEN
1349 IF( irow1.NE.irow2 ) THEN
1350 CALL sgerv2d( contxt, 1, icol2-icol1+1,
1351 $ work( irbuf+1 ), 1, up,
1352 $ mycol )
1353 t2 = t1*v2
1354 t3 = t1*v3
1355 DO 200 j = icol1, icol2
1356 sum = work( irbuf+j-icol1+1 ) +
1357 $ v2*a( ( j-1 )*lda+irow1 ) +
1358 $ v3*a( ( j-1 )*lda+irow1+1 )
1359 work( irbuf+j-icol1+1 ) = work( irbuf+
1360 $ j-icol1+1 ) - sum*t1
1361 a( ( j-1 )*lda+irow1 ) = a( ( j-1 )*
1362 $ lda+irow1 ) - sum*t2
1363 a( ( j-1 )*lda+irow1+1 ) = a( ( j-1 )*
1364 $ lda+irow1+1 ) - sum*t3
1365 200 CONTINUE
1366 CALL sgesd2d( contxt, 1, icol2-icol1+1,
1367 $ work( irbuf+1 ), 1, up,
1368 $ mycol )
1369 END IF
1370 END IF
1371 END IF
1372 210 CONTINUE
1373 END IF
1374 220 CONTINUE
1375*
1376 DO 240 ki = 1, ibulge
1377 IF( krow( ki ).GT.kp2row( ki ) )
1378 $ GO TO 240
1379 IF( ( myrow.NE.icurrow( ki ) ) .AND.
1380 $ ( down.NE.icurrow( ki ) ) )GO TO 240
1381 istart = max( k1( ki ), m )
1382 istop = min( k2( ki ), i-1 )
1383 IF( ( istart.EQ.istop ) .OR.
1384 $ ( mod( istart-1, hbl ).GE.hbl-2 ) .OR.
1385 $ ( icurrow( ki ).NE.myrow ) ) THEN
1386 DO 230 k = istart, istop
1387 v2 = work( vecsidx+( k-1 )*3+1 )
1388 v3 = work( vecsidx+( k-1 )*3+2 )
1389 t1 = work( vecsidx+( k-1 )*3+3 )
1390 nr = min( 3, i-k+1 )
1391 IF( ( nr.EQ.3 ) .AND. ( krow( ki ).LE.
1392 $ kp2row( ki ) ) ) THEN
1393 IF( ( k.LT.istop ) .AND.
1394 $ ( mod( k-1, hbl ).LT.hbl-2 ) ) THEN
1395 itmp1 = min( k2( ki )+1, i-1 ) + 1
1396 ELSE
1397 IF( mod( k-1, hbl ).LT.hbl-2 ) THEN
1398 itmp1 = min( k2( ki )+1, i-1 ) + 1
1399 END IF
1400 IF( mod( k-1, hbl ).EQ.hbl-2 ) THEN
1401 itmp1 = min( k+4, i2 ) + 1
1402 END IF
1403 IF( mod( k-1, hbl ).EQ.hbl-1 ) THEN
1404 itmp1 = min( k+3, i2 ) + 1
1405 END IF
1406 END IF
1407*
1408 irow1 = krow( ki ) + k - istart
1409 irow2 = kp2row( ki ) + k - istart
1410 CALL infog1l( itmp1, hbl, npcol, mycol, 0,
1411 $ icol1, icol2 )
1412 icol2 = numroc( i2, hbl, mycol, 0, npcol )
1413 IF( ( mod( k-1, hbl ).EQ.hbl-2 ) .AND.
1414 $ ( nprow.GT.1 ) ) THEN
1415 IF( irow1.EQ.irow2 ) THEN
1416 CALL sgerv2d( contxt, 1, icol2-icol1+1,
1417 $ a( ( icol1-1 )*lda+irow2 ),
1418 $ lda, up, mycol )
1419 END IF
1420 END IF
1421 IF( ( mod( k-1, hbl ).EQ.hbl-1 ) .AND.
1422 $ ( nprow.GT.1 ) ) THEN
1423 IF( irow1.EQ.irow2 ) THEN
1424 CALL sgerv2d( contxt, 1, icol2-icol1+1,
1425 $ a( ( icol1-1 )*lda+irow1 ),
1426 $ lda, down, mycol )
1427 END IF
1428 END IF
1429 END IF
1430 230 CONTINUE
1431 END IF
1432 240 CONTINUE
1433 250 CONTINUE
1434*
1435* Now start major set of block COL reflections
1436*
1437 DO 260 ki = 1, ibulge
1438 IF( ( mycol.NE.icurcol( ki ) ) .AND.
1439 $ ( right.NE.icurcol( ki ) ) )GO TO 260
1440 istart = max( k1( ki ), m )
1441 istop = min( k2( ki ), i-1 )
1442*
1443 IF( ( ( mod( istart-1, hbl ).LT.hbl-2 ) .OR. ( npcol.EQ.
1444 $ 1 ) ) .AND. ( icurcol( ki ).EQ.mycol ) .AND.
1445 $ ( i-istop+1.GE.3 ) ) THEN
1446 k = istart
1447 IF( ( k.LT.istop ) .AND. ( mod( k-1,
1448 $ hbl ).LT.hbl-2 ) ) THEN
1449 itmp1 = min( istart+1, i ) - 1
1450 ELSE
1451 IF( mod( k-1, hbl ).LT.hbl-2 ) THEN
1452 itmp1 = min( k+3, i )
1453 END IF
1454 IF( mod( k-1, hbl ).EQ.hbl-2 ) THEN
1455 itmp1 = max( i1, k-1 ) - 1
1456 END IF
1457 IF( mod( k-1, hbl ).EQ.hbl-1 ) THEN
1458 itmp1 = max( i1, k-2 ) - 1
1459 END IF
1460 END IF
1461*
1462 icol1 = kcol( ki )
1463 CALL infog1l( i1, hbl, nprow, myrow, 0, irow1, irow2 )
1464 irow2 = numroc( itmp1, hbl, myrow, 0, nprow )
1465 IF( irow1.LE.irow2 ) THEN
1466 itmp2 = irow2
1467 ELSE
1468 itmp2 = -1
1469 END IF
1470 CALL slaref( 'Col', a, lda, wantz, z, ldz, .true.,
1471 $ icol1, icol1, istart, istop, irow1,
1472 $ irow2, liloz, lihiz, work( vecsidx+1 ),
1473 $ v2, v3, t1, t2, t3 )
1474 k = istop
1475 IF( mod( k-1, hbl ).LT.hbl-2 ) THEN
1476*
1477* Do from ITMP1+1 to MIN(K+3,I)
1478*
1479 IF( mod( k-1, hbl ).LT.hbl-3 ) THEN
1480 irow1 = itmp2 + 1
1481 IF( mod( ( itmp1 / hbl ), nprow ).EQ.myrow )
1482 $ THEN
1483 IF( itmp2.GT.0 ) THEN
1484 irow2 = itmp2 + min( k+3, i ) - itmp1
1485 ELSE
1486 irow2 = irow1 - 1
1487 END IF
1488 ELSE
1489 irow2 = irow1 - 1
1490 END IF
1491 ELSE
1492 CALL infog1l( itmp1+1, hbl, nprow, myrow, 0,
1493 $ irow1, irow2 )
1494 irow2 = numroc( min( k+3, i ), hbl, myrow, 0,
1495 $ nprow )
1496 END IF
1497 v2 = work( vecsidx+( k-1 )*3+1 )
1498 v3 = work( vecsidx+( k-1 )*3+2 )
1499 t1 = work( vecsidx+( k-1 )*3+3 )
1500 t2 = t1*v2
1501 t3 = t1*v3
1502 icol1 = kcol( ki ) + istop - istart
1503 CALL slaref( 'Col', a, lda, .false., z, ldz,
1504 $ .false., icol1, icol1, istart, istop,
1505 $ irow1, irow2, liloz, lihiz,
1506 $ work( vecsidx+1 ), v2, v3, t1, t2,
1507 $ t3 )
1508 END IF
1509 END IF
1510 260 CONTINUE
1511*
1512 DO 320 ki = 1, ibulge
1513 IF( kcol( ki ).GT.kp2col( ki ) )
1514 $ GO TO 320
1515 IF( ( mycol.NE.icurcol( ki ) ) .AND.
1516 $ ( right.NE.icurcol( ki ) ) )GO TO 320
1517 istart = max( k1( ki ), m )
1518 istop = min( k2( ki ), i-1 )
1519 IF( mod( istart-1, hbl ).GE.hbl-2 ) THEN
1520*
1521* INFO is found in a buffer
1522*
1523 ispec = 1
1524 ELSE
1525*
1526* All INFO is local
1527*
1528 ispec = 0
1529 END IF
1530*
1531 DO 310 k = istart, istop
1532*
1533 v2 = work( vecsidx+( k-1 )*3+1 )
1534 v3 = work( vecsidx+( k-1 )*3+2 )
1535 t1 = work( vecsidx+( k-1 )*3+3 )
1536 nr = min( 3, i-k+1 )
1537 IF( ( nr.EQ.3 ) .AND. ( kcol( ki ).LE.kp2col( ki ) ) )
1538 $ THEN
1539*
1540 IF( ( k.LT.istop ) .AND.
1541 $ ( mod( k-1, hbl ).LT.hbl-2 ) ) THEN
1542 itmp1 = min( istart+1, i ) - 1
1543 ELSE
1544 IF( mod( k-1, hbl ).LT.hbl-2 ) THEN
1545 itmp1 = min( k+3, i )
1546 END IF
1547 IF( mod( k-1, hbl ).EQ.hbl-2 ) THEN
1548 itmp1 = max( i1, k-1 ) - 1
1549 END IF
1550 IF( mod( k-1, hbl ).EQ.hbl-1 ) THEN
1551 itmp1 = max( i1, k-2 ) - 1
1552 END IF
1553 END IF
1554 icol1 = kcol( ki ) + k - istart
1555 icol2 = kp2col( ki ) + k - istart
1556 CALL infog1l( i1, hbl, nprow, myrow, 0, irow1,
1557 $ irow2 )
1558 irow2 = numroc( itmp1, hbl, myrow, 0, nprow )
1559 IF( ( mod( k-1, hbl ).EQ.hbl-2 ) .AND.
1560 $ ( npcol.GT.1 ) ) THEN
1561 IF( icol1.EQ.icol2 ) THEN
1562 CALL sgesd2d( contxt, irow2-irow1+1, 1,
1563 $ a( ( icol1-1 )*lda+irow1 ),
1564 $ lda, myrow, left )
1565 CALL sgerv2d( contxt, irow2-irow1+1, 1,
1566 $ a( ( icol1-1 )*lda+irow1 ),
1567 $ lda, myrow, left )
1568 ELSE
1569 CALL sgerv2d( contxt, irow2-irow1+1, 1,
1570 $ work( icbuf+1 ), irow2-irow1+1,
1571 $ myrow, right )
1572 t2 = t1*v2
1573 t3 = t1*v3
1574 DO 270 j = irow1, irow2
1575 sum = a( ( icol1-1 )*lda+j ) +
1576 $ v2*a( icol1*lda+j ) +
1577 $ v3*work( icbuf+j-irow1+1 )
1578 a( ( icol1-1 )*lda+j ) = a( ( icol1-1 )*
1579 $ lda+j ) - sum*t1
1580 a( icol1*lda+j ) = a( icol1*lda+j ) -
1581 $ sum*t2
1582 work( icbuf+j-irow1+1 ) = work( icbuf+j-
1583 $ irow1+1 ) - sum*t3
1584 270 CONTINUE
1585 CALL sgesd2d( contxt, irow2-irow1+1, 1,
1586 $ work( icbuf+1 ), irow2-irow1+1,
1587 $ myrow, right )
1588 END IF
1589 END IF
1590 IF( ( mod( k-1, hbl ).EQ.hbl-1 ) .AND.
1591 $ ( npcol.GT.1 ) ) THEN
1592 IF( icol1.EQ.icol2 ) THEN
1593 CALL sgesd2d( contxt, irow2-irow1+1, 1,
1594 $ a( ( icol1-1 )*lda+irow1 ),
1595 $ lda, myrow, right )
1596 CALL sgerv2d( contxt, irow2-irow1+1, 1,
1597 $ a( ( icol1-1 )*lda+irow1 ),
1598 $ lda, myrow, right )
1599 ELSE
1600 CALL sgerv2d( contxt, irow2-irow1+1, 1,
1601 $ work( icbuf+1 ), irow2-irow1+1,
1602 $ myrow, left )
1603 t2 = t1*v2
1604 t3 = t1*v3
1605 DO 280 j = irow1, irow2
1606 sum = work( icbuf+j-irow1+1 ) +
1607 $ v2*a( ( icol1-1 )*lda+j ) +
1608 $ v3*a( icol1*lda+j )
1609 work( icbuf+j-irow1+1 ) = work( icbuf+j-
1610 $ irow1+1 ) - sum*t1
1611 a( ( icol1-1 )*lda+j ) = a( ( icol1-1 )*
1612 $ lda+j ) - sum*t2
1613 a( icol1*lda+j ) = a( icol1*lda+j ) -
1614 $ sum*t3
1615 280 CONTINUE
1616 CALL sgesd2d( contxt, irow2-irow1+1, 1,
1617 $ work( icbuf+1 ), irow2-irow1+1,
1618 $ myrow, left )
1619 END IF
1620 END IF
1621*
1622* If we want Z and we haven't already done any Z
1623 IF( ( wantz ) .AND. ( mod( k-1,
1624 $ hbl ).GE.hbl-2 ) .AND. ( npcol.GT.1 ) ) THEN
1625*
1626* Accumulate transformations in the matrix Z
1627*
1628 irow1 = liloz
1629 irow2 = lihiz
1630 IF( mod( k-1, hbl ).EQ.hbl-2 ) THEN
1631 IF( icol1.EQ.icol2 ) THEN
1632 CALL sgesd2d( contxt, irow2-irow1+1, 1,
1633 $ z( ( icol1-1 )*ldz+irow1 ),
1634 $ ldz, myrow, left )
1635 CALL sgerv2d( contxt, irow2-irow1+1, 1,
1636 $ z( ( icol1-1 )*ldz+irow1 ),
1637 $ ldz, myrow, left )
1638 ELSE
1639 CALL sgerv2d( contxt, irow2-irow1+1, 1,
1640 $ work( icbuf+1 ),
1641 $ irow2-irow1+1, myrow,
1642 $ right )
1643 t2 = t1*v2
1644 t3 = t1*v3
1645 icol1 = ( icol1-1 )*ldz
1646 DO 290 j = irow1, irow2
1647 sum = z( icol1+j ) +
1648 $ v2*z( icol1+j+ldz ) +
1649 $ v3*work( icbuf+j-irow1+1 )
1650 z( j+icol1 ) = z( j+icol1 ) - sum*t1
1651 z( j+icol1+ldz ) = z( j+icol1+ldz ) -
1652 $ sum*t2
1653 work( icbuf+j-irow1+1 ) = work( icbuf+
1654 $ j-irow1+1 ) - sum*t3
1655 290 CONTINUE
1656 CALL sgesd2d( contxt, irow2-irow1+1, 1,
1657 $ work( icbuf+1 ),
1658 $ irow2-irow1+1, myrow,
1659 $ right )
1660 END IF
1661 END IF
1662 IF( mod( k-1, hbl ).EQ.hbl-1 ) THEN
1663 IF( icol1.EQ.icol2 ) THEN
1664 CALL sgesd2d( contxt, irow2-irow1+1, 1,
1665 $ z( ( icol1-1 )*ldz+irow1 ),
1666 $ ldz, myrow, right )
1667 CALL sgerv2d( contxt, irow2-irow1+1, 1,
1668 $ z( ( icol1-1 )*ldz+irow1 ),
1669 $ ldz, myrow, right )
1670 ELSE
1671 CALL sgerv2d( contxt, irow2-irow1+1, 1,
1672 $ work( icbuf+1 ),
1673 $ irow2-irow1+1, myrow, left )
1674 t2 = t1*v2
1675 t3 = t1*v3
1676 icol1 = ( icol1-1 )*ldz
1677 DO 300 j = irow1, irow2
1678 sum = work( icbuf+j-irow1+1 ) +
1679 $ v2*z( j+icol1 ) +
1680 $ v3*z( j+icol1+ldz )
1681 work( icbuf+j-irow1+1 ) = work( icbuf+
1682 $ j-irow1+1 ) - sum*t1
1683 z( j+icol1 ) = z( j+icol1 ) - sum*t2
1684 z( j+icol1+ldz ) = z( j+icol1+ldz ) -
1685 $ sum*t3
1686 300 CONTINUE
1687 CALL sgesd2d( contxt, irow2-irow1+1, 1,
1688 $ work( icbuf+1 ),
1689 $ irow2-irow1+1, myrow, left )
1690 END IF
1691 END IF
1692 END IF
1693 IF( icurcol( ki ).EQ.mycol ) THEN
1694 IF( ( ispec.EQ.0 ) .OR. ( npcol.EQ.1 ) ) THEN
1695 localk2( ki ) = localk2( ki ) + 1
1696 END IF
1697 ELSE
1698 IF( ( mod( k-1, hbl ).EQ.hbl-1 ) .AND.
1699 $ ( icurcol( ki ).EQ.right ) ) THEN
1700 IF( k.GT.m ) THEN
1701 localk2( ki ) = localk2( ki ) + 2
1702 ELSE
1703 localk2( ki ) = localk2( ki ) + 1
1704 END IF
1705 END IF
1706 IF( ( mod( k-1, hbl ).EQ.hbl-2 ) .AND.
1707 $ ( i-k.EQ.2 ) .AND. ( icurcol( ki ).EQ.
1708 $ right ) ) THEN
1709 localk2( ki ) = localk2( ki ) + 2
1710 END IF
1711 END IF
1712 END IF
1713 310 CONTINUE
1714 320 CONTINUE
1715*
1716* Column work done
1717*
1718 330 CONTINUE
1719*
1720* Now do NR=2 work
1721*
1722 DO 410 ki = 1, ibulge
1723 istart = max( k1( ki ), m )
1724 istop = min( k2( ki ), i-1 )
1725 IF( mod( istart-1, hbl ).GE.hbl-2 ) THEN
1726*
1727* INFO is found in a buffer
1728*
1729 ispec = 1
1730 ELSE
1731*
1732* All INFO is local
1733*
1734 ispec = 0
1735 END IF
1736*
1737 DO 400 k = istart, istop
1738*
1739 v2 = work( vecsidx+( k-1 )*3+1 )
1740 v3 = work( vecsidx+( k-1 )*3+2 )
1741 t1 = work( vecsidx+( k-1 )*3+3 )
1742 nr = min( 3, i-k+1 )
1743 IF( nr.EQ.2 ) THEN
1744 IF ( icurrow( ki ).EQ.myrow ) THEN
1745 t2 = t1*v2
1746 END IF
1747 IF ( icurcol( ki ).EQ.mycol ) THEN
1748 t2 = t1*v2
1749 END IF
1750*
1751* Apply G from the left to transform the rows of the matrix
1752* in columns K to I2.
1753*
1754 CALL infog1l( k, hbl, npcol, mycol, 0, liloh,
1755 $ lihih )
1756 lihih = numroc( i2, hbl, mycol, 0, npcol )
1757 CALL infog1l( 1, hbl, nprow, myrow, 0, itmp2,
1758 $ itmp1 )
1759 itmp1 = numroc( k+1, hbl, myrow, 0, nprow )
1760 IF( icurrow( ki ).EQ.myrow ) THEN
1761 IF( ( ispec.EQ.0 ) .OR. ( nprow.EQ.1 ) .OR.
1762 $ ( mod( k-1, hbl ).EQ.hbl-2 ) ) THEN
1763 itmp1 = itmp1 - 1
1764 DO 340 j = ( liloh-1 )*lda,
1765 $ ( lihih-1 )*lda, lda
1766 sum = a( itmp1+j ) + v2*a( itmp1+1+j )
1767 a( itmp1+j ) = a( itmp1+j ) - sum*t1
1768 a( itmp1+1+j ) = a( itmp1+1+j ) - sum*t2
1769 340 CONTINUE
1770 ELSE
1771 IF( mod( k-1, hbl ).EQ.hbl-1 ) THEN
1772 CALL sgerv2d( contxt, 1, lihih-liloh+1,
1773 $ work( irbuf+1 ), 1, up,
1774 $ mycol )
1775 DO 350 j = liloh, lihih
1776 sum = work( irbuf+j-liloh+1 ) +
1777 $ v2*a( ( j-1 )*lda+itmp1 )
1778 work( irbuf+j-liloh+1 ) = work( irbuf+
1779 $ j-liloh+1 ) - sum*t1
1780 a( ( j-1 )*lda+itmp1 ) = a( ( j-1 )*
1781 $ lda+itmp1 ) - sum*t2
1782 350 CONTINUE
1783 CALL sgesd2d( contxt, 1, lihih-liloh+1,
1784 $ work( irbuf+1 ), 1, up,
1785 $ mycol )
1786 END IF
1787 END IF
1788 ELSE
1789 IF( ( mod( k-1, hbl ).EQ.hbl-1 ) .AND.
1790 $ ( icurrow( ki ).EQ.down ) ) THEN
1791 CALL sgesd2d( contxt, 1, lihih-liloh+1,
1792 $ a( ( liloh-1 )*lda+itmp1 ),
1793 $ lda, down, mycol )
1794 CALL sgerv2d( contxt, 1, lihih-liloh+1,
1795 $ a( ( liloh-1 )*lda+itmp1 ),
1796 $ lda, down, mycol )
1797 END IF
1798 END IF
1799*
1800* Apply G from the right to transform the columns of the
1801* matrix in rows I1 to MIN(K+3,I).
1802*
1803 CALL infog1l( i1, hbl, nprow, myrow, 0, liloh,
1804 $ lihih )
1805 lihih = numroc( i, hbl, myrow, 0, nprow )
1806*
1807 IF( icurcol( ki ).EQ.mycol ) THEN
1808* LOCAL A(LILOZ:LIHIZ,LOCALK2:LOCALK2+2)
1809 IF( ( ispec.EQ.0 ) .OR. ( npcol.EQ.1 ) .OR.
1810 $ ( mod( k-1, hbl ).EQ.hbl-2 ) ) THEN
1811 CALL infog1l( k, hbl, npcol, mycol, 0, itmp1,
1812 $ itmp2 )
1813 itmp2 = numroc( k+1, hbl, mycol, 0, npcol )
1814 DO 360 j = liloh, lihih
1815 sum = a( ( itmp1-1 )*lda+j ) +
1816 $ v2*a( itmp1*lda+j )
1817 a( ( itmp1-1 )*lda+j ) = a( ( itmp1-1 )*
1818 $ lda+j ) - sum*t1
1819 a( itmp1*lda+j ) = a( itmp1*lda+j ) -
1820 $ sum*t2
1821 360 CONTINUE
1822 ELSE
1823 itmp1 = localk2( ki )
1824 IF( mod( k-1, hbl ).EQ.hbl-1 ) THEN
1825 CALL sgerv2d( contxt, lihih-liloh+1, 1,
1826 $ work( icbuf+1 ),
1827 $ lihih-liloh+1, myrow, left )
1828 DO 370 j = liloh, lihih
1829 sum = work( icbuf+j ) +
1830 $ v2*a( ( itmp1-1 )*lda+j )
1831 work( icbuf+j ) = work( icbuf+j ) -
1832 $ sum*t1
1833 a( ( itmp1-1 )*lda+j )
1834 $ = a( ( itmp1-1 )*lda+j ) - sum*t2
1835 370 CONTINUE
1836 CALL sgesd2d( contxt, lihih-liloh+1, 1,
1837 $ work( icbuf+1 ),
1838 $ lihih-liloh+1, myrow, left )
1839 END IF
1840 END IF
1841 ELSE
1842 IF( ( mod( k-1, hbl ).EQ.hbl-1 ) .AND.
1843 $ ( icurcol( ki ).EQ.right ) ) THEN
1844 itmp1 = kcol( ki )
1845 CALL sgesd2d( contxt, lihih-liloh+1, 1,
1846 $ a( ( itmp1-1 )*lda+liloh ),
1847 $ lda, myrow, right )
1848 CALL infog1l( k, hbl, npcol, mycol, 0, itmp1,
1849 $ itmp2 )
1850 itmp2 = numroc( k+1, hbl, mycol, 0, npcol )
1851 CALL sgerv2d( contxt, lihih-liloh+1, 1,
1852 $ a( ( itmp1-1 )*lda+liloh ),
1853 $ lda, myrow, right )
1854 END IF
1855 END IF
1856*
1857 IF( wantz ) THEN
1858*
1859* Accumulate transformations in the matrix Z
1860*
1861 IF( icurcol( ki ).EQ.mycol ) THEN
1862* LOCAL Z(LILOZ:LIHIZ,LOCALK2:LOCALK2+2)
1863 IF( ( ispec.EQ.0 ) .OR. ( npcol.EQ.1 ) .OR.
1864 $ ( mod( k-1, hbl ).EQ.hbl-2 ) ) THEN
1865 itmp1 = kcol( ki ) + k - istart
1866 itmp1 = ( itmp1-1 )*ldz
1867 DO 380 j = liloz, lihiz
1868 sum = z( j+itmp1 ) +
1869 $ v2*z( j+itmp1+ldz )
1870 z( j+itmp1 ) = z( j+itmp1 ) - sum*t1
1871 z( j+itmp1+ldz ) = z( j+itmp1+ldz ) -
1872 $ sum*t2
1873 380 CONTINUE
1874 localk2( ki ) = localk2( ki ) + 1
1875 ELSE
1876 itmp1 = localk2( ki )
1877* IF WE ACTUALLY OWN COLUMN K
1878 IF( mod( k-1, hbl ).EQ.hbl-1 ) THEN
1879 CALL sgerv2d( contxt, lihiz-liloz+1, 1,
1880 $ work( icbuf+1 ), ldz,
1881 $ myrow, left )
1882 itmp1 = ( itmp1-1 )*ldz
1883 DO 390 j = liloz, lihiz
1884 sum = work( icbuf+j ) +
1885 $ v2*z( j+itmp1 )
1886 work( icbuf+j ) = work( icbuf+j ) -
1887 $ sum*t1
1888 z( j+itmp1 ) = z( j+itmp1 ) - sum*t2
1889 390 CONTINUE
1890 CALL sgesd2d( contxt, lihiz-liloz+1, 1,
1891 $ work( icbuf+1 ), ldz,
1892 $ myrow, left )
1893 localk2( ki ) = localk2( ki ) + 1
1894 END IF
1895 END IF
1896 ELSE
1897*
1898* NO WORK BUT NEED TO UPDATE ANYWAY????
1899*
1900 IF( ( mod( k-1, hbl ).EQ.hbl-1 ) .AND.
1901 $ ( icurcol( ki ).EQ.right ) ) THEN
1902 itmp1 = kcol( ki )
1903 itmp1 = ( itmp1-1 )*ldz
1904 CALL sgesd2d( contxt, lihiz-liloz+1, 1,
1905 $ z( liloz+itmp1 ), ldz,
1906 $ myrow, right )
1907 CALL sgerv2d( contxt, lihiz-liloz+1, 1,
1908 $ z( liloz+itmp1 ), ldz,
1909 $ myrow, right )
1910 localk2( ki ) = localk2( ki ) + 1
1911 END IF
1912 END IF
1913 END IF
1914 END IF
1915 400 CONTINUE
1916*
1917* Adjust local information for this bulge
1918*
1919 IF( nprow.EQ.1 ) THEN
1920 krow( ki ) = krow( ki ) + k2( ki ) - k1( ki ) + 1
1921 kp2row( ki ) = kp2row( ki ) + k2( ki ) - k1( ki ) + 1
1922 END IF
1923 IF( ( mod( k1( ki )-1, hbl ).LT.hbl-2 ) .AND.
1924 $ ( icurrow( ki ).EQ.myrow ) .AND. ( nprow.GT.1 ) )
1925 $ THEN
1926 krow( ki ) = krow( ki ) + k2( ki ) - k1( ki ) + 1
1927 END IF
1928 IF( ( mod( k2( ki ), hbl ).LT.hbl-2 ) .AND.
1929 $ ( icurrow( ki ).EQ.myrow ) .AND. ( nprow.GT.1 ) )
1930 $ THEN
1931 kp2row( ki ) = kp2row( ki ) + k2( ki ) - k1( ki ) + 1
1932 END IF
1933 IF( ( mod( k1( ki )-1, hbl ).GE.hbl-2 ) .AND.
1934 $ ( ( myrow.EQ.icurrow( ki ) ) .OR. ( down.EQ.
1935 $ icurrow( ki ) ) ) .AND. ( nprow.GT.1 ) ) THEN
1936 CALL infog1l( k2( ki )+1, hbl, nprow, myrow, 0,
1937 $ krow( ki ), itmp2 )
1938 itmp2 = numroc( n, hbl, myrow, 0, nprow )
1939 END IF
1940 IF( ( mod( k2( ki ), hbl ).GE.hbl-2 ) .AND.
1941 $ ( ( myrow.EQ.icurrow( ki ) ) .OR. ( up.EQ.
1942 $ icurrow( ki ) ) ) .AND. ( nprow.GT.1 ) ) THEN
1943 CALL infog1l( 1, hbl, nprow, myrow, 0, itmp2,
1944 $ kp2row( ki ) )
1945 kp2row( ki ) = numroc( k2( ki )+3, hbl, myrow, 0,
1946 $ nprow )
1947 END IF
1948 IF( npcol.EQ.1 ) THEN
1949 kcol( ki ) = kcol( ki ) + k2( ki ) - k1( ki ) + 1
1950 kp2col( ki ) = kp2col( ki ) + k2( ki ) - k1( ki ) + 1
1951 END IF
1952 IF( ( mod( k1( ki )-1, hbl ).LT.hbl-2 ) .AND.
1953 $ ( icurcol( ki ).EQ.mycol ) .AND. ( npcol.GT.1 ) )
1954 $ THEN
1955 kcol( ki ) = kcol( ki ) + k2( ki ) - k1( ki ) + 1
1956 END IF
1957 IF( ( mod( k2( ki ), hbl ).LT.hbl-2 ) .AND.
1958 $ ( icurcol( ki ).EQ.mycol ) .AND. ( npcol.GT.1 ) )
1959 $ THEN
1960 kp2col( ki ) = kp2col( ki ) + k2( ki ) - k1( ki ) + 1
1961 END IF
1962 IF( ( mod( k1( ki )-1, hbl ).GE.hbl-2 ) .AND.
1963 $ ( ( mycol.EQ.icurcol( ki ) ) .OR. ( right.EQ.
1964 $ icurcol( ki ) ) ) .AND. ( npcol.GT.1 ) ) THEN
1965 CALL infog1l( k2( ki )+1, hbl, npcol, mycol, 0,
1966 $ kcol( ki ), itmp2 )
1967 itmp2 = numroc( n, hbl, mycol, 0, npcol )
1968 END IF
1969 IF( ( mod( k2( ki ), hbl ).GE.hbl-2 ) .AND.
1970 $ ( ( mycol.EQ.icurcol( ki ) ) .OR. ( left.EQ.
1971 $ icurcol( ki ) ) ) .AND. ( npcol.GT.1 ) ) THEN
1972 CALL infog1l( 1, hbl, npcol, mycol, 0, itmp2,
1973 $ kp2col( ki ) )
1974 kp2col( ki ) = numroc( k2( ki )+3, hbl, mycol, 0,
1975 $ npcol )
1976 END IF
1977 k1( ki ) = k2( ki ) + 1
1978 istop = min( k1( ki )+rotn-mod( k1( ki ), rotn ), i-2 )
1979 istop = min( istop, k1( ki )+hbl-3-
1980 $ mod( k1( ki )-1, hbl ) )
1981 istop = min( istop, i2-2 )
1982 istop = max( istop, k1( ki ) )
1983* ISTOP = MIN( ISTOP , I-1 )
1984 k2( ki ) = istop
1985 IF( k1( ki ).EQ.istop ) THEN
1986 IF( ( mod( istop-1, hbl ).EQ.hbl-2 ) .AND.
1987 $ ( i-istop.GT.1 ) ) THEN
1988*
1989* Next step switches rows & cols
1990*
1991 icurrow( ki ) = mod( icurrow( ki )+1, nprow )
1992 icurcol( ki ) = mod( icurcol( ki )+1, npcol )
1993 END IF
1994 END IF
1995 410 CONTINUE
1996 IF( k2( ibulge ).LE.i-1 )
1997 $ GO TO 40
1998 END IF
1999*
2000 420 CONTINUE
2001*
2002* Failure to converge in remaining number of iterations
2003*
2004 info = i
2005 RETURN
2006*
2007 430 CONTINUE
2008*
2009 IF( l.EQ.i ) THEN
2010*
2011* H(I,I-1) is negligible: one eigenvalue has converged.
2012*
2013 CALL infog2l( i, i, desca, nprow, npcol, myrow, mycol, irow,
2014 $ icol, itmp1, itmp2 )
2015 IF( ( myrow.EQ.itmp1 ) .AND. ( mycol.EQ.itmp2 ) ) THEN
2016 wr( i ) = a( ( icol-1 )*lda+irow )
2017 ELSE
2018 wr( i ) = zero
2019 END IF
2020 wi( i ) = zero
2021 ELSE IF( l.EQ.i-1 ) THEN
2022*
2023* H(I-1,I-2) is negligible: a pair of eigenvalues have converged.
2024*
2025 CALL pselget( 'All', ' ', h11, a, l, l, desca )
2026 CALL pselget( 'All', ' ', h21, a, i, l, desca )
2027 CALL pselget( 'All', ' ', h12, a, l, i, desca )
2028 CALL pselget( 'All', ' ', h22, a, i, i, desca )
2029 CALL slanv2( h11, h12, h21, h22, wr( l ), wi( l ), wr( i ),
2030 $ wi( i ), cs, sn )
2031 IF( node .NE. 0 ) THEN
2032 wr( l ) = zero
2033 wr( i ) = zero
2034 wi( l ) = zero
2035 wi( i ) = zero
2036 ENDIF
2037 ELSE
2038*
2039* Find the eigenvalues in H(L:I,L:I), L < I-1
2040*
2041 jblk = i - l + 1
2042 IF( jblk.LE.2*iblk ) THEN
2043 CALL pslacp3( i-l+1, l, a, desca, s1, 2*iblk, 0, 0, 0 )
2044 CALL slahqr( .false., .false., jblk, 1, jblk, s1, 2*iblk,
2045 $ wr( l ), wi( l ), 1, jblk, z, ldz, ierr )
2046 IF( node.NE.0 ) THEN
2047*
2048* Erase the eigenvalues
2049*
2050 DO 440 k = l, i
2051 wr( k ) = zero
2052 wi( k ) = zero
2053 440 CONTINUE
2054 END IF
2055 END IF
2056 END IF
2057*
2058* Decrement number of remaining iterations, and return to start of
2059* the main loop with new value of I.
2060*
2061 itn = itn - its
2062 IF( m.EQ.l-10 ) THEN
2063 i = l - 1
2064 ELSE
2065 i = m
2066 END IF
2067* I = L - 1
2068 GO TO 10
2069*
2070 450 CONTINUE
2071 CALL sgsum2d( contxt, 'All', ' ', n, 1, wr, n, -1, -1 )
2072 CALL sgsum2d( contxt, 'All', ' ', n, 1, wi, n, -1, -1 )
2073 RETURN
2074*
2075* END OF PSLAHQR
2076*
2077 END
subroutine infog1l(gindx, nb, nprocs, myroc, isrcproc, lindx, rocsrc)
Definition infog1l.f:3
subroutine infog2l(grindx, gcindx, desc, nprow, npcol, myrow, mycol, lrindx, lcindx, rsrc, csrc)
Definition infog2l.f:3
#define max(A, B)
Definition pcgemr.c:180
#define min(A, B)
Definition pcgemr.c:181
subroutine pselget(scope, top, alpha, a, ia, ja, desca)
Definition pselget.f:2
subroutine pslabad(ictxt, small, large)
Definition pslabad.f:2
subroutine pslaconsb(a, desca, i, l, m, h44, h33, h43h34, buf, lwork)
Definition pslaconsb.f:3
subroutine pslacp3(m, i, a, desca, b, ldb, ii, jj, rev)
Definition pslacp3.f:2
subroutine pslahqr(wantt, wantz, n, ilo, ihi, a, desca, wr, wi, iloz, ihiz, z, descz, work, lwork, iwork, ilwork, info)
Definition pslahqr.f:4
subroutine pslasmsub(a, desca, i, l, k, smlnum, buf, lwork)
Definition pslasmsub.f:2
subroutine pslawil(ii, jj, m, a, desca, h44, h33, h43h34, v)
Definition pslawil.f:2
subroutine pxerbla(ictxt, srname, info)
Definition pxerbla.f:2
subroutine slaref(type, a, lda, wantz, z, ldz, block, irow1, icol1, istart, istop, itmp1, itmp2, liloz, lihiz, vecs, v2, v3, t1, t2, t3)
Definition slaref.f:4
subroutine slasorte(s, lds, j, out, info)
Definition slasorte.f:2