ScaLAPACK 2.1  2.1
ScaLAPACK: Scalable Linear Algebra PACKage
PB_Cpsyr2.c
Go to the documentation of this file.
1 /* ---------------------------------------------------------------------
2 *
3 * -- PBLAS auxiliary routine (version 2.0) --
4 * University of Tennessee, Knoxville, Oak Ridge National Laboratory,
5 * and University of California, Berkeley.
6 * April 1, 1998
7 *
8 * ---------------------------------------------------------------------
9 */
10 /*
11 * Include files
12 */
13 #include "../pblas.h"
14 #include "../PBpblas.h"
15 #include "../PBtools.h"
16 #include "../PBblacs.h"
17 #include "../PBblas.h"
18 
19 #ifdef __STDC__
20 void PB_Cpsyr2( PBTYP_T * TYPE, char * UPLO, int N, int K,
21  char * ALPHA, char * XC, int LDXC, char * XR, int LDXR,
22  char * YC, int LDYC, char * YR, int LDYR, char * A,
23  int IA, int JA, int * DESCA, TZSYR2_T SYR2 )
24 #else
25 void PB_Cpsyr2( TYPE, UPLO, N, K, ALPHA, XC, LDXC, XR, LDXR,
26  YC, LDYC, YR, LDYR, A, IA, JA, DESCA, SYR2 )
27 /*
28 * .. Scalar Arguments ..
29 */
30  char * UPLO;
31  int IA, JA, K, LDXC, LDXR, LDYC, LDYR, N;
32  char * ALPHA;
33  PBTYP_T * TYPE;
34  TZSYR2_T SYR2;
35 /*
36 * .. Array Arguments ..
37 */
38  int * DESCA;
39  char * A, * XC, * XR, * YC, * YR;
40 #endif
41 {
42 /*
43 * Purpose
44 * =======
45 *
46 * PB_Cpsyr2 performs a symmetric or Hermitian rank-2 update of the sub-
47 * matrix sub( A ) denoting A( IA:IA+N-1, JA:JA+N-1 ).
48 *
49 * Notes
50 * =====
51 *
52 * A description vector is associated with each 2D block-cyclicly dis-
53 * tributed matrix. This vector stores the information required to
54 * establish the mapping between a matrix entry and its corresponding
55 * process and memory location.
56 *
57 * In the following comments, the character _ should be read as
58 * "of the distributed matrix". Let A be a generic term for any 2D
59 * block cyclicly distributed matrix. Its description vector is DESC_A:
60 *
61 * NOTATION STORED IN EXPLANATION
62 * ---------------- --------------- ------------------------------------
63 * DTYPE_A (global) DESCA[ DTYPE_ ] The descriptor type.
64 * CTXT_A (global) DESCA[ CTXT_ ] The BLACS context handle, indicating
65 * the NPROW x NPCOL BLACS process grid
66 * A is distributed over. The context
67 * itself is global, but the handle
68 * (the integer value) may vary.
69 * M_A (global) DESCA[ M_ ] The number of rows in the distribu-
70 * ted matrix A, M_A >= 0.
71 * N_A (global) DESCA[ N_ ] The number of columns in the distri-
72 * buted matrix A, N_A >= 0.
73 * IMB_A (global) DESCA[ IMB_ ] The number of rows of the upper left
74 * block of the matrix A, IMB_A > 0.
75 * INB_A (global) DESCA[ INB_ ] The number of columns of the upper
76 * left block of the matrix A,
77 * INB_A > 0.
78 * MB_A (global) DESCA[ MB_ ] The blocking factor used to distri-
79 * bute the last M_A-IMB_A rows of A,
80 * MB_A > 0.
81 * NB_A (global) DESCA[ NB_ ] The blocking factor used to distri-
82 * bute the last N_A-INB_A columns of
83 * A, NB_A > 0.
84 * RSRC_A (global) DESCA[ RSRC_ ] The process row over which the first
85 * row of the matrix A is distributed,
86 * NPROW > RSRC_A >= 0.
87 * CSRC_A (global) DESCA[ CSRC_ ] The process column over which the
88 * first column of A is distributed.
89 * NPCOL > CSRC_A >= 0.
90 * LLD_A (local) DESCA[ LLD_ ] The leading dimension of the local
91 * array storing the local blocks of
92 * the distributed matrix A,
93 * IF( Lc( 1, N_A ) > 0 )
94 * LLD_A >= MAX( 1, Lr( 1, M_A ) )
95 * ELSE
96 * LLD_A >= 1.
97 *
98 * Let K be the number of rows of a matrix A starting at the global in-
99 * dex IA,i.e, A( IA:IA+K-1, : ). Lr( IA, K ) denotes the number of rows
100 * that the process of row coordinate MYROW ( 0 <= MYROW < NPROW ) would
101 * receive if these K rows were distributed over NPROW processes. If K
102 * is the number of columns of a matrix A starting at the global index
103 * JA, i.e, A( :, JA:JA+K-1, : ), Lc( JA, K ) denotes the number of co-
104 * lumns that the process MYCOL ( 0 <= MYCOL < NPCOL ) would receive if
105 * these K columns were distributed over NPCOL processes.
106 *
107 * The values of Lr() and Lc() may be determined via a call to the func-
108 * tion PB_Cnumroc:
109 * Lr( IA, K ) = PB_Cnumroc( K, IA, IMB_A, MB_A, MYROW, RSRC_A, NPROW )
110 * Lc( JA, K ) = PB_Cnumroc( K, JA, INB_A, NB_A, MYCOL, CSRC_A, NPCOL )
111 *
112 * Arguments
113 * =========
114 *
115 * TYPE (local input) pointer to a PBTYP_T structure
116 * On entry, TYPE is a pointer to a structure of type PBTYP_T,
117 * that contains type information (See pblas.h).
118 *
119 * UPLO (global input) pointer to CHAR
120 * On entry, UPLO specifies whether the local pieces of
121 * the array A containing the upper or lower triangular part
122 * of the symmetric or Hermitian submatrix sub( A ) are to be
123 * referenced as follows:
124 *
125 * UPLO = 'U' or 'u' Only the local pieces corresponding to
126 * the upper triangular part of the sym-
127 * metric or Hermitian submatrix sub( A )
128 * are to be referenced,
129 *
130 * UPLO = 'L' or 'l' Only the local pieces corresponding to
131 * the lower triangular part of the sym-
132 * metric or Hermitian submatrix sub( A )
133 * are to be referenced.
134 *
135 * N (global input) INTEGER
136 * On entry, N specifies the order of the submatrix sub( A ).
137 * N must be at least zero.
138 *
139 * K (global input) INTEGER
140 * On entry, K specifies the local number of columns of the lo-
141 * cal array XC and the local number of rows of the local array
142 * XR. K mut be at least zero.
143 *
144 * ALPHA (global input) pointer to CHAR
145 * On entry, ALPHA specifies the scalar alpha.
146 *
147 * XC (local input) pointer to CHAR
148 * On entry, XC is an array of dimension (LDXC,K). Before entry,
149 * this array contains the local entries of the matrix XC.
150 *
151 * LDXC (local input) INTEGER
152 * On entry, LDXC specifies the leading dimension of the array
153 * XC. LDXC must be at least max( 1, Lp( IA, N ) ).
154 *
155 * YC (local input) pointer to CHAR
156 * On entry, YC is an array of dimension (LDYC,K). Before entry,
157 * this array contains the local entries of the matrix YC.
158 *
159 * LDYC (local input) INTEGER
160 * On entry, LDYC specifies the leading dimension of the array
161 * YC. LDYC must be at least max( 1, Lp( IA, N ) ).
162 *
163 * XR (local input) pointer to CHAR
164 * On entry, XR is an array of dimension (LDXR,Kx), where Kx is
165 * at least Lc( JA, N ). Before entry, this array contains the
166 * local entries of the matrix XR.
167 *
168 * LDXR (local input) INTEGER
169 * On entry, LDXR specifies the leading dimension of the array
170 * XR. LDXR must be at least max( 1, K ).
171 *
172 * YR (local input) pointer to CHAR
173 * On entry, YR is an array of dimension (LDYR,Ky), where Ky is
174 * at least Lc( JA, N ). Before entry, this array contains the
175 * local entries of the matrix YR.
176 *
177 * LDYR (local input) INTEGER
178 * On entry, LDYR specifies the leading dimension of the array
179 * YR. LDYR must be at least max( 1, K ).
180 *
181 * A (local input/local output) pointer to CHAR
182 * On entry, A is an array of dimension (LLD_A, Ka), where Ka is
183 * at least Lc( 1, JA+N-1 ). Before entry, this array contains
184 * the local entries of the matrix A.
185 * Before entry with UPLO = 'U' or 'u', this array contains
186 * the local entries corresponding to the upper triangular part
187 * of the @(syhec) submatrix sub( A ), and the local entries
188 * corresponding to the strictly lower triangular of sub( A )
189 * are not referenced. On exit, the upper triangular part of
190 * sub( A ) is overwritten by the upper triangular part of the
191 * updated submatrix.
192 * Before entry with UPLO = 'L' or 'l', this array contains
193 * the local entries corresponding to the lower triangular part
194 * of the @(syhec) submatrix sub( A ), and the local entries
195 * corresponding to the strictly upper triangular of sub( A )
196 * are not referenced. On exit, the lower triangular part of
197 * sub( A ) is overwritten by the lower triangular part of the
198 * updated submatrix.
199 *
200 * IA (global input) INTEGER
201 * On entry, IA specifies A's global row index, which points to
202 * the beginning of the submatrix sub( A ).
203 *
204 * JA (global input) INTEGER
205 * On entry, JA specifies A's global column index, which points
206 * to the beginning of the submatrix sub( A ).
207 *
208 * DESCA (global and local input) INTEGER array
209 * On entry, DESCA is an integer array of dimension DLEN_. This
210 * is the array descriptor for the matrix A.
211 *
212 * SYR2 (local input) pointer to function of type TZSYR2_T
213 * On entry, SYR2 specifies the function performing the update
214 * of a single block.
215 *
216 * -- Written on April 1, 1998 by
217 * Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
218 *
219 * ---------------------------------------------------------------------
220 */
221 /*
222 * .. Local Scalars ..
223 */
224  int Acol, Arow, Aii, Aimb1, Ainb1, Ajj, Ald, Amp, Amb, Anb, Anq,
225  Aoffi, Aoffj, Arcol, Arrow, GoEast, GoSouth, IsColRepl,
226  IsRowRepl, XCinc, XRinc, Xii=0, Xjj=0, Xoffi=-1, Xoffj=-1,
227  YCinc, YRinc, iimax, ilow, imbloc, inbloc, ioffd, ioffx, iupp,
228  jjmax, joffd, joffx, lcmt, lcmt00, lmbloc, lnbloc, low, lower,
229  m1, mbloc, mblkd, mblks, mycol, myrow, n1, nbloc, nblkd,
230  nblks, npcol, nprow, pmb, qnb, size, tmp1, upp, upper;
231 /* ..
232 * .. Executable Statements ..
233 *
234 */
235 /*
236 * Quick return if possible
237 */
238  if( N <= 0 ) return;
239 /*
240 * Retrieve process grid information
241 */
242  Cblacs_gridinfo( DESCA[CTXT_], &nprow, &npcol, &myrow, &mycol );
243 /*
244 * Retrieve sub( A )'s local information: Aii, Ajj, Arow, Acol ...
245 */
246  PB_Cainfog2l( N, N, IA, JA, DESCA, nprow, npcol, myrow, mycol, &Aimb1,
247  &Ainb1, &Amp, &Anq, &Aii, &Ajj, &Arow, &Acol, &Arrow, &Arcol );
248 /*
249 * Quick return if I don't own any of sub( A ) or if sub( A ) is replicated in
250 * all processes.
251 */
252  if( ( Amp <= 0 ) || ( Anq <= 0 ) ) return;
253 
254  IsRowRepl = ( ( Arow < 0 ) || ( nprow == 1 ) );
255  IsColRepl = ( ( Acol < 0 ) || ( npcol == 1 ) );
256  Amb = DESCA[ MB_ ]; Anb = DESCA[ NB_ ]; Ald = DESCA[LLD_];
257  size = TYPE->size;
258 
259  if( IsRowRepl && IsColRepl )
260  {
261  SYR2( TYPE, UPLO, Amp, Anq, K, 0, ALPHA, XC, LDXC, YC, LDYC, XR, LDXR,
262  YR, LDYR, Mptr( A, Aii, Ajj, Ald, size ), Ald );
263  return;
264  }
265 
266  XCinc = size; XRinc = LDXR * size;
267  YCinc = size; YRinc = LDYR * size;
268  upper = ( Mupcase( UPLO[0] ) == CUPPER );
269  lower = ( Mupcase( UPLO[0] ) == CLOWER );
270 /*
271 * Initialize lcmt00, mblks, nblks, imbloc, inbloc, lmbloc, lnbloc, ilow, low,
272 * iupp, and upp.
273 */
274  PB_Cbinfo( 0, Amp, Anq, Aimb1, Ainb1, Amb, Anb, Arrow, Arcol, &lcmt00,
275  &mblks, &nblks, &imbloc, &inbloc, &lmbloc, &lnbloc, &ilow, &low,
276  &iupp, &upp );
277 
278  iimax = ( Aoffi = Aii - 1 ) + ( m1 = Amp );
279  jjmax = ( Aoffj = Ajj - 1 ) + ( n1 = Anq );
280  pmb = ( IsRowRepl ? Amb : nprow * Amb );
281  qnb = ( IsColRepl ? Anb : npcol * Anb );
282 /*
283 * Handle separately the first row and/or column of the LCM table. Update the
284 * LCM value of the curent block lcmt00, as well as the number of rows and
285 * columns mblks and nblks remaining in the LCM table.
286 */
287  GoSouth = ( lcmt00 > iupp );
288  GoEast = ( lcmt00 < ilow );
289 /*
290 * Go through the table looking for blocks owning diagonal entries.
291 */
292  if( ( !( GoSouth ) ) && ( !( GoEast ) ) )
293  {
294 /*
295 * The upper left block owns diagonal entries lcmt00 >= ilow && lcmt00 <= iupp
296 */
297  SYR2( TYPE, UPLO, imbloc, inbloc, K, lcmt00, ALPHA,
298  XC+Xii*XCinc, LDXC, YC+Xii*YCinc, LDYC,
299  XR+Xjj*XRinc, LDXR, YR+Xjj*YRinc, LDYR,
300  Mptr( A, Aii, Ajj, Ald, size ), Ald );
301 /*
302 * Decide whether one should go south or east in the table: Go east if
303 * the block below the current one only owns lower entries. If this block,
304 * however, owns diagonals, then go south.
305 */
306  GoSouth = !( GoEast = ( ( lcmt00 - ( iupp - upp + pmb ) ) < ilow ) );
307 
308  if( GoSouth )
309  {
310 /*
311 * When the upper triangular part of sub( A ) should be updated and one is
312 * planning to go south in the table, it is neccessary to take care of the
313 * remaining columns of these imbloc rows immediately.
314 */
315  if( upper && ( Anq > inbloc ) )
316  {
317  tmp1 = Anq - inbloc;
318  SYR2( TYPE, ALL, imbloc, tmp1, K, 0, ALPHA,
319  XC+Xii*XCinc, LDXC, YC+Xii*YCinc, LDYC,
320  XR+(Xjj+inbloc)*XRinc, LDXR, YR+(Xjj+inbloc)*YRinc, LDYR,
321  Mptr( A, Aii, Ajj+inbloc, Ald, size ), Ald );
322  }
323  Aii += imbloc; Xii += imbloc; m1 -= imbloc;
324  }
325  else
326  {
327 /*
328 * When the lower triangular part of sub( A ) should be updated and one is
329 * planning to go east in the table, it is neccessary to take care of the
330 * remaining rows of these inbloc columns immediately.
331 */
332  if( lower && ( Amp > imbloc ) )
333  {
334  tmp1 = Amp - imbloc;
335  SYR2( TYPE, ALL, tmp1, inbloc, K, 0, ALPHA,
336  XC+(Xii+imbloc)*XCinc, LDXC, YC+(Xii+imbloc)*YCinc, LDYC,
337  XR+Xjj*XRinc, LDXR, YR+Xjj*YRinc, LDYR,
338  Mptr( A, Aii+imbloc, Ajj, Ald, size ), Ald );
339  }
340  Ajj += inbloc; Xjj += inbloc; n1 -= inbloc;
341  }
342  }
343 
344  if( GoSouth )
345  {
346 /*
347 * Go one step south in the LCM table. Adjust the current LCM value as well as
348 * the local row indexes in A and XC.
349 */
350  lcmt00 -= ( iupp - upp + pmb ); mblks--; Aoffi += imbloc; Xoffi += imbloc;
351 /*
352 * While there are blocks remaining that own upper entries, keep going south.
353 * Adjust the current LCM value as well as the local row indexes in A and XC.
354 */
355  while( ( mblks > 0 ) && ( lcmt00 > upp ) )
356  { lcmt00 -= pmb; mblks--; Aoffi += Amb; Xoffi += Amb; }
357 /*
358 * Update the upper triangular part of sub( A ) we just skipped when necessary.
359 */
360  tmp1 = MIN( Aoffi, iimax ) - Aii + 1;
361  if( upper && ( tmp1 > 0 ) )
362  {
363  SYR2( TYPE, ALL, tmp1, n1, K, 0, ALPHA,
364  XC+Xii*XCinc, LDXC, YC+Xii*YCinc, LDYC,
365  XR+(Xoffj+1)*XRinc, LDXR, YR+(Xoffj+1)*YRinc, LDYR,
366  Mptr( A, Aii, Aoffj+1, Ald, size ), Ald );
367  Aii += tmp1; Xii += tmp1; m1 -= tmp1;
368  }
369 /*
370 * Return if no more row in the LCM table.
371 */
372  if( mblks <= 0 ) return;
373 /*
374 * lcmt00 <= upp. The current block owns either diagonals or lower entries.
375 * Save the current position in the LCM table. After this column has been
376 * completely taken care of, re-start from this row and the next column of
377 * the LCM table.
378 */
379  lcmt = lcmt00; mblkd = mblks; ioffd = Aoffi; ioffx = Xoffi;
380 
381  mbloc = Amb;
382  while( ( mblkd > 0 ) && ( lcmt >= ilow ) )
383  {
384 /*
385 * A block owning diagonals lcmt00 >= ilow && lcmt00 <= upp has been found.
386 */
387  if( mblkd == 1 ) mbloc = lmbloc;
388  SYR2( TYPE, UPLO, mbloc, inbloc, K, lcmt, ALPHA,
389  XC+(ioffx+1)*XCinc, LDXC, YC+(ioffx+1)*YCinc, LDYC,
390  XR+(Xoffj+1)*XRinc, LDXR, YR+(Xoffj+1)*YRinc, LDYR,
391  Mptr( A, ioffd+1, Aoffj+1, Ald, size ), Ald );
392  lcmt00 = lcmt; lcmt -= pmb;
393  mblks = mblkd; mblkd--;
394  Aoffi = ioffd; ioffd += mbloc;
395  Xoffi = ioffx; ioffx += mbloc;
396  }
397 /*
398 * Update the lower triangular part of sub( A ).
399 */
400  tmp1 = m1 - ioffd + Aii - 1;
401  if( lower && ( tmp1 > 0 ) )
402  SYR2( TYPE, ALL, tmp1, inbloc, K, 0, ALPHA,
403  XC+(ioffx+1)*XCinc, LDXC, YC+(ioffx+1)*YCinc, LDYC,
404  XR+(Xoffj+1)*XRinc, LDXR, YR+(Xoffj+1)*YRinc, LDYR,
405  Mptr( A, ioffd+1, Aoffj+1, Ald, size ), Ald );
406 
407  tmp1 = Aoffi - Aii + 1;
408  m1 -= tmp1;
409  n1 -= inbloc;
410  lcmt00 += low - ilow + qnb;
411  nblks--;
412  Aoffj += inbloc;
413  Xoffj += inbloc;
414 /*
415 * Update the upper triangular part of sub( A ).
416 */
417  if( upper && ( tmp1 > 0 ) && ( n1 > 0 ) )
418  SYR2( TYPE, ALL, tmp1, n1, K, 0, ALPHA,
419  XC+Xii*XCinc, LDXC, YC+Xii*YCinc, LDYC,
420  XR+(Xoffj+1)*XRinc, LDXR, YR+(Xoffj+1)*YRinc, LDYR,
421  Mptr( A, Aii, Aoffj+1, Ald, size ), Ald );
422  Aii = Aoffi + 1; Ajj = Aoffj + 1;
423  Xii = Xoffi + 1; Xjj = Xoffj + 1;
424  }
425  else if( GoEast )
426  {
427 /*
428 * Go one step east in the LCM table. Adjust the current LCM value as well as
429 * the local column index in A and XR.
430 */
431  lcmt00 += low - ilow + qnb; nblks--; Aoffj += inbloc; Xoffj += inbloc;
432 /*
433 * While there are blocks remaining that own lower entries, keep going east.
434 * Adjust the current LCM value as well as the local column index in A and XR.
435 */
436  while( ( nblks > 0 ) && ( lcmt00 < low ) )
437  { lcmt00 += qnb; nblks--; Aoffj += Anb; Xoffj += Anb; }
438 /*
439 * Update the lower triangular part of sub( A ).
440 */
441  tmp1 = MIN( Aoffj, jjmax ) - Ajj + 1;
442  if( lower && ( tmp1 > 0 ) )
443  {
444  SYR2( TYPE, ALL, m1, tmp1, K, 0, ALPHA,
445  XC+Xii*XCinc, LDXC, YC+Xii*YCinc, LDYC,
446  XR+Xjj*XRinc, LDXR, YR+Xjj*YRinc, LDYR,
447  Mptr( A, Aii, Ajj, Ald, size ), Ald );
448  Ajj += tmp1; Xjj += tmp1; n1 -= tmp1;
449  }
450 /*
451 * Return if no more column in the LCM table.
452 */
453  if( nblks <= 0 ) return;
454 /*
455 * lcmt00 >= low. The current block owns either diagonals or upper entries.
456 * Save the current position in the LCM table. After this row has been
457 * completely taken care of, re-start from this column and the next row of
458 * the LCM table.
459 */
460  lcmt = lcmt00; nblkd = nblks; joffd = Aoffj; joffx = Xoffj;
461 
462  nbloc = Anb;
463  while( ( nblkd > 0 ) && ( lcmt <= iupp ) )
464  {
465 /*
466 * A block owning diagonals lcmt00 >= low && lcmt00 <= iupp has been found.
467 */
468  if( nblkd == 1 ) nbloc = lnbloc;
469  SYR2( TYPE, UPLO, imbloc, nbloc, K, lcmt, ALPHA,
470  XC+Xii*XCinc, LDXC, YC+Xii*YCinc, LDYC,
471  XR+(joffx+1)*XRinc, LDXR, YR+(joffx+1)*YRinc, LDYR,
472  Mptr( A, Aii, joffd+1, Ald, size ), Ald );
473  lcmt00 = lcmt; lcmt += qnb;
474  nblks = nblkd; nblkd--;
475  Aoffj = joffd; joffd += nbloc;
476  Xoffj = joffx; joffx += nbloc;
477  }
478 /*
479 * Update the upper triangular part of sub( A ).
480 */
481  tmp1 = n1 - joffd + Ajj - 1;
482  if( upper && ( tmp1 > 0 ) )
483  SYR2( TYPE, ALL, imbloc, tmp1, K, 0, ALPHA,
484  XC+Xii*XCinc, LDXC, YC+Xii*YCinc, LDYC,
485  XR+(joffx+1)*XRinc, LDXR, YR+(joffx+1)*YRinc, LDYR,
486  Mptr( A, Aii, joffd+1, Ald, size ), Ald );
487 
488  tmp1 = Aoffj - Ajj + 1;
489  m1 -= imbloc;
490  n1 -= tmp1;
491  lcmt00 -= ( iupp - upp + pmb );
492  mblks--;
493  Aoffi += imbloc;
494  Xoffi += imbloc;
495 /*
496 * Update the lower triangular part of sub( A ).
497 */
498  if( lower && ( m1 > 0 ) && ( tmp1 > 0 ) )
499  SYR2( TYPE, ALL, m1, tmp1, K, 0, ALPHA,
500  XC+(Xoffi+1)*XCinc, LDXC, YC+(Xoffi+1)*YCinc, LDYC,
501  XR+Xjj*XRinc, LDXR, YR+Xjj*YRinc, LDYR,
502  Mptr( A, Aoffi+1, Ajj, Ald, size ), Ald );
503  Aii = Aoffi + 1; Ajj = Aoffj + 1;
504  Xii = Xoffi + 1; Xjj = Xoffj + 1;
505  }
506 /*
507 * Loop over the remaining columns of the LCM table.
508 */
509  nbloc = Anb;
510  while( nblks > 0 )
511  {
512  if( nblks == 1 ) nbloc = lnbloc;
513 /*
514 * While there are blocks remaining that own upper entries, keep going south.
515 * Adjust the current LCM value as well as the local row index in A and XC.
516 */
517  while( ( mblks > 0 ) && ( lcmt00 > upp ) )
518  { lcmt00 -= pmb; mblks--; Aoffi += Amb; Xoffi += Amb; }
519 /*
520 * Update the upper triangular part of sub( A ).
521 */
522  tmp1 = MIN( Aoffi, iimax ) - Aii + 1;
523  if( upper && ( tmp1 > 0 ) )
524  {
525  SYR2( TYPE, ALL, tmp1, n1, K, 0, ALPHA,
526  XC+Xii*XCinc, LDXC, YC+Xii*YCinc, LDYC,
527  XR+(Xoffj+1)*XRinc, LDXR, YR+(Xoffj+1)*YRinc, LDYR,
528  Mptr( A, Aii, Aoffj+1, Ald, size ), Ald );
529  Aii += tmp1;
530  Xii += tmp1;
531  m1 -= tmp1;
532  }
533 /*
534 * Return if no more row in the LCM table.
535 */
536  if( mblks <= 0 ) return;
537 /*
538 * lcmt00 <= upp. The current block owns either diagonals or lower entries.
539 * Save the current position in the LCM table. After this column has been
540 * completely taken care of, re-start from this row and the next column of
541 * the LCM table.
542 */
543  lcmt = lcmt00; mblkd = mblks; ioffd = Aoffi; ioffx = Xoffi;
544 
545  mbloc = Amb;
546  while( ( mblkd > 0 ) && ( lcmt >= low ) )
547  {
548 /*
549 * A block owning diagonals lcmt00 >= low && lcmt00 <= upp has been found.
550 */
551  if( mblkd == 1 ) mbloc = lmbloc;
552  SYR2( TYPE, UPLO, mbloc, nbloc, K, lcmt, ALPHA,
553  XC+(ioffx+1)*XCinc, LDXC, YC+(ioffx+1)*YCinc, LDYC,
554  XR+(Xoffj+1)*XRinc, LDXR, YR+(Xoffj+1)*YRinc, LDYR,
555  Mptr( A, ioffd+1, Aoffj+1, Ald, size ), Ald );
556  lcmt00 = lcmt; lcmt -= pmb;
557  mblks = mblkd; mblkd--;
558  Aoffi = ioffd; Xoffi = ioffx;
559  ioffd += mbloc; ioffx += mbloc;
560  }
561 /*
562 * Update the lower triangular part of sub( A ).
563 */
564  tmp1 = m1 - ioffd + Aii - 1;
565  if( lower && ( tmp1 > 0 ) )
566  SYR2( TYPE, ALL, tmp1, nbloc, K, 0, ALPHA,
567  XC+(ioffx+1)*XCinc, LDXC, YC+(ioffx+1)*YCinc, LDYC,
568  XR+(Xoffj+1)*XRinc, LDXR, YR+(Xoffj+1)*YRinc, LDYR,
569  Mptr( A, ioffd+1, Aoffj+1, Ald, size ), Ald );
570 
571  tmp1 = MIN( Aoffi, iimax ) - Aii + 1;
572  m1 -= tmp1;
573  n1 -= nbloc;
574  lcmt00 += qnb;
575  nblks--;
576  Aoffj += nbloc;
577  Xoffj += nbloc;
578 /*
579 * Update the upper triangular part of sub( A ).
580 */
581  if( upper && ( tmp1 > 0 ) && ( n1 > 0 ) )
582  SYR2( TYPE, ALL, tmp1, n1, K, 0, ALPHA,
583  XC+Xii*XCinc, LDXC, YC+Xii*YCinc, LDYC,
584  XR+(Xoffj+1)*XRinc, LDXR, YR+(Xoffj+1)*YRinc, LDYR,
585  Mptr( A, Aii, Aoffj+1, Ald, size ), Ald );
586  Aii = Aoffi + 1; Ajj = Aoffj + 1;
587  Xii = Xoffi + 1; Xjj = Xoffj + 1;
588  }
589 /*
590 * End of PB_Cpsyr2
591 */
592 }
TYPE
#define TYPE
Definition: clamov.c:7
MB_
#define MB_
Definition: PBtools.h:43
NB_
#define NB_
Definition: PBtools.h:44
PB_Cainfog2l
void PB_Cainfog2l()
LLD_
#define LLD_
Definition: PBtools.h:47
PB_Cbinfo
void PB_Cbinfo()
CLOWER
#define CLOWER
Definition: PBblas.h:25
ALL
#define ALL
Definition: PBblas.h:50
MIN
#define MIN(a_, b_)
Definition: PBtools.h:76
PB_Cpsyr2
void PB_Cpsyr2(PBTYP_T *TYPE, char *UPLO, int N, int K, char *ALPHA, char *XC, int LDXC, char *XR, int LDXR, char *YC, int LDYC, char *YR, int LDYR, char *A, int IA, int JA, int *DESCA, TZSYR2_T SYR2)
Definition: PB_Cpsyr2.c:25
TZSYR2_T
void(* TZSYR2_T)()
Definition: pblas.h:426
Cblacs_gridinfo
void Cblacs_gridinfo()
PBTYP_T
Definition: pblas.h:325
Mupcase
#define Mupcase(C)
Definition: PBtools.h:83
CUPPER
#define CUPPER
Definition: PBblas.h:26
Mptr
#define Mptr(a_, i_, j_, lda_, siz_)
Definition: PBtools.h:132
CTXT_
#define CTXT_
Definition: PBtools.h:38