Routine PZLARFG() Documented By Universal Report

Routine: PZLARFG()

# lines:	291
# code:	291
# comment:	0
# blank:	0
# Variables:	42
# Callers:	11
# Callings:	0
# Words:	119
# Keywords:	67

..
     .. Array Arguments ..
     ..
  Purpose
  =======
  PZLARFG generates a complex elementary reflector H of order n, such
  that
     H * sub( X ) = H * ( x(iax,jax) ) = ( alpha ),   H' * H = I.
                        (      x     )   (   0   )
  where alpha is a real scalar, and sub( X ) is an (N-1)-element
  complex distributed vector X(IX:IX+N-2,JX) if INCX = 1 and
  X(IX,JX:JX+N-2) if INCX = DESCX(M_).  H is represented in the form
        H = I - tau * ( 1 ) * ( 1 v' ) ,
                      ( v )
  where tau is a complex scalar and v is a complex (N-1)-element
  vector. Note that H is not Hermitian.
  If the elements of sub( X ) are all zero and X(IAX,JAX) is real,
  then tau = 0 and H is taken to be the unit matrix.
  Otherwise  1 <= real(tau) <= 2 and abs(tau-1) <= 1.
  Notes
  =====
  Each global data object is described by an associated description
  vector.  This vector stores the information required to establish
  the mapping between an object element and its corresponding process
  and memory location.
  Let A be a generic term for any 2D block cyclicly distributed array.
  Such a global array has an associated description vector DESCA.
  In the following comments, the character _ should be read as
  "of the global array".
  NOTATION        STORED IN      EXPLANATION
  --------------- -------------- --------------------------------------
  DTYPE_A(global) DESCA( DTYPE_ )The descriptor type.  In this case,
                                 DTYPE_A = 1.
  CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
                                 the BLACS process grid A is distribu-
                                 ted over. The context itself is glo-
                                 bal, but the handle (the integer
                                 value) may vary.
  M_A    (global) DESCA( M_ )    The number of rows in the global
                                 array A.
  N_A    (global) DESCA( N_ )    The number of columns in the global
                                 array A.
  MB_A   (global) DESCA( MB_ )   The blocking factor used to distribute
                                 the rows of the array.
  NB_A   (global) DESCA( NB_ )   The blocking factor used to distribute
                                 the columns of the array.
  RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
                                 row of the array A is distributed.
  CSRC_A (global) DESCA( CSRC_ ) The process column over which the
                                 first column of the array A is
                                 distributed.
  LLD_A  (local)  DESCA( LLD_ )  The leading dimension of the local
                                 array.  LLD_A >= MAX(1,LOCr(M_A)).
  Let K be the number of rows or columns of a distributed matrix,
  and assume that its process grid has dimension p x q.
  LOCr( K ) denotes the number of elements of K that a process
  would receive if K were distributed over the p processes of its
  process column.
  Similarly, LOCc( K ) denotes the number of elements of K that a
  process would receive if K were distributed over the q processes of
  its process row.
  The values of LOCr() and LOCc() may be determined via a call to the
  ScaLAPACK tool function, NUMROC:
          LOCr( M ) = NUMROC( M, MB_A, MYROW, RSRC_A, NPROW ),
          LOCc( N ) = NUMROC( N, NB_A, MYCOL, CSRC_A, NPCOL ).
  An upper bound for these quantities may be computed by:
          LOCr( M ) <= ceil( ceil(M/MB_A)/NPROW )*MB_A
          LOCc( N ) <= ceil( ceil(N/NB_A)/NPCOL )*NB_A
  Because vectors may be viewed as a subclass of matrices, a
  distributed vector is considered to be a distributed matrix.
  Arguments
  =========
  N       (global input) INTEGER
          The global order of the elementary reflector. N >= 0.
  ALPHA   (local output) COMPLEX*16
          On exit, alpha is computed in the process scope having the
          vector sub( X ).
  IAX     (global input) INTEGER
          The global row index in X of X(IAX,JAX).
  JAX     (global input) INTEGER
          The global column index in X of X(IAX,JAX).
  X       (local input/local output) COMPLEX*16, pointer into the
          local memory to an array of dimension (LLD_X,*). This array
          contains the local pieces of the distributed vector sub( X ).
          Before entry, the incremented array sub( X ) must contain
          the vector x. On exit, it is overwritten with the vector v.
  IX      (global input) INTEGER
          The row index in the global array X indicating the first
          row of sub( X ).
  JX      (global input) INTEGER
          The column index in the global array X indicating the
          first column of sub( X ).
  DESCX   (global and local input) INTEGER array of dimension DLEN_.
          The array descriptor for the distributed matrix X.
  INCX    (global input) INTEGER
          The global increment for the elements of X. Only two values
          of INCX are supported in this version, namely 1 and M_X.
          INCX must not be zero.
  TAU     (local output) COMPLEX*16, array, dimension  LOCc(JX)
          if INCX = 1, and LOCr(IX) otherwise. This array contains the
          Householder scalars related to the Householder vectors.
          TAU is tied to the distributed matrix X.
  =====================================================================
     .. Parameters ..

Display dynamic version

FontSize: - +

001        SUBROUTINE PZLARFG( N , ALPHA , IAX , JAX , X , IX , JX , DESCX , INCX ,
002       $TAU )
003  
004  *     -- ScaLAPACK auxiliary routine(version 1.7) --
005  *     University of Tennessee , Knoxville , Oak Ridge National Laboratory ,
006  *     and University of California , Berkeley.
007  *     May 1 , 1997
008  
009  *     .. Scalar Arguments ..
010        INTEGER IAX , INCX , IX , JAX , JX , N
011        COMPLEX*16 ALPHA
012        INTEGER BLOCK_CYCLIC_2D , CSRC_ , CTXT_ , DLEN_ , DTYPE_ ,
013       $LLD_ , MB_ , M_ , NB_ , N_ , RSRC_
014        PARAMETER( BLOCK_CYCLIC_2D = 1 , DLEN_ = 9 , DTYPE_ = 1 ,
015       $CTXT_ = 2 , M_ = 3 , N_ = 4 , MB_ = 5 , NB_ = 6 ,
016       $RSRC_ = 7 , CSRC_ = 8 , LLD_ = 9 )
017        DOUBLE PRECISION ONE , ZERO
018        PARAMETER( ONE = 1.0D + 0 , ZERO = 0.0D + 0 )
019  *     ..
020  *     .. Local Scalars ..
021        INTEGER ICTXT , IIAX , INDXTAU , IXCOL , IXROW , J , JJAX ,
022       $KNT , MYCOL , MYROW , NPCOL , NPROW
023        DOUBLE PRECISION ALPHI , ALPHR , BETA , RSAFMN , SAFMIN , XNORM
024  *     ..
025  *     .. External Subroutines ..
026        EXTERNAL BLACS_GRIDINFO , INFOG2L , PDZNRM2 ,
027       $ZGEBR2D , ZGEBS2D , PZSCAL ,
028       $PZDSCAL
029  *     ..
030  *     .. External Functions ..
031        DOUBLE PRECISION DLAMCH , DLAPY3
032        COMPLEX*16 ZLADIV
033        EXTERNAL DLAMCH , DLAPY3 , ZLADIV
034  *     ..
035  *     .. Intrinsic Functions ..
036        INTRINSIC ABS , DBLE , DCMPLX , DIMAG , SIGN
037  *     ..
038  *     .. Executable Statements ..
039  
040  *     Get grid parameters.
041  
042        ICTXT = DESCX( CTXT_ )
043        CALL BLACS_GRIDINFO( ICTXT , NPROW , NPCOL , MYROW , MYCOL )
044  
045        IF( INCX.EQ.DESCX( M_ ) ) THEN
046  
047  *         sub( X ) is distributed across a process row.
048  
048               
049            CALL INFOG2L( IX , JAX , DESCX , NPROW , NPCOL , MYROW , MYCOL ,
050       $    IIAX , JJAX , IXROW , IXCOL )
051  
052            IF( MYROW.NE.IXROW )
052                   
053       $        RETURN
054  
055  *             Broadcast X(IAX , JAX) across the process row.
056  
057                IF( MYCOL.EQ.IXCOL ) THEN
057                       
058                    J = IIAX + (JJAX - 1)*DESCX( LLD_ )
059                    CALL ZGEBS2D( ICTXT , 'Rowwise' , ' ' , 1 , 1 , X( J ) , 1 )
060                    ALPHA = X( J )
061                ELSE
061                       
062                    CALL ZGEBR2D( ICTXT , 'Rowwise' , ' ' , 1 , 1 , ALPHA , 1 ,
063       $            MYROW , IXCOL )
064                END IF
065  
066                INDXTAU = IIAX
067  
068            ELSE
069  
070  *             sub( X ) is distributed across a process column.
071  
071                   
072                CALL INFOG2L( IAX , JX , DESCX , NPROW , NPCOL , MYROW , MYCOL ,
073       $        IIAX , JJAX , IXROW , IXCOL )
074  
075                IF( MYCOL.NE.IXCOL )
075                       
076       $            RETURN
077  
078  *                 Broadcast X(IAX , JAX) across the process column.
079  
080                    IF( MYROW.EQ.IXROW ) THEN
080                           
081                        J = IIAX + (JJAX - 1)*DESCX( LLD_ )
082                        CALL ZGEBS2D( ICTXT , 'Columnwise' , ' ' , 1 , 1 , X( J ) , 1 )
083                        ALPHA = X( J )
084                    ELSE
084                           
085                        CALL ZGEBR2D( ICTXT , 'Columnwise' , ' ' , 1 , 1 , ALPHA , 1 ,
086       $                IXROW , MYCOL )
087                    END IF
088  
089                    INDXTAU = JJAX
090  
091                END IF
092  
093                IF( N.LE.0 ) THEN
093                       
094                    TAU( INDXTAU ) = ZERO
095                    RETURN
096                END IF
097  
098                CALL PDZNRM2( N - 1 , XNORM , X , IX , JX , DESCX , INCX )
099                ALPHR = DBLE( ALPHA )
100                ALPHI = DIMAG( ALPHA )
101  
102                IF( XNORM.EQ.ZERO .AND. ALPHI.EQ.ZERO ) THEN
103  
104  *                 H = I
105  
105                       
106                    TAU( INDXTAU ) = ZERO
107  
108                ELSE
109  
110  *                 General case
111  
111                       
112                    BETA = - SIGN( DLAPY3( ALPHR , ALPHI , XNORM ) , ALPHR )
113                    SAFMIN = DLAMCH( 'S' )
114                    RSAFMN = ONE / SAFMIN
115                    IF( ABS( BETA ).LT.SAFMIN ) THEN
116  
117  *                     XNORM , BETA may be inaccurate ; scale X and recompute them
118  
118                           
119                        KNT = 0
120     10 CONTINUE
121        KNT = KNT + 1
122        CALL PZDSCAL( N - 1 , RSAFMN , X , IX , JX , DESCX , INCX )
123        BETA = BETA*RSAFMN
124        ALPHI = ALPHI*RSAFMN
125        ALPHR = ALPHR*RSAFMN
126        IF( ABS( BETA ).LT.SAFMIN )
126               
127       $    GO TO 10
128  
129  *         New BETA is at most 1 , at least SAFMIN
130  
131            CALL PDZNRM2( N - 1 , XNORM , X , IX , JX , DESCX , INCX )
132            ALPHA = DCMPLX( ALPHR , ALPHI )
133            BETA = - SIGN( DLAPY3( ALPHR , ALPHI , XNORM ) , ALPHR )
134            TAU( INDXTAU ) = DCMPLX(( BETA - ALPHR ) / BETA ,
135       $    - ALPHI / BETA )
136            ALPHA = ZLADIV( DCMPLX( ONE ) , ALPHA - BETA )
137            CALL PZSCAL( N - 1 , ALPHA , X , IX , JX , DESCX , INCX )
138  
139  *         If ALPHA is subnormal , it may lose relative accuracy
140  
141            ALPHA = BETA
142            DO 20 J = 1 , KNT
142                   
143                ALPHA = ALPHA*SAFMIN
144     20     CONTINUE
145        ELSE
145               
146            TAU( INDXTAU ) = DCMPLX(( BETA - ALPHR ) / BETA ,
147       $    - ALPHI / BETA )
148            ALPHA = ZLADIV( DCMPLX( ONE ) , ALPHA - BETA )
149            CALL PZSCAL( N - 1 , ALPHA , X , IX , JX , DESCX , INCX )
150            ALPHA = BETA
151        END IF
152        END IF
153  
154        RETURN
155  
156  *     End of PZLARFG
157  
158        END2636

Variables in Routine PZLARFG()

Summary Report
Data Type	Quantity	Size(byte)
*COMPLEX16**	2	?
DOUBLE PRECISION	10	40
INTEGER	29	116
REAL	1	4
TOTAL	42	160

List of Variables

COMPLEX*16

ALPHA ZLADIV

DOUBLE PRECISION

ALPHI ALPHR BETA DLAMCH DLAPY3
ONE RSAFMN SAFMIN XNORM ZERO

INTEGER

BLOCK_CYCLIC_2D CSRC_ CTXT_ DLEN_ DTYPE_
IAX ICTXT IIAX INCX INDXTAU
IX IXCOL IXROW J JAX
JJAX JX KNT LLD_ M_
MB_ MYCOL MYROW N N_
NB_ NPCOL NPROW RSRC_

REAL

TAU

Variables Dependence Graph
Put the mouse over a right hand side variable to display the corresponding line of the dependence
-
- -
ALPHA <--- ALPHA{2ALPHA = ALPHA*SAFMIN, 3}, ALPHIALPHA = DCMPLX( ALPHR, ALPHI ), JALPHA = X( J ){2ALPHA = X( J )}, ALPHRALPHA = DCMPLX( ALPHR, ALPHI ), ONE{2}, SAFMINALPHA = ALPHA*SAFMIN, BETA{2ALPHA = BETA, 3, 4ALPHA = BETA}, ZLADIV{2}

ALPHI <--- ALPHAALPHI = DIMAG( ALPHA ), ALPHIALPHI = ALPHI*RSAFMN, RSAFMNALPHI = ALPHI*RSAFMN

ALPHR <--- ALPHAALPHR = DBLE( ALPHA ), ALPHRALPHR = ALPHR*RSAFMN, RSAFMNALPHR = ALPHR*RSAFMN

BETA <--- ALPHI{2}, ALPHR{2}, RSAFMNBETA = BETA*RSAFMN, BETABETA = BETA*RSAFMN, XNORM{2}, DLAPY3{2}

ICTXT <--- CTXT_ICTXT = DESCX( CTXT_ )

INDXTAU <--- IIAXINDXTAU = IIAX, JJAXINDXTAU = JJAX

J <--- IIAXJ = IIAX+(JJAX-1)*DESCX( LLD_ ){2J = IIAX+(JJAX-1)*DESCX( LLD_ )}, JJAXJ = IIAX+(JJAX-1)*DESCX( LLD_ ){2J = IIAX+(JJAX-1)*DESCX( LLD_ )}, KNTDO 20 J = 1, KNT, LLD_J = IIAX+(JJAX-1)*DESCX( LLD_ ){2J = IIAX+(JJAX-1)*DESCX( LLD_ )}

KNT <--- KNTKNT = KNT + 1

RSAFMN <--- ONERSAFMN = ONE / SAFMIN, SAFMINRSAFMN = ONE / SAFMIN

SAFMIN <--- DLAMCHSAFMIN = DLAMCH( 'S' )

TAU <--- ALPHI{2}, ALPHR{2}, BETA{2}, ZEROTAU( INDXTAU ) = ZERO{2TAU( INDXTAU ) = ZERO}