OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
pdlahqr.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine pdlahqr (wantt, wantz, n, ilo, ihi, a, desca, wr, wi, iloz, ihiz, z, descz, work, lwork, iwork, ilwork, info)

Function/Subroutine Documentation

◆ pdlahqr()

subroutine pdlahqr ( logical wantt,
logical wantz,
integer n,
integer ilo,
integer ihi,
double precision, dimension( * ) a,
integer, dimension( * ) desca,
double precision, dimension( * ) wr,
double precision, dimension( * ) wi,
integer iloz,
integer ihiz,
double precision, dimension( * ) z,
integer, dimension( * ) descz,
double precision, dimension( * ) work,
integer lwork,
integer, dimension( * ) iwork,
integer ilwork,
integer info )

Definition at line 1 of file pdlahqr.f.

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 DOUBLE PRECISION A( * ), WI( * ), WORK( * ), WR( * ), Z( * )
16* ..
17*
18* Purpose
19* =======
20*
21* PDLAHQR 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). PDLAHQR 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) DOUBLE PRECISION 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) DOUBLE PRECISION array,
115* dimension (N)
116* WI (global replicated output) DOUBLE PRECISION 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) DOUBLE PRECISION 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) DOUBLE PRECISION 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: PDLAHQR 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* PDLACONSB -> To determine where to start each iteration
188* PDLAWIL -> Given the shift, get the transformation
189* DLASORTE -> Pair up eigenvalues so that reals are paired.
190* PDLACP3 -> Parallel array to local replicated array copy &
191* back.
192* DLAREF -> Row/column reflector applier. Core routine
193* here.
194* PDLASMSUB -> 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 DOUBLE PRECISION ZERO, ONE, HALF
238 parameter( zero = 0.0d+0, one = 1.0d+0, half = 0.5d+0 )
239 DOUBLE PRECISION CONST
240 parameter( const = 1.50d+0 )
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 DOUBLE PRECISION 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 DOUBLE PRECISION S1( 2*IBLK, 2*IBLK ), SMALLA( 6, 6, IBLK ),
265 $ VCOPY( 3 )
266* ..
267* .. External Functions ..
268 INTEGER ILCM, NUMROC
269 DOUBLE PRECISION PDLAMCH
270 EXTERNAL ilcm, numroc, pdlamch
271* ..
272* .. External Subroutines ..
273 EXTERNAL blacs_gridinfo, dcopy, dgebr2d, dgebs2d,
274 $ dgerv2d, dgesd2d, dgsum2d, dlahqr, dlaref,
275 $ dlarfg, dlasorte, igamn2d, infog1l, infog2l,
276 $ pdlabad, pdlaconsb, pdlacp3, pdlasmsub,
277 $ pdlawil, pxerbla, dlanv2
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, 'PDLAHQR', -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 = pdlamch( contxt, 'SAFE MINIMUM' )
389 ovfl = one / unfl
390 CALL pdlabad( contxt, unfl, ovfl )
391 ulp = pdlamch( 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 pdlasmsub( 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 pdlacp3( 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 dlahqr( .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* PDLACONSB is the routine that does this.
521*
522c CALL PDLACONSB( 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 PDLACONSB 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 dlasorte( 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 pdlawil( 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 pdlawil( 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 pdlacp3( 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 pdlacp3( 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* DLAHQR 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 dcopy( 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 dlarfg( 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 dcopy( 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 dlarfg( 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 dlarfg( 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 dlarfg( 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 dgerv2d( 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 dgesd2d( 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 dgebr2d( 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 dgebs2d( 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 dgebr2d( 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 dgebs2d( 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 dgebr2d( 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 DGEBS2D( 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.EQ..AND..GT..AND. IF( ( MYCOLICURCOL( KI ) ) ( NPROW1 )
1012.LE. $ ( ISTARTISTOP ) ) THEN
1013.NE. IF( MYROWICURROW( KI ) ) THEN
1014 CALL DGEBR2D( 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 DGEBS2D( 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.EQ..AND. IF( ( MYROWICURROW( KI ) )
1036.EQ..AND. $ ( MYCOLICURCOL( KI ) )
1037.EQ..AND..LT. $ ( MODKM1HBL-2 ) ( ISTARTI-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.EQ. IF( NR3 ) 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.EQ..AND. IF( ( MOD( ISTART-1, HBL )HBL-1 )
1075.LE..AND. $ ( ISTARTISTOP )
1076.EQ..AND. $ ( MYROWICURROW( KI ) )
1077.EQ. $ ( MYCOLICURCOL( 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.EQ. IF( NR3 ) 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.EQ..AND. IF( ( MYROWICURROW( KI ) )
1116.EQ..AND. $ ( MYCOLICURCOL( KI ) )
1117.EQ..AND..EQ. $ ( ( ( MODKM1HBL-2 ) ( ISTARTI-
1118.OR..LT..AND..LE. $ 1 ) ) ( ( MODKM1HBL-2 ) ( ISTARTI-
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.LT. IF( KISTOP ) THEN
1134*
1135* Do some work so next step is ready...
1136*
1137 T2 = T1*V2
1138 T3 = T1*V3
1139 CALL DLAREF( '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.EQ..AND. IF( ( NR3 ) ( MOD( K-1,
1149.LT. $ HBL )HBL-2 ) ) THEN
1150 T2 = T1*V2
1151 T3 = T1*V3
1152 CALL DLAREF( '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.GE..AND..LE. IF( ( MODKM1HBL-2 ) ( KI-1 ) ) THEN
1168.EQ..AND..LT. IF( ( MODKM1HBL-2 ) ( KI-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 PDLACP3( MIN( 6, N-K+2 ), K-1, A, DESCA,
1175 $ SMALLA( 1, 1, KI ), 6, ITMP1, ITMP2,
1176 $ 1 )
1177*
1178 END IF
1179.EQ. IF( MODKM1HBL-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 PDLACP3( 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.NE..AND. IF( ( MYROWICURROW( KI ) )
1197.NE. $ ( DOWNICURROW( KI ) ) )GO TO 160
1198 ISTART = MAX( K1( KI ), M )
1199 ISTOP = MIN( K2( KI ), I-1 )
1200*
1201.GT..AND. IF( ( ISTOPISTART )
1202.LT..AND. $ ( MOD( ISTART-1, HBL )HBL-2 )
1203.EQ. $ ( ICURROW( KI )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 DLAREF( '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.GT. IF( KROW( KI )KP2ROW( KI ) )
1218 $ GO TO 180
1219.NE..AND. IF( ( MYROWICURROW( KI ) )
1220.NE. $ ( DOWNICURROW( KI ) ) )GO TO 180
1221 ISTART = MAX( K1( KI ), M )
1222 ISTOP = MIN( K2( KI ), I-1 )
1223.EQ..OR. IF( ( ISTARTISTOP )
1224.GE..OR. $ ( MOD( ISTART-1, HBL )HBL-2 )
1225.NE. $ ( ICURROW( KI )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.EQ..AND..LE. IF( ( NR3 ) ( KROW( KI )
1232 $ KP2ROW( KI ) ) ) THEN
1233.LT..AND. IF( ( KISTOP )
1234.LT. $ ( MOD( K-1, HBL )HBL-2 ) ) THEN
1235 ITMP1 = MIN( K2( KI )+1, I-1 ) + 1
1236 ELSE
1237.LT. IF( MOD( K-1, HBL )HBL-2 ) THEN
1238 ITMP1 = MIN( K2( KI )+1, I-1 ) + 1
1239 END IF
1240.EQ. IF( MOD( K-1, HBL )HBL-2 ) THEN
1241 ITMP1 = MIN( K+4, I2 ) + 1
1242 END IF
1243.EQ. IF( MOD( K-1, HBL )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.LT..OR. IF( ( MOD( K-1, HBL )HBL-2 )
1256.EQ. $ ( NPROW1 ) ) THEN
1257 T2 = T1*V2
1258 T3 = T1*V3
1259 CALL DLAREF( '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.EQ..AND. IF( ( MOD( K-1, HBL )HBL-2 )
1266.GT. $ ( NPROW1 ) ) THEN
1267.EQ. IF( IROW1IROW2 ) THEN
1268 CALL DGESD2D( CONTXT, 1, ICOL2-ICOL1+1,
1269 $ A( ( ICOL1-1 )*LDA+IROW2 ),
1270 $ LDA, UP, MYCOL )
1271 END IF
1272 END IF
1273.EQ..AND. IF( ( MOD( K-1, HBL )HBL-1 )
1274.GT. $ ( NPROW1 ) ) THEN
1275.EQ. IF( IROW1IROW2 ) THEN
1276 CALL DGESD2D( 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.GT. IF( KROW( KI )KP2ROW( KI ) )
1288 $ GO TO 220
1289.NE..AND. IF( ( MYROWICURROW( KI ) )
1290.NE. $ ( DOWNICURROW( KI ) ) )GO TO 220
1291 ISTART = MAX( K1( KI ), M )
1292 ISTOP = MIN( K2( KI ), I-1 )
1293.EQ..OR. IF( ( ISTARTISTOP )
1294.GE..OR. $ ( MOD( ISTART-1, HBL )HBL-2 )
1295.NE. $ ( ICURROW( KI )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.EQ..AND..LE. IF( ( NR3 ) ( KROW( KI )
1302 $ KP2ROW( KI ) ) ) THEN
1303.LT..AND. IF( ( KISTOP )
1304.LT. $ ( MOD( K-1, HBL )HBL-2 ) ) THEN
1305 ITMP1 = MIN( K2( KI )+1, I-1 ) + 1
1306 ELSE
1307.LT. IF( MOD( K-1, HBL )HBL-2 ) THEN
1308 ITMP1 = MIN( K2( KI )+1, I-1 ) + 1
1309 END IF
1310.EQ. IF( MOD( K-1, HBL )HBL-2 ) THEN
1311 ITMP1 = MIN( K+4, I2 ) + 1
1312 END IF
1313.EQ. IF( MOD( K-1, HBL )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.EQ..AND. IF( ( MOD( K-1, HBL )HBL-2 )
1324.GT. $ ( NPROW1 ) ) THEN
1325.NE. IF( IROW1IROW2 ) THEN
1326 CALL DGERV2D( 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 DGESD2D( CONTXT, 1, ICOL2-ICOL1+1,
1343 $ WORK( IRBUF+1 ), 1, DOWN,
1344 $ MYCOL )
1345 END IF
1346 END IF
1347.EQ..AND. IF( ( MOD( K-1, HBL )HBL-1 )
1348.GT. $ ( NPROW1 ) ) THEN
1349.NE. IF( IROW1IROW2 ) THEN
1350 CALL DGERV2D( 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 DGESD2D( 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.GT. IF( KROW( KI )KP2ROW( KI ) )
1378 $ GO TO 240
1379.NE..AND. IF( ( MYROWICURROW( KI ) )
1380.NE. $ ( DOWNICURROW( KI ) ) )GO TO 240
1381 ISTART = MAX( K1( KI ), M )
1382 ISTOP = MIN( K2( KI ), I-1 )
1383.EQ..OR. IF( ( ISTARTISTOP )
1384.GE..OR. $ ( MOD( ISTART-1, HBL )HBL-2 )
1385.NE. $ ( ICURROW( KI )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.EQ..AND..LE. IF( ( NR3 ) ( KROW( KI )
1392 $ KP2ROW( KI ) ) ) THEN
1393.LT..AND. IF( ( KISTOP )
1394.LT. $ ( MOD( K-1, HBL )HBL-2 ) ) THEN
1395 ITMP1 = MIN( K2( KI )+1, I-1 ) + 1
1396 ELSE
1397.LT. IF( MOD( K-1, HBL )HBL-2 ) THEN
1398 ITMP1 = MIN( K2( KI )+1, I-1 ) + 1
1399 END IF
1400.EQ. IF( MOD( K-1, HBL )HBL-2 ) THEN
1401 ITMP1 = MIN( K+4, I2 ) + 1
1402 END IF
1403.EQ. IF( MOD( K-1, HBL )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.EQ..AND. IF( ( MOD( K-1, HBL )HBL-2 )
1414.GT. $ ( NPROW1 ) ) THEN
1415.EQ. IF( IROW1IROW2 ) THEN
1416 CALL DGERV2D( CONTXT, 1, ICOL2-ICOL1+1,
1417 $ A( ( ICOL1-1 )*LDA+IROW2 ),
1418 $ LDA, UP, MYCOL )
1419 END IF
1420 END IF
1421.EQ..AND. IF( ( MOD( K-1, HBL )HBL-1 )
1422.GT. $ ( NPROW1 ) ) THEN
1423.EQ. IF( IROW1IROW2 ) THEN
1424 CALL DGERV2D( 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.NE..AND. IF( ( MYCOLICURCOL( KI ) )
1439.NE. $ ( RIGHTICURCOL( KI ) ) )GO TO 260
1440 ISTART = MAX( K1( KI ), M )
1441 ISTOP = MIN( K2( KI ), I-1 )
1442*
1443.LT..OR..EQ. IF( ( ( MOD( ISTART-1, HBL )HBL-2 ) ( NPCOL
1444.AND..EQ..AND. $ 1 ) ) ( ICURCOL( KI )MYCOL )
1445.GE. $ ( I-ISTOP+13 ) ) THEN
1446 K = ISTART
1447.LT..AND. IF( ( KISTOP ) ( MOD( K-1,
1448.LT. $ HBL )HBL-2 ) ) THEN
1449 ITMP1 = MIN( ISTART+1, I ) - 1
1450 ELSE
1451.LT. IF( MOD( K-1, HBL )HBL-2 ) THEN
1452 ITMP1 = MIN( K+3, I )
1453 END IF
1454.EQ. IF( MOD( K-1, HBL )HBL-2 ) THEN
1455 ITMP1 = MAX( I1, K-1 ) - 1
1456 END IF
1457.EQ. IF( MOD( K-1, HBL )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.LE. IF( IROW1IROW2 ) THEN
1466 ITMP2 = IROW2
1467 ELSE
1468 ITMP2 = -1
1469 END IF
1470 CALL DLAREF( '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.LT. IF( MOD( K-1, HBL )HBL-2 ) THEN
1476*
1477* Do from ITMP1+1 to MIN(K+3,I)
1478*
1479.LT. IF( MOD( K-1, HBL )HBL-3 ) THEN
1480 IROW1 = ITMP2 + 1
1481.EQ. IF( MOD( ( ITMP1 / HBL ), NPROW )MYROW )
1482 $ THEN
1483.GT. IF( ITMP20 ) 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 DLAREF( '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.GT. IF( KCOL( KI )KP2COL( KI ) )
1514 $ GO TO 320
1515.NE..AND. IF( ( MYCOLICURCOL( KI ) )
1516.NE. $ ( RIGHTICURCOL( KI ) ) )GO TO 320
1517 ISTART = MAX( K1( KI ), M )
1518 ISTOP = MIN( K2( KI ), I-1 )
1519.GE. IF( MOD( ISTART-1, HBL )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.EQ..AND..LE. IF( ( NR3 ) ( KCOL( KI )KP2COL( KI ) ) )
1538 $ THEN
1539*
1540.LT..AND. IF( ( KISTOP )
1541.LT. $ ( MOD( K-1, HBL )HBL-2 ) ) THEN
1542 ITMP1 = MIN( ISTART+1, I ) - 1
1543 ELSE
1544.LT. IF( MOD( K-1, HBL )HBL-2 ) THEN
1545 ITMP1 = MIN( K+3, I )
1546 END IF
1547.EQ. IF( MOD( K-1, HBL )HBL-2 ) THEN
1548 ITMP1 = MAX( I1, K-1 ) - 1
1549 END IF
1550.EQ. IF( MOD( K-1, HBL )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.EQ..AND. IF( ( MOD( K-1, HBL )HBL-2 )
1560.GT. $ ( NPCOL1 ) ) THEN
1561.EQ. IF( ICOL1ICOL2 ) THEN
1562 CALL DGESD2D( CONTXT, IROW2-IROW1+1, 1,
1563 $ A( ( ICOL1-1 )*LDA+IROW1 ),
1564 $ LDA, MYROW, LEFT )
1565 CALL DGERV2D( CONTXT, IROW2-IROW1+1, 1,
1566 $ A( ( ICOL1-1 )*LDA+IROW1 ),
1567 $ LDA, MYROW, LEFT )
1568 ELSE
1569 CALL DGERV2D( 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 DGESD2D( CONTXT, IROW2-IROW1+1, 1,
1586 $ WORK( ICBUF+1 ), IROW2-IROW1+1,
1587 $ MYROW, RIGHT )
1588 END IF
1589 END IF
1590.EQ..AND. IF( ( MOD( K-1, HBL )HBL-1 )
1591.GT. $ ( NPCOL1 ) ) THEN
1592.EQ. IF( ICOL1ICOL2 ) THEN
1593 CALL DGESD2D( CONTXT, IROW2-IROW1+1, 1,
1594 $ A( ( ICOL1-1 )*LDA+IROW1 ),
1595 $ LDA, MYROW, RIGHT )
1596 CALL DGERV2D( CONTXT, IROW2-IROW1+1, 1,
1597 $ A( ( ICOL1-1 )*LDA+IROW1 ),
1598 $ LDA, MYROW, RIGHT )
1599 ELSE
1600 CALL DGERV2D( 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 DGESD2D( 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.AND. IF( ( WANTZ ) ( MOD( K-1,
1624.GE..AND..GT. $ HBL )HBL-2 ) ( NPCOL1 ) ) THEN
1625*
1626* Accumulate transformations in the matrix Z
1627*
1628 IROW1 = LILOZ
1629 IROW2 = LIHIZ
1630.EQ. IF( MOD( K-1, HBL )HBL-2 ) THEN
1631.EQ. IF( ICOL1ICOL2 ) THEN
1632 CALL DGESD2D( CONTXT, IROW2-IROW1+1, 1,
1633 $ Z( ( ICOL1-1 )*LDZ+IROW1 ),
1634 $ LDZ, MYROW, LEFT )
1635 CALL DGERV2D( CONTXT, IROW2-IROW1+1, 1,
1636 $ Z( ( ICOL1-1 )*LDZ+IROW1 ),
1637 $ LDZ, MYROW, LEFT )
1638 ELSE
1639 CALL DGERV2D( 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 DGESD2D( CONTXT, IROW2-IROW1+1, 1,
1657 $ WORK( ICBUF+1 ),
1658 $ IROW2-IROW1+1, MYROW,
1659 $ RIGHT )
1660 END IF
1661 END IF
1662.EQ. IF( MOD( K-1, HBL )HBL-1 ) THEN
1663.EQ. IF( ICOL1ICOL2 ) THEN
1664 CALL DGESD2D( CONTXT, IROW2-IROW1+1, 1,
1665 $ Z( ( ICOL1-1 )*LDZ+IROW1 ),
1666 $ LDZ, MYROW, RIGHT )
1667 CALL DGERV2D( CONTXT, IROW2-IROW1+1, 1,
1668 $ Z( ( ICOL1-1 )*LDZ+IROW1 ),
1669 $ LDZ, MYROW, RIGHT )
1670 ELSE
1671 CALL DGERV2D( 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 DGESD2D( 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.EQ. IF( ICURCOL( KI )MYCOL ) THEN
1694.EQ..OR..EQ. IF( ( ISPEC0 ) ( NPCOL1 ) ) THEN
1695 LOCALK2( KI ) = LOCALK2( KI ) + 1
1696 END IF
1697 ELSE
1698.EQ..AND. IF( ( MOD( K-1, HBL )HBL-1 )
1699.EQ. $ ( ICURCOL( KI )RIGHT ) ) THEN
1700.GT. IF( KM ) THEN
1701 LOCALK2( KI ) = LOCALK2( KI ) + 2
1702 ELSE
1703 LOCALK2( KI ) = LOCALK2( KI ) + 1
1704 END IF
1705 END IF
1706.EQ..AND. IF( ( MOD( K-1, HBL )HBL-2 )
1707.EQ..AND..EQ. $ ( I-K2 ) ( ICURCOL( KI )
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.GE. IF( MOD( ISTART-1, HBL )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.EQ. IF( NR2 ) THEN
1744.EQ. IF ( ICURROW( KI )MYROW ) THEN
1745 T2 = T1*V2
1746 END IF
1747.EQ. IF ( ICURCOL( KI )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.EQ. IF( ICURROW( KI )MYROW ) THEN
1761.EQ..OR..EQ..OR. IF( ( ISPEC0 ) ( NPROW1 )
1762.EQ. $ ( MOD( K-1, HBL )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.EQ. IF( MOD( K-1, HBL )HBL-1 ) THEN
1772 CALL DGERV2D( 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 DGESD2D( CONTXT, 1, LIHIH-LILOH+1,
1784 $ WORK( IRBUF+1 ), 1, UP,
1785 $ MYCOL )
1786 END IF
1787 END IF
1788 ELSE
1789.EQ..AND. IF( ( MOD( K-1, HBL )HBL-1 )
1790.EQ. $ ( ICURROW( KI )DOWN ) ) THEN
1791 CALL DGESD2D( CONTXT, 1, LIHIH-LILOH+1,
1792 $ A( ( LILOH-1 )*LDA+ITMP1 ),
1793 $ LDA, DOWN, MYCOL )
1794 CALL DGERV2D( 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.EQ. IF( ICURCOL( KI )MYCOL ) THEN
1808* LOCAL A(LILOZ:LIHIZ,LOCALK2:LOCALK2+2)
1809.EQ..OR..EQ..OR. IF( ( ISPEC0 ) ( NPCOL1 )
1810.EQ. $ ( MOD( K-1, HBL )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.EQ. IF( MOD( K-1, HBL )HBL-1 ) THEN
1825 CALL DGERV2D( 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 DGESD2D( CONTXT, LIHIH-LILOH+1, 1,
1837 $ WORK( ICBUF+1 ),
1838 $ LIHIH-LILOH+1, MYROW, LEFT )
1839 END IF
1840 END IF
1841 ELSE
1842.EQ..AND. IF( ( MOD( K-1, HBL )HBL-1 )
1843.EQ. $ ( ICURCOL( KI )RIGHT ) ) THEN
1844 ITMP1 = KCOL( KI )
1845 CALL DGESD2D( 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 DGERV2D( 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.EQ. IF( ICURCOL( KI )MYCOL ) THEN
1862* LOCAL Z(LILOZ:LIHIZ,LOCALK2:LOCALK2+2)
1863.EQ..OR..EQ..OR. IF( ( ISPEC0 ) ( NPCOL1 )
1864.EQ. $ ( MOD( K-1, HBL )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.EQ. IF( MOD( K-1, HBL )HBL-1 ) THEN
1879 CALL DGERV2D( 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 DGESD2D( 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.EQ..AND. IF( ( MOD( K-1, HBL )HBL-1 )
1901.EQ. $ ( ICURCOL( KI )RIGHT ) ) THEN
1902 ITMP1 = KCOL( KI )
1903 ITMP1 = ( ITMP1-1 )*LDZ
1904 CALL DGESD2D( CONTXT, LIHIZ-LILOZ+1, 1,
1905 $ Z( LILOZ+ITMP1 ), LDZ,
1906 $ MYROW, RIGHT )
1907 CALL DGERV2D( 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.EQ. IF( NPROW1 ) 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.LT..AND. IF( ( MOD( K1( KI )-1, HBL )HBL-2 )
1924.EQ..AND..GT. $ ( ICURROW( KI )MYROW ) ( NPROW1 ) )
1925 $ THEN
1926 KROW( KI ) = KROW( KI ) + K2( KI ) - K1( KI ) + 1
1927 END IF
1928.LT..AND. IF( ( MOD( K2( KI ), HBL )HBL-2 )
1929.EQ..AND..GT. $ ( ICURROW( KI )MYROW ) ( NPROW1 ) )
1930 $ THEN
1931 KP2ROW( KI ) = KP2ROW( KI ) + K2( KI ) - K1( KI ) + 1
1932 END IF
1933.GE..AND. IF( ( MOD( K1( KI )-1, HBL )HBL-2 )
1934.EQ..OR..EQ. $ ( ( MYROWICURROW( KI ) ) ( DOWN
1935.AND..GT. $ ICURROW( KI ) ) ) ( NPROW1 ) ) 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.GE..AND. IF( ( MOD( K2( KI ), HBL )HBL-2 )
1941.EQ..OR..EQ. $ ( ( MYROWICURROW( KI ) ) ( UP
1942.AND..GT. $ ICURROW( KI ) ) ) ( NPROW1 ) ) 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.EQ. IF( NPCOL1 ) 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.LT..AND. IF( ( MOD( K1( KI )-1, HBL )HBL-2 )
1953.EQ..AND..GT. $ ( ICURCOL( KI )MYCOL ) ( NPCOL1 ) )
1954 $ THEN
1955 KCOL( KI ) = KCOL( KI ) + K2( KI ) - K1( KI ) + 1
1956 END IF
1957.LT..AND. IF( ( MOD( K2( KI ), HBL )HBL-2 )
1958.EQ..AND..GT. $ ( ICURCOL( KI )MYCOL ) ( NPCOL1 ) )
1959 $ THEN
1960 KP2COL( KI ) = KP2COL( KI ) + K2( KI ) - K1( KI ) + 1
1961 END IF
1962.GE..AND. IF( ( MOD( K1( KI )-1, HBL )HBL-2 )
1963.EQ..OR..EQ. $ ( ( MYCOLICURCOL( KI ) ) ( RIGHT
1964.AND..GT. $ ICURCOL( KI ) ) ) ( NPCOL1 ) ) 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.GE..AND. IF( ( MOD( K2( KI ), HBL )HBL-2 )
1970.EQ..OR..EQ. $ ( ( MYCOLICURCOL( KI ) ) ( LEFT
1971.AND..GT. $ ICURCOL( KI ) ) ) ( NPCOL1 ) ) 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.EQ. IF( K1( KI )ISTOP ) THEN
1986.EQ..AND. IF( ( MOD( ISTOP-1, HBL )HBL-2 )
1987.GT. $ ( I-ISTOP1 ) ) 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.LE. IF( K2( IBULGE )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.EQ. IF( LI ) 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.EQ..AND..EQ. IF( ( MYROWITMP1 ) ( MYCOLITMP2 ) ) THEN
2016 WR( I ) = A( ( ICOL-1 )*LDA+IROW )
2017 ELSE
2018 WR( I ) = ZERO
2019 END IF
2020 WI( I ) = ZERO
2021.EQ. ELSE IF( LI-1 ) THEN
2022*
2023* H(I-1,I-2) is negligible: a pair of eigenvalues have converged.
2024*
2025 CALL PDELGET( 'all', ' ', H11, A, L, L, DESCA )
2026 CALL PDELGET( 'all', ' ', H21, A, I, L, DESCA )
2027 CALL PDELGET( 'all', ' ', H12, A, L, I, DESCA )
2028 CALL PDELGET( 'all', ' ', H22, A, I, I, DESCA )
2029 CALL DLANV2( H11, H12, H21, H22, WR( L ), WI( L ), WR( I ),
2030 $ WI( I ), CS, SN )
2031.NE. IF( NODE 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.LE. IF( JBLK2*IBLK ) THEN
2043 CALL PDLACP3( I-L+1, L, A, DESCA, S1, 2*IBLK, 0, 0, 0 )
2044 CALL DLAHQR( .FALSE., .FALSE., JBLK, 1, JBLK, S1, 2*IBLK,
2045 $ WR( L ), WI( L ), 1, JBLK, Z, LDZ, IERR )
2046.NE. IF( NODE0 ) 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.EQ. IF( ML-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 DGSUM2D( CONTXT, 'all', ' ', N, 1, WR, N, -1, -1 )
2072 CALL DGSUM2D( CONTXT, 'all', ' ', N, 1, WI, N, -1, -1 )
2073 RETURN
2074*
2075* END OF PDLAHQR
2076*
dlaref
subroutine dlaref(type, a, lda, wantz, z, ldz, block, irow1, icol1, istart, istop, itmp1, itmp2, liloz, lihiz, vecs, v2, v3, t1, t2, t3)
Definition dlaref.f:4
dlasorte
subroutine dlasorte(s, lds, j, out, info)
Definition dlasorte.f:2
dlahqr
subroutine dlahqr(wantt, wantz, n, ilo, ihi, h, ldh, wr, wi, iloz, ihiz, z, ldz, info)
DLAHQR computes the eigenvalues and Schur factorization of an upper Hessenberg matrix,...
Definition dlahqr.f:207
dlanv2
subroutine dlanv2(a, b, c, d, rt1r, rt1i, rt2r, rt2i, cs, sn)
DLANV2 computes the Schur factorization of a real 2-by-2 nonsymmetric matrix in standard form.
Definition dlanv2.f:127
dlarfg
subroutine dlarfg(n, alpha, x, incx, tau)
DLARFG generates an elementary reflector (Householder matrix).
Definition dlarfg.f:106
dcopy
subroutine dcopy(n, dx, incx, dy, incy)
DCOPY
Definition dcopy.f:82
ilcm
integer function ilcm(m, n)
Definition ilcm.f:2
infog1l
subroutine infog1l(gindx, nb, nprocs, myroc, isrcproc, lindx, rocsrc)
Definition infog1l.f:3
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
dgebs2d
subroutine dgebs2d(contxt, scope, top, m, n, a, lda)
Definition mpi.f:1082
pxerbla
subroutine pxerbla(contxt, srname, info)
Definition mpi.f:1600
dgebr2d
subroutine dgebr2d(contxt, scope, top, m, n, a, lda)
Definition mpi.f:1123
infog2l
subroutine infog2l(grindx, gcindx, desc, nprow, npcol, myrow, mycol, lrindx, lcindx, rsrc, csrc)
Definition mpi.f:937
blacs_gridinfo
subroutine blacs_gridinfo(cntxt, nprow, npcol, myrow, mycol)
Definition mpi.f:754
numroc
integer function numroc(n, nb, iproc, isrcproc, nprocs)
Definition mpi.f:786
pdlamch
double precision function pdlamch(ictxt, cmach)
Definition pdblastst.f:6769
pdlabad
subroutine pdlabad(ictxt, small, large)
Definition pdlabad.f:2
pdlaconsb
subroutine pdlaconsb(a, desca, i, l, m, h44, h33, h43h34, buf, lwork)
Definition pdlaconsb.f:3
pdlacp3
subroutine pdlacp3(m, i, a, desca, b, ldb, ii, jj, rev)
Definition pdlacp3.f:2
pdlasmsub
subroutine pdlasmsub(a, desca, i, l, k, smlnum, buf, lwork)
Definition pdlasmsub.f:2
pdlawil
subroutine pdlawil(ii, jj, m, a, desca, h44, h33, h43h34, v)
Definition pdlawil.f:2