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

Go to the source code of this file.

Functions/Subroutines

subroutine pctrevc (side, howmny, select, n, t, desct, vl, descvl, vr, descvr, mm, m, work, rwork, info)

Function/Subroutine Documentation

◆ pctrevc()

subroutine pctrevc ( character side,
character howmny,
logical, dimension( * ) select,
integer n,
complex, dimension( * ) t,
integer, dimension( * ) desct,
complex, dimension( * ) vl,
integer, dimension( * ) descvl,
complex, dimension( * ) vr,
integer, dimension( * ) descvr,
integer mm,
integer m,
complex, dimension( * ) work,
real, dimension( * ) rwork,
integer info )

Definition at line 1 of file pctrevc.f.

3*
4* -- ScaLAPACK routine (version 1.7) --
5* University of Tennessee, Knoxville, Oak Ridge National Laboratory,
6* and University of California, Berkeley.
7* July 31, 2001
8*
9* .. Scalar Arguments ..
10 CHARACTER HOWMNY, SIDE
11 INTEGER INFO, M, MM, N
12* ..
13* .. Array Arguments ..
14 LOGICAL SELECT( * )
15 INTEGER DESCT( * ), DESCVL( * ), DESCVR( * )
16 REAL RWORK( * )
17 COMPLEX T( * ), VL( * ), VR( * ), WORK( * )
18* ..
19*
20* Purpose
21* =======
22*
23* PCTREVC computes some or all of the right and/or left eigenvectors of
24* a complex upper triangular matrix T in parallel.
25*
26* The right eigenvector x and the left eigenvector y of T corresponding
27* to an eigenvalue w are defined by:
28*
29* T*x = w*x, y'*T = w*y'
30*
31* where y' denotes the conjugate transpose of the vector y.
32*
33* If all eigenvectors are requested, the routine may either return the
34* matrices X and/or Y of right or left eigenvectors of T, or the
35* products Q*X and/or Q*Y, where Q is an input unitary
36* matrix. If T was obtained from the Schur factorization of an
37* original matrix A = Q*T*Q', then Q*X and Q*Y are the matrices of
38* right or left eigenvectors of A.
39*
40* Notes
41* =====
42*
43* Each global data object is described by an associated description
44* vector. This vector stores the information required to establish
45* the mapping between an object element and its corresponding process
46* and memory location.
47*
48* Let A be a generic term for any 2D block cyclicly distributed array.
49* Such a global array has an associated description vector DESCA.
50* In the following comments, the character _ should be read as
51* "of the global array".
52*
53* NOTATION STORED IN EXPLANATION
54* --------------- -------------- --------------------------------------
55* DTYPE_A(global) DESCA( DTYPE_ )The descriptor type. In this case,
56* DTYPE_A = 1.
57* CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
58* the BLACS process grid A is distribu-
59* ted over. The context itself is glo-
60* bal, but the handle (the integer
61* value) may vary.
62* M_A (global) DESCA( M_ ) The number of rows in the global
63* array A.
64* N_A (global) DESCA( N_ ) The number of columns in the global
65* array A.
66* MB_A (global) DESCA( MB_ ) The blocking factor used to distribute
67* the rows of the array.
68* NB_A (global) DESCA( NB_ ) The blocking factor used to distribute
69* the columns of the array.
70* RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
71* row of the array A is distributed.
72* CSRC_A (global) DESCA( CSRC_ ) The process column over which the
73* first column of the array A is
74* distributed.
75* LLD_A (local) DESCA( LLD_ ) The leading dimension of the local
76* array. LLD_A >= MAX(1,LOCr(M_A)).
77*
78* Let K be the number of rows or columns of a distributed matrix,
79* and assume that its process grid has dimension r x c.
80* LOCr( K ) denotes the number of elements of K that a process
81* would receive if K were distributed over the r processes of its
82* process column.
83* Similarly, LOCc( K ) denotes the number of elements of K that a
84* process would receive if K were distributed over the c processes of
85* its process row.
86* The values of LOCr() and LOCc() may be determined via a call to the
87* ScaLAPACK tool function, NUMROC:
88* LOCr( M ) = NUMROC( M, MB_A, MYROW, RSRC_A, NPROW ),
89* LOCc( N ) = NUMROC( N, NB_A, MYCOL, CSRC_A, NPCOL ).
90* An upper bound for these quantities may be computed by:
91* LOCr( M ) <= ceil( ceil(M/MB_A)/NPROW )*MB_A
92* LOCc( N ) <= ceil( ceil(N/NB_A)/NPCOL )*NB_A
93*
94* Arguments
95* =========
96*
97* SIDE (global input) CHARACTER*1
98* = 'R': compute right eigenvectors only;
99* = 'L': compute left eigenvectors only;
100* = 'B': compute both right and left eigenvectors.
101*
102* HOWMNY (global input) CHARACTER*1
103* = 'A': compute all right and/or left eigenvectors;
104* = 'B': compute all right and/or left eigenvectors,
105* and backtransform them using the input matrices
106* supplied in VR and/or VL;
107* = 'S': compute selected right and/or left eigenvectors,
108* specified by the logical array SELECT.
109*
110* SELECT (global input) LOGICAL array, dimension (N)
111* If HOWMNY = 'S', SELECT specifies the eigenvectors to be
112* computed.
113* If HOWMNY = 'A' or 'B', SELECT is not referenced.
114* To select the eigenvector corresponding to the j-th
115* eigenvalue, SELECT(j) must be set to .TRUE..
116*
117* N (global input) INTEGER
118* The order of the matrix T. N >= 0.
119*
120* T (global input/output) COMPLEX array, dimension
121* (DESCT(LLD_),*)
122* The upper triangular matrix T. T is modified, but restored
123* on exit.
124*
125* DESCT (global and local input) INTEGER array of dimension DLEN_.
126* The array descriptor for the distributed matrix T.
127*
128* VL (global input/output) COMPLEX array, dimension
129* (DESCVL(LLD_),MM)
130* On entry, if SIDE = 'L' or 'B' and HOWMNY = 'B', VL must
131* contain an N-by-N matrix Q (usually the unitary matrix Q of
132* Schur vectors returned by CHSEQR).
133* On exit, if SIDE = 'L' or 'B', VL contains:
134* if HOWMNY = 'A', the matrix Y of left eigenvectors of T;
135* if HOWMNY = 'B', the matrix Q*Y;
136* if HOWMNY = 'S', the left eigenvectors of T specified by
137* SELECT, stored consecutively in the columns
138* of VL, in the same order as their
139* eigenvalues.
140* If SIDE = 'R', VL is not referenced.
141*
142* DESCVL (global and local input) INTEGER array of dimension DLEN_.
143* The array descriptor for the distributed matrix VL.
144*
145* VR (global input/output) COMPLEX array, dimension
146* (DESCVR(LLD_),MM)
147* On entry, if SIDE = 'R' or 'B' and HOWMNY = 'B', VR must
148* contain an N-by-N matrix Q (usually the unitary matrix Q of
149* Schur vectors returned by CHSEQR).
150* On exit, if SIDE = 'R' or 'B', VR contains:
151* if HOWMNY = 'A', the matrix X of right eigenvectors of T;
152* if HOWMNY = 'B', the matrix Q*X;
153* if HOWMNY = 'S', the right eigenvectors of T specified by
154* SELECT, stored consecutively in the columns
155* of VR, in the same order as their
156* eigenvalues.
157* If SIDE = 'L', VR is not referenced.
158*
159* DESCVR (global and local input) INTEGER array of dimension DLEN_.
160* The array descriptor for the distributed matrix VR.
161*
162* MM (global input) INTEGER
163* The number of columns in the arrays VL and/or VR. MM >= M.
164*
165* M (global output) INTEGER
166* The number of columns in the arrays VL and/or VR actually
167* used to store the eigenvectors. If HOWMNY = 'A' or 'B', M
168* is set to N. Each selected eigenvector occupies one
169* column.
170*
171* WORK (local workspace) COMPLEX array,
172* dimension ( 2*DESCT(LLD_) )
173* Additional workspace may be required if PCLATTRS is updated
174* to use WORK.
175*
176* RWORK (local workspace) REAL array,
177* dimension ( DESCT(LLD_) )
178*
179* INFO (global output) INTEGER
180* = 0: successful exit
181* < 0: if INFO = -i, the i-th argument had an illegal value
182*
183* Further Details
184* ===============
185*
186* The algorithm used in this program is basically backward (forward)
187* substitution. It is the hope that scaling would be used to make the
188* the code robust against possible overflow. But scaling has not yet
189* been implemented in PCLATTRS which is called by this routine to solve
190* the triangular systems. PCLATTRS just calls PCTRSV.
191*
192* Each eigenvector is normalized so that the element of largest
193* magnitude has magnitude 1; here the magnitude of a complex number
194* (x,y) is taken to be |x| + |y|.
195*
196* Further Details
197* ===============
198*
199* Implemented by Mark R. Fahey, June, 2000
200*
201* =====================================================================
202*
203* .. Parameters ..
204 REAL ZERO, ONE
205 parameter( zero = 0.0e+0, one = 1.0e+0 )
206 COMPLEX CZERO, CONE
207 parameter( czero = ( 0.0e+0, 0.0e+0 ),
208 $ cone = ( 1.0e+0, 0.0e+0 ) )
209 INTEGER BLOCK_CYCLIC_2D, DLEN_, DTYPE_, CTXT_, M_, N_,
210 $ MB_, NB_, RSRC_, CSRC_, LLD_
211 parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
212 $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
213 $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
214* ..
215* .. Local Scalars ..
216 LOGICAL ALLV, BOTHV, LEFTV, OVER, RIGHTV, SOMEV
217 INTEGER CONTXT, CSRC, I, ICOL, II, IROW, IS, ITMP1,
218 $ ITMP2, J, K, KI, LDT, LDVL, LDVR, LDW, MB,
219 $ MYCOL, MYROW, NB, NPCOL, NPROW, RSRC
220 REAL SELF
221 REAL OVFL, REMAXD, SCALE, SMLNUM, ULP, UNFL
222 COMPLEX CDUM, REMAXC, SHIFT
223* ..
224* .. Local Arrays ..
225 INTEGER DESCW( DLEN_ )
226 REAL SMIN( 1 )
227* ..
228* .. External Functions ..
229 LOGICAL LSAME
230 REAL PSLAMCH
231 EXTERNAL lsame, pslamch
232* ..
233* .. External Subroutines ..
234 EXTERNAL blacs_gridinfo, descinit, sgsum2d, igamn2d,
235 $ infog2l, pslabad, pscasum, pxerbla, pcamax,
236 $ pccopy, pcsscal, pcgemv, pclaset, pclattrs,
237 $ cgsum2d
238* ..
239* .. Intrinsic Functions ..
240 INTRINSIC abs, real, cmplx, conjg, aimag, max
241* ..
242* .. Statement Functions ..
243 REAL CABS1
244* ..
245* .. Statement Function definitions ..
246 cabs1( cdum ) = abs( real( cdum ) ) + abs( aimag( cdum ) )
247* ..
248* .. Executable Statements ..
249*
250* This is just to keep ftnchek happy
251 IF( block_cyclic_2d*csrc_*ctxt_*dlen_*dtype_*lld_*mb_*m_*nb_*n_*
252 $ rsrc_.LT.0 )RETURN
253*
254 contxt = desct( ctxt_ )
255 rsrc = desct( rsrc_ )
256 csrc = desct( csrc_ )
257 mb = desct( mb_ )
258 nb = desct( nb_ )
259 ldt = desct( lld_ )
260 ldw = ldt
261 ldvr = descvr( lld_ )
262 ldvl = descvl( lld_ )
263*
264 CALL blacs_gridinfo( contxt, nprow, npcol, myrow, mycol )
265 self = myrow*npcol + mycol
266*
267* Decode and test the input parameters
268*
269 bothv = lsame( side, 'B' )
270 rightv = lsame( side, 'R' ) .OR. bothv
271 leftv = lsame( side, 'l.OR.' ) BOTHV
272*
273 ALLV = LSAME( HOWMNY, 'a' )
274 OVER = LSAME( HOWMNY, 'b.OR.' ) LSAME( HOWMNY, 'o' )
275 SOMEV = LSAME( HOWMNY, 's' )
276*
277* Set M to the number of columns required to store the selected
278* eigenvectors.
279*
280 IF( SOMEV ) THEN
281 M = 0
282 DO 10 J = 1, N
283 IF( SELECT( J ) )
284 $ M = M + 1
285 10 CONTINUE
286 ELSE
287 M = N
288 END IF
289*
290 INFO = 0
291.NOT..AND..NOT. IF( RIGHTV LEFTV ) THEN
292 INFO = -1
293.NOT..AND..NOT..AND..NOT. ELSE IF( ALLV OVER SOMEV ) THEN
294 INFO = -2
295.LT. ELSE IF( N0 ) THEN
296 INFO = -4
297.LT. ELSE IF( MMM ) THEN
298 INFO = -11
299 END IF
300 CALL IGAMN2D( CONTXT, 'all', ' ', 1, 1, INFO, 1, ITMP1, ITMP2, -1,
301 $ -1, -1 )
302.LT. IF( INFO0 ) THEN
303 CALL PXERBLA( CONTXT, 'pctrevc', -INFO )
304 RETURN
305 END IF
306*
307* Quick return if possible.
308*
309.EQ. IF( N0 )
310 $ RETURN
311*
312* Set the constants to control overflow.
313*
314 UNFL = PSLAMCH( CONTXT, 'safe minimum' )
315 OVFL = ONE / UNFL
316 CALL PSLABAD( CONTXT, UNFL, OVFL )
317 ULP = PSLAMCH( CONTXT, 'precision' )
318 SMLNUM = UNFL*( N / ULP )
319*
320* Store the diagonal elements of T in working array WORK( LDW+1 ).
321*
322 DO 20 I = 1, N
323 CALL INFOG2L( I, I, DESCT, NPROW, NPCOL, MYROW, MYCOL, IROW,
324 $ ICOL, ITMP1, ITMP2 )
325.EQ..AND..EQ. IF( ( MYROWITMP1 ) ( MYCOLITMP2 ) ) THEN
326 WORK( LDW+IROW ) = T( ( ICOL-1 )*LDT+IROW )
327 END IF
328 20 CONTINUE
329*
330* Compute 1-norm of each column of strictly upper triangular
331* part of T to control overflow in triangular solver. Computed,
332* but not used. For use in PCLATTRS.
333*
334 RWORK( 1 ) = ZERO
335 DO 30 J = 2, N
336 CALL PSCASUM( J-1, RWORK( J ), T, 1, J, DESCT, 1 )
337 30 CONTINUE
338* I replicate the norms in RWORK. Should they be distributed
339* over the process rows?
340 CALL SGSUM2D( CONTXT, 'row', ' ', N, 1, RWORK, N, -1, -1 )
341*
342 IF( RIGHTV ) THEN
343*
344* Compute right eigenvectors.
345*
346* Need to set the distribution pattern of WORK
347*
348 CALL DESCINIT( DESCW, N, 1, NB, 1, RSRC, CSRC, CONTXT, LDW,
349 $ INFO )
350*
351 IS = M
352 DO 70 KI = N, 1, -1
353*
354 IF( SOMEV ) THEN
355.NOT. IF( SELECT( KI ) )
356 $ GO TO 70
357 END IF
358*
359 SMIN( 1 ) = ZERO
360 SHIFT = CZERO
361 CALL INFOG2L( KI, KI, DESCT, NPROW, NPCOL, MYROW, MYCOL,
362 $ IROW, ICOL, ITMP1, ITMP2 )
363.EQ..AND..EQ. IF( ( MYROWITMP1 ) ( MYCOLITMP2 ) ) THEN
364 SHIFT = T( ( ICOL-1 )*LDT+IROW )
365 SMIN( 1 ) = MAX( ULP*( CABS1( SHIFT ) ), SMLNUM )
366 END IF
367 CALL SGSUM2D( CONTXT, 'all', ' ', 1, 1, SMIN, 1, -1, -1 )
368 CALL CGSUM2D( CONTXT, 'all', ' ', 1, 1, SHIFT, 1, -1, -1 )
369*
370 CALL INFOG2L( 1, 1, DESCW, NPROW, NPCOL, MYROW, MYCOL, IROW,
371 $ ICOL, ITMP1, ITMP2 )
372.EQ..AND..EQ. IF( ( MYROWITMP1 ) ( MYCOLITMP2 ) ) THEN
373 WORK( 1 ) = CONE
374 END IF
375*
376* Form right-hand side. Distribute rhs onto first column
377* of processor grid.
378*
379.GT. IF( KI1 ) THEN
380 CALL PCCOPY( KI-1, T, 1, KI, DESCT, 1, WORK, 1, 1, DESCW,
381 $ 1 )
382 END IF
383 DO 40 K = 1, KI - 1
384 CALL INFOG2L( K, 1, DESCW, NPROW, NPCOL, MYROW, MYCOL,
385 $ IROW, ICOL, ITMP1, ITMP2 )
386.EQ..AND..EQ. IF( MYROWITMP1 MYCOLITMP2 ) THEN
387 WORK( IROW ) = -WORK( IROW )
388 END IF
389 40 CONTINUE
390*
391* Solve the triangular system:
392* (T(1:KI-1,1:KI-1) - T(KI,KI))*X = SCALE*WORK.
393*
394 DO 50 K = 1, KI - 1
395 CALL INFOG2L( K, K, DESCT, NPROW, NPCOL, MYROW, MYCOL,
396 $ IROW, ICOL, ITMP1, ITMP2 )
397.EQ..AND..EQ. IF( ( MYROWITMP1 ) ( MYCOLITMP2 ) ) THEN
398 T( ( ICOL-1 )*LDT+IROW ) = T( ( ICOL-1 )*LDT+IROW ) -
399 $ SHIFT
400.LT. IF( CABS1( T( ( ICOL-1 )*LDT+IROW ) )SMIN( 1 ) )
401 $ THEN
402 T( ( ICOL-1 )*LDT+IROW ) = CMPLX( SMIN( 1 ) )
403 END IF
404 END IF
405 50 CONTINUE
406*
407.GT. IF( KI1 ) THEN
408 CALL PCLATTRS( 'upper', 'no transpose', 'non-unit', 'y',
409 $ KI-1, T, 1, 1, DESCT, WORK, 1, 1, DESCW,
410 $ SCALE, RWORK, INFO )
411 CALL INFOG2L( KI, 1, DESCW, NPROW, NPCOL, MYROW, MYCOL,
412 $ IROW, ICOL, ITMP1, ITMP2 )
413.EQ..AND..EQ. IF( MYROWITMP1 MYCOLITMP2 ) THEN
414 WORK( IROW ) = CMPLX( SCALE )
415 END IF
416 END IF
417*
418* Copy the vector x or Q*x to VR and normalize.
419*
420.NOT. IF( OVER ) THEN
421 CALL PCCOPY( KI, WORK, 1, 1, DESCW, 1, VR, 1, IS, DESCVR,
422 $ 1 )
423*
424 CALL PCAMAX( KI, REMAXC, II, VR, 1, IS, DESCVR, 1 )
425 REMAXD = ONE / MAX( CABS1( REMAXC ), UNFL )
426 CALL PCSSCAL( KI, REMAXD, VR, 1, IS, DESCVR, 1 )
427*
428 CALL PCLASET( ' ', N-KI, 1, CZERO, CZERO, VR, KI+1, IS,
429 $ DESCVR )
430 ELSE
431.GT. IF( KI1 )
432 $ CALL PCGEMV( 'n', N, KI-1, CONE, VR, 1, 1, DESCVR,
433 $ WORK, 1, 1, DESCW, 1, CMPLX( SCALE ),
434 $ VR, 1, KI, DESCVR, 1 )
435*
436 CALL PCAMAX( N, REMAXC, II, VR, 1, KI, DESCVR, 1 )
437 REMAXD = ONE / MAX( CABS1( REMAXC ), UNFL )
438 CALL PCSSCAL( N, REMAXD, VR, 1, KI, DESCVR, 1 )
439 END IF
440*
441* Set back the original diagonal elements of T.
442*
443 DO 60 K = 1, KI - 1
444 CALL INFOG2L( K, K, DESCT, NPROW, NPCOL, MYROW, MYCOL,
445 $ IROW, ICOL, ITMP1, ITMP2 )
446.EQ..AND..EQ. IF( ( MYROWITMP1 ) ( MYCOLITMP2 ) ) THEN
447 T( ( ICOL-1 )*LDT+IROW ) = WORK( LDW+IROW )
448 END IF
449 60 CONTINUE
450*
451 IS = IS - 1
452 70 CONTINUE
453 END IF
454*
455 IF( LEFTV ) THEN
456*
457* Compute left eigenvectors.
458*
459* Need to set the distribution pattern of WORK
460*
461 CALL DESCINIT( DESCW, N, 1, MB, 1, RSRC, CSRC, CONTXT, LDW,
462 $ INFO )
463*
464 IS = 1
465 DO 110 KI = 1, N
466*
467 IF( SOMEV ) THEN
468.NOT. IF( SELECT( KI ) )
469 $ GO TO 110
470 END IF
471*
472 SMIN( 1 ) = ZERO
473 SHIFT = CZERO
474 CALL INFOG2L( KI, KI, DESCT, NPROW, NPCOL, MYROW, MYCOL,
475 $ IROW, ICOL, ITMP1, ITMP2 )
476.EQ..AND..EQ. IF( ( MYROWITMP1 ) ( MYCOLITMP2 ) ) THEN
477 SHIFT = T( ( ICOL-1 )*LDT+IROW )
478 SMIN( 1 ) = MAX( ULP*( CABS1( SHIFT ) ), SMLNUM )
479 END IF
480 CALL SGSUM2D( CONTXT, 'all', ' ', 1, 1, SMIN, 1, -1, -1 )
481 CALL CGSUM2D( CONTXT, 'all', ' ', 1, 1, SHIFT, 1, -1, -1 )
482*
483 CALL INFOG2L( N, 1, DESCW, NPROW, NPCOL, MYROW, MYCOL, IROW,
484 $ ICOL, ITMP1, ITMP2 )
485.EQ..AND..EQ. IF( ( MYROWITMP1 ) ( MYCOLITMP2 ) ) THEN
486 WORK( IROW ) = CONE
487 END IF
488*
489* Form right-hand side.
490*
491.LT. IF( KIN ) THEN
492 CALL PCCOPY( N-KI, T, KI, KI+1, DESCT, N, WORK, KI+1, 1,
493 $ DESCW, 1 )
494 END IF
495 DO 80 K = KI + 1, N
496 CALL INFOG2L( K, 1, DESCW, NPROW, NPCOL, MYROW, MYCOL,
497 $ IROW, ICOL, ITMP1, ITMP2 )
498.EQ..AND..EQ. IF( MYROWITMP1 MYCOLITMP2 ) THEN
499 WORK( IROW ) = -CONJG( WORK( IROW ) )
500 END IF
501 80 CONTINUE
502*
503* Solve the triangular system:
504* (T(KI+1:N,KI+1:N) - T(KI,KI))'*X = SCALE*WORK.
505*
506 DO 90 K = KI + 1, N
507 CALL INFOG2L( K, K, DESCT, NPROW, NPCOL, MYROW, MYCOL,
508 $ IROW, ICOL, ITMP1, ITMP2 )
509.EQ..AND..EQ. IF( ( MYROWITMP1 ) ( MYCOLITMP2 ) ) THEN
510 T( ( ICOL-1 )*LDT+IROW ) = T( ( ICOL-1 )*LDT+IROW ) -
511 $ SHIFT
512.LT. IF( CABS1( T( ( ICOL-1 )*LDT+IROW ) )SMIN( 1 ) )
513 $ T( ( ICOL-1 )*LDT+IROW ) = CMPLX( SMIN( 1 ) )
514 END IF
515 90 CONTINUE
516*
517.LT. IF( KIN ) THEN
518 CALL PCLATTRS( 'upper', 'conjugate transpose', 'nonunit',
519 $ 'y', N-KI, T, KI+1, KI+1, DESCT, WORK,
520 $ KI+1, 1, DESCW, SCALE, RWORK, INFO )
521 CALL INFOG2L( KI, 1, DESCW, NPROW, NPCOL, MYROW, MYCOL,
522 $ IROW, ICOL, ITMP1, ITMP2 )
523.EQ..AND..EQ. IF( MYROWITMP1 MYCOLITMP2 ) THEN
524 WORK( IROW ) = CMPLX( SCALE )
525 END IF
526 END IF
527*
528* Copy the vector x or Q*x to VL and normalize.
529*
530.NOT. IF( OVER ) THEN
531 CALL PCCOPY( N-KI+1, WORK, KI, 1, DESCW, 1, VL, KI, IS,
532 $ DESCVL, 1 )
533*
534 CALL PCAMAX( N-KI+1, REMAXC, II, VL, KI, IS, DESCVL, 1 )
535 REMAXD = ONE / MAX( CABS1( REMAXC ), UNFL )
536 CALL PCSSCAL( N-KI+1, REMAXD, VL, KI, IS, DESCVL, 1 )
537*
538 CALL PCLASET( ' ', KI-1, 1, CZERO, CZERO, VL, 1, IS,
539 $ DESCVL )
540 ELSE
541.LT. IF( KIN )
542 $ CALL PCGEMV( 'n', N, N-KI, CONE, VL, 1, KI+1, DESCVL,
543 $ WORK, KI+1, 1, DESCW, 1, CMPLX( SCALE ),
544 $ VL, 1, KI, DESCVL, 1 )
545*
546 CALL PCAMAX( N, REMAXC, II, VL, 1, KI, DESCVL, 1 )
547 REMAXD = ONE / MAX( CABS1( REMAXC ), UNFL )
548 CALL PCSSCAL( N, REMAXD, VL, 1, KI, DESCVL, 1 )
549 END IF
550*
551* Set back the original diagonal elements of T.
552*
553 DO 100 K = KI + 1, N
554 CALL INFOG2L( K, K, DESCT, NPROW, NPCOL, MYROW, MYCOL,
555 $ IROW, ICOL, ITMP1, ITMP2 )
556.EQ..AND..EQ. IF( ( MYROWITMP1 ) ( MYCOLITMP2 ) ) THEN
557 T( ( ICOL-1 )*LDT+IROW ) = WORK( LDW+IROW )
558 END IF
559 100 CONTINUE
560*
561 IS = IS + 1
562 110 CONTINUE
563 END IF
564*
565 RETURN
566*
567* End of PCTREVC
568*
cmplx
float cmplx[2]
Definition pblas.h:136
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:53
max
#define max(a, b)
Definition macros.h:21
pxerbla
subroutine pxerbla(contxt, srname, info)
Definition mpi.f:1600
infog2l
subroutine infog2l(grindx, gcindx, desc, nprow, npcol, myrow, mycol, lrindx, lcindx, rsrc, csrc)
Definition mpi.f:937
descinit
subroutine descinit(desc, m, n, mb, nb, irsrc, icsrc, ictxt, lld, info)
Definition mpi.f:777
blacs_gridinfo
subroutine blacs_gridinfo(cntxt, nprow, npcol, myrow, mycol)
Definition mpi.f:754
pclaset
subroutine pclaset(uplo, m, n, alpha, beta, a, ia, ja, desca)
Definition pcblastst.f:7508
pslamch
real function pslamch(ictxt, cmach)
Definition pcblastst.f:7455
pclattrs
subroutine pclattrs(uplo, trans, diag, normin, n, a, ia, ja, desca, x, ix, jx, descx, scale, cnorm, info)
Definition pclattrs.f:3
pctrevc
subroutine pctrevc(side, howmny, select, n, t, desct, vl, descvl, vr, descvr, mm, m, work, rwork, info)
Definition pctrevc.f:3
pslabad
subroutine pslabad(ictxt, small, large)
Definition pslabad.f:2