ScaLAPACK 2.1  2.1
ScaLAPACK: Scalable Linear Algebra PACKage
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 *
522 c CALL PSLACONSB( A, DESCA, I, L, M, H44, H33, H43H34,
523 c $ 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
slaref
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
max
#define max(A, B)
Definition: pcgemr.c:180
pslacp3
subroutine pslacp3(M, I, A, DESCA, B, LDB, II, JJ, REV)
Definition: pslacp3.f:2
infog1l
subroutine infog1l(GINDX, NB, NPROCS, MYROC, ISRCPROC, LINDX, ROCSRC)
Definition: infog1l.f:3
pslabad
subroutine pslabad(ICTXT, SMALL, LARGE)
Definition: pslabad.f:2
infog2l
subroutine infog2l(GRINDX, GCINDX, DESC, NPROW, NPCOL, MYROW, MYCOL, LRINDX, LCINDX, RSRC, CSRC)
Definition: infog2l.f:3
pslaconsb
subroutine pslaconsb(A, DESCA, I, L, M, H44, H33, H43H34, BUF, LWORK)
Definition: pslaconsb.f:3
slasorte
subroutine slasorte(S, LDS, J, OUT, INFO)
Definition: slasorte.f:2
pselget
subroutine pselget(SCOPE, TOP, ALPHA, A, IA, JA, DESCA)
Definition: pselget.f:2
pslawil
subroutine pslawil(II, JJ, M, A, DESCA, H44, H33, H43H34, V)
Definition: pslawil.f:2
pslahqr
subroutine pslahqr(WANTT, WANTZ, N, ILO, IHI, A, DESCA, WR, WI, ILOZ, IHIZ, Z, DESCZ, WORK, LWORK, IWORK, ILWORK, INFO)
Definition: pslahqr.f:4
pslasmsub
subroutine pslasmsub(A, DESCA, I, L, K, SMLNUM, BUF, LWORK)
Definition: pslasmsub.f:2
pxerbla
subroutine pxerbla(ICTXT, SRNAME, INFO)
Definition: pxerbla.f:2
min
#define min(A, B)
Definition: pcgemr.c:181