Routine: PZLABRD()  File: SRC\pzlabrd.f

 
 
# lines: 515
  # code: 515
  # comment: 0
  # blank:0
# Variables:45
# Callers:1
# Callings:2
# Words:228
# Keywords:107
 

 

 
..
     .. Array Arguments ..
     ..
  Purpose
  =======
  PZLABRD reduces the first NB rows and columns of a complex general
  M-by-N distributed matrix sub( A ) = A(IA:IA+M-1,JA:JA+N-1) to upper
  or lower bidiagonal form by an unitary transformation Q' * A * P, and
  returns the matrices X and Y which are needed to apply the transfor-
  mation to the unreduced part of sub( A ).
  If M >= N, sub( A ) is reduced to upper bidiagonal form; if M < N, to
  lower bidiagonal form.
  This is an auxiliary routine called by PZGEBRD.
  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
  Arguments
  =========
  M       (global input) INTEGER
          The number of rows to be operated on, i.e. the number of rows
          of the distributed submatrix sub( A ). M >= 0.
  N       (global input) INTEGER
          The number of columns to be operated on, i.e. the number of
          columns of the distributed submatrix sub( A ). N >= 0.
  NB      (global input) INTEGER
          The number of leading rows and columns of sub( A ) to be
          reduced.
  A       (local input/local output) COMPLEX*16 pointer into the
          local memory to an array of dimension (LLD_A,LOCc(JA+N-1)).
          On entry, this array contains the local pieces of the
          general distributed matrix sub( A ) to be reduced. On exit,
          the first NB rows and columns of the matrix are overwritten;
          the rest of the distributed matrix sub( A ) is unchanged.
          If m >= n, elements on and below the diagonal in the first NB
            columns, with the array TAUQ, represent the unitary
            matrix Q as a product of elementary reflectors; and
            elements above the diagonal in the first NB rows, with the
            array TAUP, represent the unitary matrix P as a product
            of elementary reflectors.
          If m < n, elements below the diagonal in the first NB
            columns, with the array TAUQ, represent the unitary
            matrix Q as a product of elementary reflectors, and
            elements on and above the diagonal in the first NB rows,
            with the array TAUP, represent the unitary matrix P as
            a product of elementary reflectors.
          See Further Details.
  IA      (global input) INTEGER
          The row index in the global array A indicating the first
          row of sub( A ).
  JA      (global input) INTEGER
          The column index in the global array A indicating the
          first column of sub( A ).
  DESCA   (global and local input) INTEGER array of dimension DLEN_.
          The array descriptor for the distributed matrix A.
  D       (local output) DOUBLE PRECISION array, dimension
          LOCr(IA+MIN(M,N)-1) if M >= N; LOCc(JA+MIN(M,N)-1) otherwise.
          The distributed diagonal elements of the bidiagonal matrix
          B: D(i) = A(ia+i-1,ja+i-1). D is tied to the distributed
          matrix A.
  E       (local output) DOUBLE PRECISION array, dimension
          LOCr(IA+MIN(M,N)-1) if M >= N; LOCc(JA+MIN(M,N)-2) otherwise.
          The distributed off-diagonal elements of the bidiagonal
          distributed matrix B:
          if m >= n, E(i) = A(ia+i-1,ja+i) for i = 1,2,...,n-1;
          if m < n, E(i) = A(ia+i,ja+i-1) for i = 1,2,...,m-1.
          E is tied to the distributed matrix A.
  TAUQ    (local output) COMPLEX*16 array dimension
          LOCc(JA+MIN(M,N)-1). The scalar factors of the elementary
          reflectors which represent the unitary matrix Q. TAUQ is
          tied to the distributed matrix A. See Further Details.
  TAUP    (local output) COMPLEX*16 array, dimension
          LOCr(IA+MIN(M,N)-1). The scalar factors of the elementary
          reflectors which represent the unitary matrix P. TAUP is
          tied to the distributed matrix A. See Further Details.
  X       (local output) COMPLEX*16 pointer into the local memory
          to an array of dimension (LLD_X,NB). On exit, the local
          pieces of the distributed M-by-NB matrix
          X(IX:IX+M-1,JX:JX+NB-1) required to update the unreduced
          part of sub( A ).
  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.
  Y       (local output) COMPLEX*16 pointer into the local memory
          to an array of dimension (LLD_Y,NB).  On exit, the local
          pieces of the distributed N-by-NB matrix
          Y(IY:IY+N-1,JY:JY+NB-1) required to update the unreduced
          part of sub( A ).
  IY      (global input) INTEGER
          The row index in the global array Y indicating the first
          row of sub( Y ).
  JY      (global input) INTEGER
          The column index in the global array Y indicating the
          first column of sub( Y ).
  DESCY   (global and local input) INTEGER array of dimension DLEN_.
          The array descriptor for the distributed matrix Y.
  WORK    (local workspace) COMPLEX*16 array, dimension (LWORK)
          LWORK >= NB_A + NQ, with
          NQ = NUMROC( N+MOD( IA-1, NB_Y ), NB_Y, MYCOL, IACOL, NPCOL )
          IACOL = INDXG2P( JA, NB_A, MYCOL, CSRC_A, NPCOL )
          INDXG2P and NUMROC are ScaLAPACK tool functions;
          MYROW, MYCOL, NPROW and NPCOL can be determined by calling
          the subroutine BLACS_GRIDINFO.
  Further Details
  ===============
  The matrices Q and P are represented as products of elementary
  reflectors:
     Q = H(1) H(2) . . . H(nb)  and  P = G(1) G(2) . . . G(nb)
  Each H(i) and G(i) has the form:
     H(i) = I - tauq * v * v'  and G(i) = I - taup * u * u'
  where tauq and taup are complex scalars, and v and u are complex
  vectors.
  If m >= n, v(1:i-1) = 0, v(i) = 1, and v(i:m) is stored on exit in
  A(ia+i-1:ia+m-1,ja+i-1); u(1:i) = 0, u(i+1) = 1, and u(i+1:n) is
  stored on exit in A(ia+i-1,ja+i:ja+n-1); tauq is stored in
  TAUQ(ja+i-1) and taup in TAUP(ia+i-1).
  If m < n, v(1:i) = 0, v(i+1) = 1, and v(i+1:m) is stored on exit in
  A(ia+i+1:ia+m-1,ja+i-1); u(1:i-1) = 0, u(i) = 1, and u(i:n) is
  stored on exit in A(ia+i-1,ja+i:ja+n-1); tauq is stored in
  TAUQ(ja+i-1) and taup in TAUP(ia+i-1).
  The elements of the vectors v and u together form the m-by-nb matrix
  V and the nb-by-n matrix U' which are needed, with X and Y, to apply
  the transformation to the unreduced part of the matrix, using a block
  update of the form:  sub( A ) := sub( A ) - V*Y' - X*U'.
  The contents of sub( A ) on exit are illustrated by the following
  examples with nb = 2:
  m = 6 and n = 5 (m > n):          m = 5 and n = 6 (m < n):
    (  1   1   u1  u1  u1 )           (  1   u1  u1  u1  u1  u1 )
    (  v1  1   1   u2  u2 )           (  1   1   u2  u2  u2  u2 )
    (  v1  v2  a   a   a  )           (  v1  1   a   a   a   a  )
    (  v1  v2  a   a   a  )           (  v1  v2  a   a   a   a  )
    (  v1  v2  a   a   a  )           (  v1  v2  a   a   a   a  )
    (  v1  v2  a   a   a  )
  where a denotes an element of the original matrix which is unchanged,
  vi denotes an element of the vector defining H(i), and ui an element
  of the vector defining G(i).
  =====================================================================
     .. Parameters ..
 
