Routine: PCLABRD()  File: SRC\pclabrd.f

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

 

 
..
     .. Array Arguments ..
     ..
  Purpose
  =======
  PCLABRD 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 PCGEBRD.
  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 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) REAL 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) REAL 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 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 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 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 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 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 PCLABRD( 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 ONE , ZERO
017        PARAMETER( ONE =( 1.0E + 0 , 0.0E + 0 ) ,
018       $ZERO =( 0.0E + 0 , 0.0E + 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 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 , PCCOPY ,
030       $PCELGET , PCELSET , PCGEMV , PCLACGV ,
031       $PCLARFG , PCSCAL , PSELSET
032  *     ..
033  *     .. Intrinsic Functions ..
034        INTRINSIC CMPLX , 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 PCGEMV( '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 PCGEMV( '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 PCELSET( 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 PCLARFG ( M - K + 1 , ALPHA , I , J , A , I + 1 , J , DESCA , 1 ,
090       $            TAUQ )
091                    CALL PSELSET( D , 1 , J , DESCD , REAL( ALPHA ) )
092                    CALL PCELSET( A , I , J , DESCA , ONE )
093  
094  *                 Compute Y(IA + I : IA + N - 1 , J)
095  
096                    CALL PCGEMV( '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 PCGEMV( '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 PCGEMV( '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 PCGEMV( '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 PCGEMV( '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 PCELGET( 'Rowwise' , ' ' , TAU , TAUQ , 1 , J , DESCTQ )
113                    CALL PCSCAL( N - K , TAU , WORK( IPY ) , 1 , JWY , DESCWY ,
114       $            DESCWY( M_ ) )
115                    CALL PCLACGV ( N - K , WORK( IPY ) , 1 , JWY , DESCWY ,
116       $            DESCWY( M_ ) )
117                    CALL PCCOPY( 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 PCLACGV ( N - K , A , I , J + 1 , DESCA , DESCA( M_ ) )
123                    CALL PCLACGV ( K , A , I , JA , DESCA , DESCA( M_ ) )
124                    CALL PCGEMV( '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 PCLACGV ( K , A , I , JA , DESCA , DESCA( M_ ) )
128                    CALL PCLACGV ( K - 1 , X , IX + K - 1 , JX , DESCX , DESCX( M_ ) )
129                    CALL PCGEMV( '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 PCLACGV ( K - 1 , X , IX + K - 1 , JX , DESCX , DESCX( M_ ) )
133                    CALL PCELSET( A , I , J , DESCA , CMPLX( REAL( ALPHA ) ) )
134  
135  *                 Generate reflection P(i) to annihilate A(i , j + 2 : ja + n - 1)
136  
137                    CALL PCLARFG ( N - K , ALPHA , I , J + 1 , A , I ,
138       $            MIN( J + 2 , N + JA - 1 ) , DESCA , DESCA( M_ ) , TAUP )
139                    CALL PSELSET( E , I , 1 , DESCE , REAL( ALPHA ) )
140                    CALL PCELSET( A , I , J + 1 , DESCA , ONE )
141  
142  *                 Compute X(I + 1 : IA + M - 1 , J)
143  
144                    CALL PCGEMV( '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 PCGEMV( '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 PCGEMV( '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 PCGEMV( '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 PCGEMV( '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 PCELGET( 'Columnwise' , ' ' , TAU , TAUP , I , 1 , DESCTP )
161                    CALL PCSCAL( M - K , TAU , X , IX + K , JX + K - 1 , DESCX , 1 )
162                    CALL PCLACGV ( 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 PCLACGV ( N - K + 1 , A , I , J , DESCA , DESCA( M_ ) )
182                    IF( K.GT.1 ) THEN
183                        CALL PCLACGV ( K - 1 , A , I , JA , DESCA , DESCA( M_ ) )
184                        CALL PCGEMV( '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 PCLACGV ( K - 1 , A , I , JA , DESCA , DESCA( M_ ) )
189                        CALL PCLACGV ( K - 1 , X , IX + K - 1 , JX , DESCX , DESCX( M_ ) )
190                        CALL PCGEMV( '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 PCLACGV ( K - 1 , X , IX + K - 1 , JX , DESCX , DESCX( M_ ) )
195                        CALL PCELSET( A , I , J - 1 , DESCA , CMPLX( REAL( ALPHA ) ) )
196                    END IF
197  
198  *                 Generate reflection P(i) to annihilate A(i , j + 1 : ja + n - 1)
199  
200                    CALL PCLARFG ( N - K + 1 , ALPHA , I , J , A , I , J + 1 , DESCA ,
201       $            DESCA( M_ ) , TAUP )
202                    CALL PSELSET( D , I , 1 , DESCD , REAL( ALPHA ) )
203                    CALL PCELSET( A , I , J , DESCA , ONE )
204  
205  *                 Compute X(i + 1 : ia + m - 1 , j)
206  
207                    CALL PCGEMV( '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 PCGEMV( '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 PCGEMV( '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 PCGEMV( '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 PCGEMV( '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 PCELGET( 'Columnwise' , ' ' , TAU , TAUP , I , 1 , DESCTP )
224                    CALL PCSCAL( M - K , TAU , X , IX + K , JX + K - 1 , DESCX , 1 )
225                    CALL PCLACGV ( N - K + 1 , A , I , J , DESCA , DESCA( M_ ) )
226  
227  *                 Update A(i + 1 : ia + m - 1 , j)
228  
229                    CALL PCGEMV( '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 PCGEMV( '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 PCELSET( A , I , J , DESCA , ALPHA )
236  
237  *                 Generate reflection Q(i) to annihilate A(i + 2 : ia + m - 1 , j)
238  
239                    CALL PCLARFG ( M - K , ALPHA , I + 1 , J , A , MIN( I + 2 , M + IA - 1 ) ,
240       $            J , DESCA , 1 , TAUQ )
241                    CALL PSELSET( E , 1 , J , DESCE , REAL( ALPHA ) )
242                    CALL PCELSET( A , I + 1 , J , DESCA , ONE )
243  
244  *                 Compute Y(ia + i : ia + n - 1 , j)
245  
246                    CALL PCGEMV( '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 PCGEMV( '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 PCGEMV( '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 PCGEMV( '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 PCGEMV( '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 PCELGET( 'Rowwise' , ' ' , TAU , TAUQ , 1 , J , DESCTQ )
263                    CALL PCSCAL( N - K , TAU , WORK( IPY ) , 1 , JWY , DESCWY ,
264       $            DESCWY( M_ ) )
265                    CALL PCLACGV ( N - K , WORK( IPY ) , 1 , JWY , DESCWY ,
266       $            DESCWY( M_ ) )
267                    CALL PCCOPY( 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 PCLABRD
275  
276        END