ScaLAPACK 2.1  2.1
ScaLAPACK: Scalable Linear Algebra PACKage
pdlase2.f
Go to the documentation of this file.
1  SUBROUTINE pdlase2( UPLO, M, N, ALPHA, BETA, A, IA, JA, DESCA )
2 *
3 * -- ScaLAPACK auxiliary routine (version 1.7) --
4 * University of Tennessee, Knoxville, Oak Ridge National Laboratory,
5 * and University of California, Berkeley.
6 * May 1, 1997
7 *
8 * .. Scalar Arguments ..
9  CHARACTER uplo
10  INTEGER ia, ja, m, n
11  DOUBLE PRECISION alpha, beta
12 * ..
13 * .. Array Arguments ..
14  INTEGER desca( * )
15  DOUBLE PRECISION A( * )
16 * ..
17 *
18 * Purpose
19 * =======
20 *
21 * PDLASE2 initializes an M-by-N distributed matrix sub( A ) denoting
22 * A(IA:IA+M-1,JA:JA+N-1) to BETA on the diagonal and ALPHA on the
23 * offdiagonals. PDLASE2 requires that only dimension of the matrix
24 * operand is distributed.
25 *
26 * Notes
27 * =====
28 *
29 * Each global data object is described by an associated description
30 * vector. This vector stores the information required to establish
31 * the mapping between an object element and its corresponding process
32 * and memory location.
33 *
34 * Let A be a generic term for any 2D block cyclicly distributed array.
35 * Such a global array has an associated description vector DESCA.
36 * In the following comments, the character _ should be read as
37 * "of the global array".
38 *
39 * NOTATION STORED IN EXPLANATION
40 * --------------- -------------- --------------------------------------
41 * DTYPE_A(global) DESCA( DTYPE_ )The descriptor type. In this case,
42 * DTYPE_A = 1.
43 * CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
44 * the BLACS process grid A is distribu-
45 * ted over. The context itself is glo-
46 * bal, but the handle (the integer
47 * value) may vary.
48 * M_A (global) DESCA( M_ ) The number of rows in the global
49 * array A.
50 * N_A (global) DESCA( N_ ) The number of columns in the global
51 * array A.
52 * MB_A (global) DESCA( MB_ ) The blocking factor used to distribute
53 * the rows of the array.
54 * NB_A (global) DESCA( NB_ ) The blocking factor used to distribute
55 * the columns of the array.
56 * RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
57 * row of the array A is distributed.
58 * CSRC_A (global) DESCA( CSRC_ ) The process column over which the
59 * first column of the array A is
60 * distributed.
61 * LLD_A (local) DESCA( LLD_ ) The leading dimension of the local
62 * array. LLD_A >= MAX(1,LOCr(M_A)).
63 *
64 * Let K be the number of rows or columns of a distributed matrix,
65 * and assume that its process grid has dimension p x q.
66 * LOCr( K ) denotes the number of elements of K that a process
67 * would receive if K were distributed over the p processes of its
68 * process column.
69 * Similarly, LOCc( K ) denotes the number of elements of K that a
70 * process would receive if K were distributed over the q processes of
71 * its process row.
72 * The values of LOCr() and LOCc() may be determined via a call to the
73 * ScaLAPACK tool function, NUMROC:
74 * LOCr( M ) = NUMROC( M, MB_A, MYROW, RSRC_A, NPROW ),
75 * LOCc( N ) = NUMROC( N, NB_A, MYCOL, CSRC_A, NPCOL ).
76 * An upper bound for these quantities may be computed by:
77 * LOCr( M ) <= ceil( ceil(M/MB_A)/NPROW )*MB_A
78 * LOCc( N ) <= ceil( ceil(N/NB_A)/NPCOL )*NB_A
79 *
80 * Arguments
81 * =========
82 *
83 * UPLO (global input) CHARACTER
84 * Specifies the part of the distributed matrix sub( A ) to be
85 * set:
86 * = 'U': Upper triangular part is set; the strictly lower
87 * triangular part of sub( A ) is not changed;
88 * = 'L': Lower triangular part is set; the strictly upper
89 * triangular part of sub( A ) is not changed;
90 * Otherwise: All of the matrix sub( A ) is set.
91 *
92 * M (global input) INTEGER
93 * The number of rows to be operated on i.e the number of rows
94 * of the distributed submatrix sub( A ). M >= 0.
95 *
96 * N (global input) INTEGER
97 * The number of columns to be operated on i.e the number of
98 * columns of the distributed submatrix sub( A ). N >= 0.
99 *
100 * ALPHA (global input) DOUBLE PRECISION
101 * The constant to which the offdiagonal elements are to be
102 * set.
103 *
104 * BETA (global input) DOUBLE PRECISION
105 * The constant to which the diagonal elements are to be set.
106 *
107 * A (local output) DOUBLE PRECISION pointer into the local memory
108 * to an array of dimension (LLD_A,LOCc(JA+N-1)). This array
109 * contains the local pieces of the distributed matrix sub( A )
110 * to be set. On exit, the leading M-by-N submatrix sub( A )
111 * is set as follows:
112 *
113 * if UPLO = 'U', A(IA+i-1,JA+j-1) = ALPHA, 1<=i<=j-1, 1<=j<=N,
114 * if UPLO = 'L', A(IA+i-1,JA+j-1) = ALPHA, j+1<=i<=M, 1<=j<=N,
115 * otherwise, A(IA+i-1,JA+j-1) = ALPHA, 1<=i<=M, 1<=j<=N,
116 * IA+i.NE.JA+j,
117 * and, for all UPLO, A(IA+i-1,JA+i-1) = BETA, 1<=i<=min(M,N).
118 *
119 * IA (global input) INTEGER
120 * The row index in the global array A indicating the first
121 * row of sub( A ).
122 *
123 * JA (global input) INTEGER
124 * The column index in the global array A indicating the
125 * first column of sub( A ).
126 *
127 * DESCA (global and local input) INTEGER array of dimension DLEN_.
128 * The array descriptor for the distributed matrix A.
129 *
130 * =====================================================================
131 *
132 * .. Parameters ..
133  INTEGER BLOCK_CYCLIC_2D, CSRC_, CTXT_, DLEN_, DTYPE_,
134  $ LLD_, MB_, M_, NB_, N_, RSRC_
135  parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
136  $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
137  $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
138 * ..
139 * .. Local Scalars ..
140  INTEGER HEIGHT, IACOL, IAROW, IBASE, ICOFFA, II, IIA,
141  $ IIBEG, IIEND, IINXT, ILEFT, IRIGHT, IROFFA,
142  $ ITOP, JJ, JJA, JJBEG, JJEND, JJNXT, LDA, MBA,
143  $ MP, MPA, MYCOL, MYDIST, MYROW, NBA, NPCOL,
144  $ NPROW, NQ, NQA, WIDE
145 * ..
146 * .. External Subroutines ..
147  EXTERNAL blacs_gridinfo, dlaset, infog2l
148 * ..
149 * .. External Functions ..
150  LOGICAL LSAME
151  INTEGER ICEIL, NUMROC
152  EXTERNAL iceil, lsame, numroc
153 * ..
154 * .. Intrinsic Functions ..
155  INTRINSIC max, min, mod
156 * ..
157 * .. Executable Statements ..
158 *
159  IF( m.EQ.0 .OR. n.EQ.0 )
160  $ RETURN
161 *
162 * Get grid parameters
163 *
164  CALL blacs_gridinfo( desca( ctxt_ ), nprow, npcol, myrow, mycol )
165 *
166  CALL infog2l( ia, ja, desca, nprow, npcol, myrow, mycol, iia, jja,
167  $ iarow, iacol )
168  mba = desca( mb_ )
169  nba = desca( nb_ )
170  lda = desca( lld_ )
171  iroffa = mod( ia-1, mba )
172  icoffa = mod( ja-1, nba )
173 *
174  IF( n.LE.( nba-icoffa ) ) THEN
175 *
176 * It is assumed that the local columns JJA:JJA+N-1 of the matrix
177 * A are in the same process column (IACOL).
178 *
179 * N
180 * JJA JJA+N-1
181 * / --------------------- \
182 * IROFFA| | | |
183 * \ |...................| | ( IAROW )
184 * IIA |x | | MB_A
185 * | x | |
186 * |--x----------------| /
187 * | x |
188 * | x | ITOP
189 * | x | |
190 * | x | /-------\
191 * |-------x-----------| |-------x-----------|
192 * | x | | x |
193 * | x | | x |
194 * | x | | x |
195 * | x | | x |
196 * |------------x------| |------------x------|
197 * | x | \____________/
198 * | x | |
199 * | x | IBASE
200 * | x |
201 * |-----------------x-| Local picture
202 * | x|
203 * | |
204 * | |
205 * | |
206 * |-------------------|
207 * | |
208 * . .
209 * . .
210 * . (IACOL) .
211 *
212  IF( mycol.EQ.iacol ) THEN
213 *
214  mpa = numroc( m+iroffa, mba, myrow, iarow, nprow )
215  IF( mpa.LE.0 )
216  $ RETURN
217  IF( myrow.EQ.iarow )
218  $ mpa = mpa - iroffa
219  mydist = mod( myrow-iarow+nprow, nprow )
220  itop = mydist * mba - iroffa
221 *
222  IF( lsame( uplo, 'U' ) ) THEN
223 *
224  itop = max( 0, itop )
225  iibeg = iia
226  iiend = iia + mpa - 1
227  iinxt = min( iceil( iibeg, mba ) * mba, iiend )
228 *
229  10 CONTINUE
230  IF( ( n-itop ).GT.0 ) THEN
231  CALL dlaset( uplo, iinxt-iibeg+1, n-itop, alpha, beta,
232  $ a( iibeg+(jja+itop-1)*lda ), lda )
233  mydist = mydist + nprow
234  itop = mydist * mba - iroffa
235  iibeg = iinxt +1
236  iinxt = min( iinxt+mba, iiend )
237  GO TO 10
238  END IF
239 *
240  ELSE IF( lsame( uplo, 'L' ) ) THEN
241 *
242  ii = iia
243  jj = jja
244  mp = mpa
245  ibase = min( itop+mba, n )
246  itop = min( max( 0, itop ), n )
247 *
248  20 CONTINUE
249  IF( jj.LE.( jja+n-1 ) ) THEN
250  height = ibase - itop
251  CALL dlaset( 'All', mp, itop-jj+jja, alpha, alpha,
252  $ a( ii+(jj-1)*lda ), lda )
253  CALL dlaset( uplo, mp, height, alpha, beta,
254  $ a( ii+(jja+itop-1)*lda ), lda )
255  mp = max( 0, mp - height )
256  ii = ii + height
257  jj = jja + ibase
258  mydist = mydist + nprow
259  itop = mydist * mba - iroffa
260  ibase = min( itop + mba, n )
261  itop = min( itop, n )
262  GO TO 20
263  END IF
264 *
265  ELSE
266 *
267  ii = iia
268  jj = jja
269  mp = mpa
270  ibase = min( itop+mba, n )
271  itop = min( max( 0, itop ), n )
272 *
273  30 CONTINUE
274  IF( jj.LE.( jja+n-1 ) ) THEN
275  height = ibase - itop
276  CALL dlaset( 'All', mpa, itop-jj+jja, alpha, alpha,
277  $ a( iia+(jj-1)*lda ), lda )
278  CALL dlaset( 'All', mpa-mp, height, alpha, alpha,
279  $ a( iia+(jja+itop-1)*lda ), lda )
280  CALL dlaset( 'All', mp, height, alpha, beta,
281  $ a( ii+(jja+itop-1)*lda ), lda )
282  mp = max( 0, mp - height )
283  ii = ii + height
284  jj = jja + ibase
285  mydist = mydist + nprow
286  itop = mydist * mba - iroffa
287  ibase = min( itop + mba, n )
288  itop = min( itop, n )
289  GO TO 30
290  END IF
291 *
292  END IF
293 *
294  END IF
295 *
296  ELSE IF( m.LE.( mba-iroffa ) ) THEN
297 *
298 * It is assumed that the local rows IIA:IIA+M-1 of the matrix A
299 * are in the same process row (IAROW).
300 *
301 * ICOFFA
302 * / \JJA
303 * IIA ------------------ .... --------
304 * | .x | | | / | | \
305 * | . x | | | ILEFT| | | |
306 * | . x | | | | | |
307 * | . x | | \ x | |
308 * | . |x | | |x | | IRIGHT
309 * | . | x | | | x | |
310 * (IAROW) | . | x | | | x | |
311 * | . | x| | | x| |
312 * | . | x | | x /
313 * | . | |x | | |
314 * | . | | x | | |
315 * | . | | x | | |
316 * | . | | x| | |
317 * IIA+M-1 ------------------ .... -------
318 * NB_A
319 * (IACOL) Local picture
320 *
321  IF( myrow.EQ.iarow ) THEN
322 *
323  nqa = numroc( n+icoffa, nba, mycol, iacol, npcol )
324  IF( nqa.LE.0 )
325  $ RETURN
326  IF( mycol.EQ.iacol )
327  $ nqa = nqa - icoffa
328  mydist = mod( mycol-iacol+npcol, npcol )
329  ileft = mydist * nba - icoffa
330 *
331  IF( lsame( uplo, 'L' ) ) THEN
332 *
333  ileft = max( 0, ileft )
334  jjbeg = jja
335  jjend = jja + nqa - 1
336  jjnxt = min( iceil( jjbeg, nba ) * nba, jjend )
337 *
338  40 CONTINUE
339  IF( ( m-ileft ).GT.0 ) THEN
340  CALL dlaset( uplo, m-ileft, jjnxt-jjbeg+1, alpha,
341  $ beta, a( iia+ileft+(jjbeg-1)*lda ), lda )
342  mydist = mydist + npcol
343  ileft = mydist * nba - icoffa
344  jjbeg = jjnxt +1
345  jjnxt = min( jjnxt+nba, jjend )
346  GO TO 40
347  END IF
348 *
349  ELSE IF( lsame( uplo, 'U' ) ) THEN
350 *
351  ii = iia
352  jj = jja
353  nq = nqa
354  iright = min( ileft+nba, m )
355  ileft = min( max( 0, ileft ), m )
356 *
357  50 CONTINUE
358  IF( ii.LE.( iia+m-1 ) ) THEN
359  wide = iright - ileft
360  CALL dlaset( 'All', ileft-ii+iia, nq, alpha, alpha,
361  $ a( ii+(jj-1)*lda ), lda )
362  CALL dlaset( uplo, wide, nq, alpha, beta,
363  $ a( iia+ileft+(jj-1)*lda ), lda )
364  nq = max( 0, nq - wide )
365  ii = iia + iright
366  jj = jj + wide
367  mydist = mydist + npcol
368  ileft = mydist * nba - icoffa
369  iright = min( ileft + nba, m )
370  ileft = min( ileft, m )
371  GO TO 50
372  END IF
373 *
374  ELSE
375 *
376  ii = iia
377  jj = jja
378  nq = nqa
379  iright = min( ileft+nba, m )
380  ileft = min( max( 0, ileft ), m )
381 *
382  60 CONTINUE
383  IF( ii.LE.( iia+m-1 ) ) THEN
384  wide = iright - ileft
385  CALL dlaset( 'All', ileft-ii+iia, nqa, alpha, alpha,
386  $ a( ii+(jja-1)*lda ), lda )
387  CALL dlaset( 'All', wide, nqa-nq, alpha, alpha,
388  $ a( iia+ileft+(jja-1)*lda ), lda )
389  CALL dlaset( 'All', wide, nq, alpha, beta,
390  $ a( iia+ileft+(jj-1)*lda ), lda )
391  nq = max( 0, nq - wide )
392  ii = iia + iright
393  jj = jj + wide
394  mydist = mydist + npcol
395  ileft = mydist * nba - icoffa
396  iright = min( ileft + nba, m )
397  ileft = min( ileft, m )
398  GO TO 60
399  END IF
400 *
401  END IF
402 *
403  END IF
404 *
405  END IF
406 *
407  RETURN
408 *
409 * End of PDLASE2
410 *
411  END
pdlase2
subroutine pdlase2(UPLO, M, N, ALPHA, BETA, A, IA, JA, DESCA)
Definition: pdlase2.f:2
max
#define max(A, B)
Definition: pcgemr.c:180
infog2l
subroutine infog2l(GRINDX, GCINDX, DESC, NPROW, NPCOL, MYROW, MYCOL, LRINDX, LCINDX, RSRC, CSRC)
Definition: infog2l.f:3
min
#define min(A, B)
Definition: pcgemr.c:181