Display dynamic version Find AutoScroll Reload FontSize: - + Hide Comments Hide Blanks Frame FullScreen MailPrint

  
001        SUBROUTINE PZLABRD( M , N , NB , A , IA , JA , DESCA , D , E , TAUQ , TAUP ,
002       $X , IX , JX , DESCX , Y , IY , JY , DESCY , WORK )
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 IA , IX , IY , JA , JX , JY , M , N , NB
011        INTEGER BLOCK_CYCLIC_2D , CSRC_ , CTXT_ , DLEN_ , DTYPE_ ,
012       $LLD_ , MB_ , M_ , NB_ , N_ , RSRC_
013        PARAMETER( BLOCK_CYCLIC_2D = 1 , DLEN_ = 9 , DTYPE_ = 1 ,
014       $CTXT_ = 2 , M_ = 3 , N_ = 4 , MB_ = 5 , NB_ = 6 ,
015       $RSRC_ = 7 , CSRC_ = 8 , LLD_ = 9 )
016        COMPLEX*16 ONE , ZERO
017        PARAMETER( ONE =( 1.0D + 0 , 0.0D + 0 ) ,
018       $ZERO =( 0.0D + 0 , 0.0D + 0 ) )
019  *     ..
020  *     .. Local Scalars ..
021        INTEGER I , IACOL , IAROW , ICTXT , II , IPY , IW , J , JJ ,
022       $JWY , K , MYCOL , MYROW , NPCOL , NPROW
023        COMPLEX*16 ALPHA , TAU
024        INTEGER DESCD( DLEN_ ) , DESCE( DLEN_ ) ,
025       $DESCTP( DLEN_ ) , DESCTQ( DLEN_ ) ,
026       $DESCW( DLEN_ ) , DESCWY( DLEN_ )
027  *     ..
028  *     .. External Subroutines ..
029        EXTERNAL BLACS_GRIDINFO , DESCSET , INFOG2L , PDELSET ,
030       $PZCOPY , PZELGET , PZELSET , PZGEMV ,
031       $PZLACGV , PZLARFG , PZSCAL
032  *     ..
033  *     .. Intrinsic Functions ..
034        INTRINSIC DCMPLX , MIN , MOD
035  *     ..
036  *     .. Executable Statements ..
037  
038  *     Quick return if possible
039  
040        IF( M.LE.0 .OR. N.LE.0 )
041       $    RETURN
042  
043            ICTXT = DESCA( CTXT_ )
044            CALL BLACS_GRIDINFO( ICTXT , NPROW , NPCOL , MYROW , MYCOL )
045            CALL INFOG2L( IA , JA , DESCA , NPROW , NPCOL , MYROW , MYCOL , II , JJ ,
046       $    IAROW , IACOL )
047            IPY = DESCA( MB_ ) + 1
048            IW = MOD( IA - 1 , DESCA( NB_ ) ) + 1
049            ALPHA = ZERO
050  
051            CALL DESCSET( DESCWY , 1 , N + MOD( IA - 1 , DESCY( NB_ ) ) , 1 ,
052       $    DESCA( NB_ ) , IAROW , IACOL , ICTXT , 1 )
053            CALL DESCSET( DESCW , DESCA( MB_ ) , 1 , DESCA( MB_ ) , 1 , IAROW ,
054       $    IACOL , ICTXT , DESCA( MB_ ) )
055            CALL DESCSET( DESCTQ , 1 , JA + MIN(M , N) - 1 , 1 , DESCA( NB_ ) , IAROW ,
056       $    DESCA( CSRC_ ) , DESCA( CTXT_ ) , 1 )
057            CALL DESCSET( DESCTP , IA + MIN(M , N) - 1 , 1 , DESCA( MB_ ) , 1 ,
058       $    DESCA( RSRC_ ) , IACOL , DESCA( CTXT_ ) ,
059       $    DESCA( LLD_ ) )
060  
061            IF( M.GE.N ) THEN
062  
063  *             Reduce to upper bidiagonal form
064  
065                CALL DESCSET( DESCD , 1 , JA + MIN(M , N) - 1 , 1 , DESCA( NB_ ) , MYROW ,
066       $        DESCA( CSRC_ ) , DESCA( CTXT_ ) , 1 )
067                CALL DESCSET( DESCE , IA + MIN(M , N) - 1 , 1 , DESCA( MB_ ) , 1 ,
068       $        DESCA( RSRC_ ) , MYCOL , DESCA( CTXT_ ) ,
069       $        DESCA( LLD_ ) )
070                DO 10 K = 1 , NB
071                    I = IA + K - 1
072                    J = JA + K - 1
073                    JWY = IW + K
074  
075  *                 Update A(i : ia + m - 1 , j)
076  
077                    IF( K.GT.1 ) THEN
078                        CALL PZGEMV( 'No transpose' , M - K + 1 , K - 1 , - ONE , A , I , JA ,
079       $                DESCA , Y , IY , JY + K - 1 , DESCY , 1 , ONE , A , I ,
080       $                J , DESCA , 1 )
081                        CALL PZGEMV( 'No transpose' , M - K + 1 , K - 1 , - ONE , X , IX + K - 1 ,
082       $                JX , DESCX , A , IA , J , DESCA , 1 , ONE , A , I , J ,
083       $                DESCA , 1 )
084                        CALL PZELSET( A , I - 1 , J , DESCA , ALPHA )
085                    END IF
086  
087  *                 Generate reflection Q(i) to annihilate A(i + 1 : ia + m - 1 , j)
088  
089                    CALL PZLARFG ( M - K + 1 , ALPHA , I , J , A , I + 1 , J , DESCA , 1 ,
090       $            TAUQ )
091                    CALL PDELSET( D , 1 , J , DESCD , DBLE( ALPHA ) )
092                    CALL PZELSET( A , I , J , DESCA , ONE )
093  
094  *                 Compute Y(IA + I : IA + N - 1 , J)
095  
096                    CALL PZGEMV( 'Conjugate transpose' , M - K + 1 , N - K , ONE , A , I ,
097       $            J + 1 , DESCA , A , I , J , DESCA , 1 , ZERO ,
098       $            WORK( IPY ) , 1 , JWY , DESCWY , DESCWY( M_ ) )
099                    CALL PZGEMV( 'Conjugate transpose' , M - K + 1 , K - 1 , ONE , A , I ,
100       $            JA , DESCA , A , I , J , DESCA , 1 , ZERO , WORK , IW ,
101       $            1 , DESCW , 1 )
102                    CALL PZGEMV( 'Conjugate transpose' , K - 1 , N - K , - ONE , Y , IY ,
103       $            JY + K , DESCY , WORK , IW , 1 , DESCW , 1 , ONE ,
104       $            WORK( IPY ) , 1 , JWY , DESCWY , DESCWY( M_ ) )
105                    CALL PZGEMV( 'Conjugate transpose' , M - K + 1 , K - 1 , ONE , X ,
106       $            IX + K - 1 , JX , DESCX , A , I , J , DESCA , 1 , ZERO ,
107       $            WORK , IW , 1 , DESCW , 1 )
108                    CALL PZGEMV( 'Conjugate transpose' , K - 1 , N - K , - ONE , A , IA ,
109       $            J + 1 , DESCA , WORK , IW , 1 , DESCW , 1 , ONE ,
110       $            WORK( IPY ) , 1 , JWY , DESCWY , DESCWY( M_ ) )
111  
112                    CALL PZELGET( 'Rowwise' , ' ' , TAU , TAUQ , 1 , J , DESCTQ )
113                    CALL PZSCAL( N - K , TAU , WORK( IPY ) , 1 , JWY , DESCWY ,
114       $            DESCWY( M_ ) )
115                    CALL PZLACGV ( N - K , WORK( IPY ) , 1 , JWY , DESCWY ,
116       $            DESCWY( M_ ) )
117                    CALL PZCOPY( N - K , WORK( IPY ) , 1 , JWY , DESCWY , DESCWY( M_ ) ,
118       $            Y , IY + K - 1 , JY + K , DESCY , DESCY( M_ ) )
119  
120  *                 Update A(i , j + 1 : ja + n - 1)
121  
122                    CALL PZLACGV ( N - K , A , I , J + 1 , DESCA , DESCA( M_ ) )
123                    CALL PZLACGV ( K , A , I , JA , DESCA , DESCA( M_ ) )
124                    CALL PZGEMV( 'Conjugate transpose' , K , N - K , - ONE , Y , IY ,
125       $            JY + K , DESCY , A , I , JA , DESCA , DESCA( M_ ) , ONE ,
126       $            A , I , J + 1 , DESCA , DESCA( M_ ) )
127                    CALL PZLACGV ( K , A , I , JA , DESCA , DESCA( M_ ) )
128                    CALL PZLACGV ( K - 1 , X , IX + K - 1 , JX , DESCX , DESCX( M_ ) )
129                    CALL PZGEMV( 'Conjugate transpose' , K - 1 , N - K , - ONE , A , IA ,
130       $            J + 1 , DESCA , X , IX + K - 1 , JX , DESCX , DESCX( M_ ) ,
131       $            ONE , A , I , J + 1 , DESCA , DESCA( M_ ) )
132                    CALL PZLACGV ( K - 1 , X , IX + K - 1 , JX , DESCX , DESCX( M_ ) )
133                    CALL PZELSET( A , I , J , DESCA , DCMPLX( DBLE( ALPHA ) ) )
134  
135  *                 Generate reflection P(i) to annihilate A(i , j + 2 : ja + n - 1)
136  
137                    CALL PZLARFG ( N - K , ALPHA , I , J + 1 , A , I ,
138       $            MIN( J + 2 , N + JA - 1 ) , DESCA , DESCA( M_ ) , TAUP )
139                    CALL PDELSET( E , I , 1 , DESCE , DBLE( ALPHA ) )
140                    CALL PZELSET( A , I , J + 1 , DESCA , ONE )
141  
142  *                 Compute X(I + 1 : IA + M - 1 , J)
143  
144                    CALL PZGEMV( 'No transpose' , M - K , N - K , ONE , A , I + 1 , J + 1 ,
145       $            DESCA , A , I , J + 1 , DESCA , DESCA( M_ ) , ZERO , X ,
146       $            IX + K , JX + K - 1 , DESCX , 1 )
147                    CALL PZGEMV( 'No transpose' , K , N - K , ONE , Y , IY , JY + K ,
148       $            DESCY , A , I , J + 1 , DESCA , DESCA( M_ ) , ZERO ,
149       $            WORK , IW , 1 , DESCW , 1 )
150                    CALL PZGEMV( 'No transpose' , M - K , K , - ONE , A , I + 1 , JA ,
151       $            DESCA , WORK , IW , 1 , DESCW , 1 , ONE , X , IX + K ,
152       $            JX + K - 1 , DESCX , 1 )
153                    CALL PZGEMV( 'No transpose' , K - 1 , N - K , ONE , A , IA , J + 1 ,
154       $            DESCA , A , I , J + 1 , DESCA , DESCA( M_ ) , ZERO ,
155       $            WORK , IW , 1 , DESCW , 1 )
156                    CALL PZGEMV( 'No transpose' , M - K , K - 1 , - ONE , X , IX + K , JX ,
157       $            DESCX , WORK , IW , 1 , DESCW , 1 , ONE , X , IX + K ,
158       $            JX + K - 1 , DESCX , 1 )
159  
160                    CALL PZELGET( 'Columnwise' , ' ' , TAU , TAUP , I , 1 , DESCTP )
161                    CALL PZSCAL( M - K , TAU , X , IX + K , JX + K - 1 , DESCX , 1 )
162                    CALL PZLACGV ( N - K , A , I , J + 1 , DESCA , DESCA( M_ ) )
163     10         CONTINUE
164  
165            ELSE
166  
167  *             Reduce to lower bidiagonal form
168  
169                CALL DESCSET( DESCD , IA + MIN(M , N) - 1 , 1 , DESCA( MB_ ) , 1 ,
170       $        DESCA( RSRC_ ) , MYCOL , DESCA( CTXT_ ) ,
171       $        DESCA( LLD_ ) )
172                CALL DESCSET( DESCE , 1 , JA + MIN(M , N) - 1 , 1 , DESCA( NB_ ) , MYROW ,
173       $        DESCA( CSRC_ ) , DESCA( CTXT_ ) , 1 )
174                DO 20 K = 1 , NB
175                    I = IA + K - 1
176                    J = JA + K - 1
177                    JWY = IW + K
178  
179  *                 Update A(i , j : ja + n - 1)
180  
181                    CALL PZLACGV ( N - K + 1 , A , I , J , DESCA , DESCA( M_ ) )
182                    IF( K.GT.1 ) THEN
183                        CALL PZLACGV ( K - 1 , A , I , JA , DESCA , DESCA( M_ ) )
184                        CALL PZGEMV( 'Conjugate transpose' , K - 1 , N - K + 1 , - ONE , Y ,
185       $                IY , JY + K - 1 , DESCY , A , I , JA , DESCA ,
186       $                DESCA( M_ ) , ONE , A , I , J , DESCA ,
187       $                DESCA( M_ ) )
188                        CALL PZLACGV ( K - 1 , A , I , JA , DESCA , DESCA( M_ ) )
189                        CALL PZLACGV ( K - 1 , X , IX + K - 1 , JX , DESCX , DESCX( M_ ) )
190                        CALL PZGEMV( 'Conjugate transpose' , K - 1 , N - K + 1 , - ONE , A ,
191       $                IA , J , DESCA , X , IX + K - 1 , JX , DESCX ,
192       $                DESCX( M_ ) , ONE , A , I , J , DESCA ,
193       $                DESCA( M_ ) )
194                        CALL PZLACGV ( K - 1 , X , IX + K - 1 , JX , DESCX , DESCX( M_ ) )
195                        CALL PZELSET( A , I , J - 1 , DESCA , DCMPLX( DBLE( ALPHA ) ) )
196                    END IF
197  
198  *                 Generate reflection P(i) to annihilate A(i , j + 1 : ja + n - 1)
199  
200                    CALL PZLARFG ( N - K + 1 , ALPHA , I , J , A , I , J + 1 , DESCA ,
201       $            DESCA( M_ ) , TAUP )
202                    CALL PDELSET( D , I , 1 , DESCD , DBLE( ALPHA ) )
203                    CALL PZELSET( A , I , J , DESCA , ONE )
204  
205  *                 Compute X(i + 1 : ia + m - 1 , j)
206  
207                    CALL PZGEMV( 'No transpose' , M - K , N - K + 1 , ONE , A , I + 1 , J ,
208       $            DESCA , A , I , J , DESCA , DESCA( M_ ) , ZERO , X ,
209       $            IX + K , JX + K - 1 , DESCX , 1 )
210                    CALL PZGEMV( 'No transpose' , K - 1 , N - K + 1 , ONE , Y , IY , JY + K - 1 ,
211       $            DESCY , A , I , J , DESCA , DESCA( M_ ) , ZERO ,
212       $            WORK , IW , 1 , DESCW , 1 )
213                    CALL PZGEMV( 'No transpose' , M - K , K - 1 , - ONE , A , I + 1 , JA ,
214       $            DESCA , WORK , IW , 1 , DESCW , 1 , ONE , X , IX + K ,
215       $            JX + K - 1 , DESCX , 1 )
216                    CALL PZGEMV( 'No transpose' , K - 1 , N - K + 1 , ONE , A , IA , J ,
217       $            DESCA , A , I , J , DESCA , DESCA( M_ ) , ZERO ,
218       $            WORK , IW , 1 , DESCW , 1 )
219                    CALL PZGEMV( 'No transpose' , M - K , K - 1 , - ONE , X , IX + K , JX ,
220       $            DESCX , WORK , IW , 1 , DESCW , 1 , ONE , X , IX + K ,
221       $            JX + K - 1 , DESCX , 1 )
222  
223                    CALL PZELGET( 'Columnwise' , ' ' , TAU , TAUP , I , 1 , DESCTP )
224                    CALL PZSCAL( M - K , TAU , X , IX + K , JX + K - 1 , DESCX , 1 )
225                    CALL PZLACGV ( N - K + 1 , A , I , J , DESCA , DESCA( M_ ) )
226  
227  *                 Update A(i + 1 : ia + m - 1 , j)
228  
229                    CALL PZGEMV( 'No transpose' , M - K , K - 1 , - ONE , A , I + 1 , JA ,
230       $            DESCA , Y , IY , JY + K - 1 , DESCY , 1 , ONE , A , I + 1 , J ,
231       $            DESCA , 1 )
232                    CALL PZGEMV( 'No transpose' , M - K , K , - ONE , X , IX + K , JX ,
233       $            DESCX , A , IA , J , DESCA , 1 , ONE , A , I + 1 , J ,
234       $            DESCA , 1 )
235                    CALL PZELSET( A , I , J , DESCA , ALPHA )
236  
237  *                 Generate reflection Q(i) to annihilate A(i + 2 : ia + m - 1 , j)
238  
239                    CALL PZLARFG ( M - K , ALPHA , I + 1 , J , A , MIN( I + 2 , M + IA - 1 ) ,
240       $            J , DESCA , 1 , TAUQ )
241                    CALL PDELSET( E , 1 , J , DESCE , DBLE( ALPHA ) )
242                    CALL PZELSET( A , I + 1 , J , DESCA , ONE )
243  
244  *                 Compute Y(ia + i : ia + n - 1 , j)
245  
246                    CALL PZGEMV( 'Conjugate transpose' , M - K , N - K , ONE , A , I + 1 ,
247       $            J + 1 , DESCA , A , I + 1 , J , DESCA , 1 , ZERO ,
248       $            WORK( IPY ) , 1 , JWY , DESCWY , DESCWY( M_ ) )
249                    CALL PZGEMV( 'Conjugate transpose' , M - K , K - 1 , ONE , A , I + 1 ,
250       $            JA , DESCA , A , I + 1 , J , DESCA , 1 , ZERO , WORK , IW ,
251       $            1 , DESCW , 1 )
252                    CALL PZGEMV( 'Conjugate transpose' , K - 1 , N - K , - ONE , Y , IY ,
253       $            JY + K , DESCY , WORK , IW , 1 , DESCW , 1 , ONE ,
254       $            WORK( IPY ) , 1 , JWY , DESCWY , DESCWY( M_ ) )
255                    CALL PZGEMV( 'Conjugate transpose' , M - K , K , ONE , X , IX + K ,
256       $            JX , DESCX , A , I + 1 , J , DESCA , 1 , ZERO , WORK , IW ,
257       $            1 , DESCW , 1 )
258                    CALL PZGEMV( 'Conjugate transpose' , K , N - K , - ONE , A , IA ,
259       $            J + 1 , DESCA , WORK , IW , 1 , DESCW , 1 , ONE ,
260       $            WORK( IPY ) , 1 , JWY , DESCWY , DESCWY( M_ ) )
261  
262                    CALL PZELGET( 'Rowwise' , ' ' , TAU , TAUQ , 1 , J , DESCTQ )
263                    CALL PZSCAL( N - K , TAU , WORK( IPY ) , 1 , JWY , DESCWY ,
264       $            DESCWY( M_ ) )
265                    CALL PZLACGV ( N - K , WORK( IPY ) , 1 , JWY , DESCWY ,
266       $            DESCWY( M_ ) )
267                    CALL PZCOPY( N - K , WORK( IPY ) , 1 , JWY , DESCWY , DESCWY( M_ ) ,
268       $            Y , IY + K - 1 , JY + K , DESCY , DESCY( M_ ) )
269     20         CONTINUE
270            END IF
271  
272            RETURN
273  
274  *         End of PZLABRD
275  
276        END