ScaLAPACK 2.1  2.1
ScaLAPACK: Scalable Linear Algebra PACKage
pzlange.f
Go to the documentation of this file.
1  DOUBLE PRECISION FUNCTION pzlange( NORM, M, N, A, IA, JA, DESCA,
2  $ WORK )
3  IMPLICIT NONE
4 *
5 * -- ScaLAPACK auxiliary routine (version 1.7) --
6 * University of Tennessee, Knoxville, Oak Ridge National Laboratory,
7 * and University of California, Berkeley.
8 * May 1, 1997
9 *
10 * .. Scalar Arguments ..
11  CHARACTER norm
12  INTEGER ia, ja, m, n
13 * ..
14 * .. Array Arguments ..
15  INTEGER desca( * )
16  DOUBLE PRECISION work( * )
17  COMPLEX*16 a( * )
18 * ..
19 *
20 * Purpose
21 * =======
22 *
23 * PZLANGE returns the value of the one norm, or the Frobenius norm,
24 * or the infinity norm, or the element of largest absolute value of a
25 * distributed matrix sub( A ) = A(IA:IA+M-1, JA:JA+N-1).
26 *
27 * PZLANGE returns the value
28 *
29 * ( max(abs(A(i,j))), NORM = 'M' or 'm' with IA <= i <= IA+M-1,
30 * ( and JA <= j <= JA+N-1,
31 * (
32 * ( norm1( sub( A ) ), NORM = '1', 'O' or 'o'
33 * (
34 * ( normI( sub( A ) ), NORM = 'I' or 'i'
35 * (
36 * ( normF( sub( A ) ), NORM = 'F', 'f', 'E' or 'e'
37 *
38 * where norm1 denotes the one norm of a matrix (maximum column sum),
39 * normI denotes the infinity norm of a matrix (maximum row sum) and
40 * normF denotes the Frobenius norm of a matrix (square root of sum of
41 * squares). Note that max(abs(A(i,j))) is not a matrix norm.
42 *
43 * Notes
44 * =====
45 *
46 * Each global data object is described by an associated description
47 * vector. This vector stores the information required to establish
48 * the mapping between an object element and its corresponding process
49 * and memory location.
50 *
51 * Let A be a generic term for any 2D block cyclicly distributed array.
52 * Such a global array has an associated description vector DESCA.
53 * In the following comments, the character _ should be read as
54 * "of the global array".
55 *
56 * NOTATION STORED IN EXPLANATION
57 * --------------- -------------- --------------------------------------
58 * DTYPE_A(global) DESCA( DTYPE_ )The descriptor type. In this case,
59 * DTYPE_A = 1.
60 * CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
61 * the BLACS process grid A is distribu-
62 * ted over. The context itself is glo-
63 * bal, but the handle (the integer
64 * value) may vary.
65 * M_A (global) DESCA( M_ ) The number of rows in the global
66 * array A.
67 * N_A (global) DESCA( N_ ) The number of columns in the global
68 * array A.
69 * MB_A (global) DESCA( MB_ ) The blocking factor used to distribute
70 * the rows of the array.
71 * NB_A (global) DESCA( NB_ ) The blocking factor used to distribute
72 * the columns of the array.
73 * RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
74 * row of the array A is distributed.
75 * CSRC_A (global) DESCA( CSRC_ ) The process column over which the
76 * first column of the array A is
77 * distributed.
78 * LLD_A (local) DESCA( LLD_ ) The leading dimension of the local
79 * array. LLD_A >= MAX(1,LOCr(M_A)).
80 *
81 * Let K be the number of rows or columns of a distributed matrix,
82 * and assume that its process grid has dimension p x q.
83 * LOCr( K ) denotes the number of elements of K that a process
84 * would receive if K were distributed over the p processes of its
85 * process column.
86 * Similarly, LOCc( K ) denotes the number of elements of K that a
87 * process would receive if K were distributed over the q processes of
88 * its process row.
89 * The values of LOCr() and LOCc() may be determined via a call to the
90 * ScaLAPACK tool function, NUMROC:
91 * LOCr( M ) = NUMROC( M, MB_A, MYROW, RSRC_A, NPROW ),
92 * LOCc( N ) = NUMROC( N, NB_A, MYCOL, CSRC_A, NPCOL ).
93 * An upper bound for these quantities may be computed by:
94 * LOCr( M ) <= ceil( ceil(M/MB_A)/NPROW )*MB_A
95 * LOCc( N ) <= ceil( ceil(N/NB_A)/NPCOL )*NB_A
96 *
97 * Arguments
98 * =========
99 *
100 * NORM (global input) CHARACTER
101 * Specifies the value to be returned in PZLANGE as described
102 * above.
103 *
104 * M (global input) INTEGER
105 * The number of rows to be operated on i.e the number of rows
106 * of the distributed submatrix sub( A ). When M = 0, PZLANGE
107 * is set to zero. M >= 0.
108 *
109 * N (global input) INTEGER
110 * The number of columns to be operated on i.e the number of
111 * columns of the distributed submatrix sub( A ). When N = 0,
112 * PZLANGE is set to zero. N >= 0.
113 *
114 * A (local input) COMPLEX*16 pointer into the local memory
115 * to an array of dimension (LLD_A, LOCc(JA+N-1)) containing the
116 * local pieces of the distributed matrix sub( A ).
117 *
118 * IA (global input) INTEGER
119 * The row index in the global array A indicating the first
120 * row of sub( A ).
121 *
122 * JA (global input) INTEGER
123 * The column index in the global array A indicating the
124 * first column of sub( A ).
125 *
126 * DESCA (global and local input) INTEGER array of dimension DLEN_.
127 * The array descriptor for the distributed matrix A.
128 *
129 * WORK (local workspace) DOUBLE PRECISION array dimension (LWORK)
130 * LWORK >= 0 if NORM = 'M' or 'm' (not referenced),
131 * Nq0 if NORM = '1', 'O' or 'o',
132 * Mp0 if NORM = 'I' or 'i',
133 * 0 if NORM = 'F', 'f', 'E' or 'e' (not referenced),
134 * where
135 *
136 * IROFFA = MOD( IA-1, MB_A ), ICOFFA = MOD( JA-1, NB_A ),
137 * IAROW = INDXG2P( IA, MB_A, MYROW, RSRC_A, NPROW ),
138 * IACOL = INDXG2P( JA, NB_A, MYCOL, CSRC_A, NPCOL ),
139 * Mp0 = NUMROC( M+IROFFA, MB_A, MYROW, IAROW, NPROW ),
140 * Nq0 = NUMROC( N+ICOFFA, NB_A, MYCOL, IACOL, NPCOL ),
141 *
142 * INDXG2P and NUMROC are ScaLAPACK tool functions; MYROW,
143 * MYCOL, NPROW and NPCOL can be determined by calling the
144 * subroutine BLACS_GRIDINFO.
145 *
146 * =====================================================================
147 *
148 * .. Parameters ..
149  INTEGER block_cyclic_2d, csrc_, ctxt_, dlen_, dtype_,
150  $ lld_, mb_, m_, nb_, n_, rsrc_
151  parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
152  $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
153  $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
154  DOUBLE PRECISION one, zero
155  parameter( one = 1.0d+0, zero = 0.0d+0 )
156 * ..
157 * .. Local Scalars ..
158  INTEGER i, iacol, iarow, ictxt, ii, icoff, ioffa,
159  $ iroff, j, jj, lda, mp, mycol, myrow, npcol,
160  $ nprow, nq
161  DOUBLE PRECISION sum, value
162 * ..
163 * .. Local Arrays ..
164  DOUBLE PRECISION ssq( 2 ), colssq( 2 )
165 * ..
166 * .. External Subroutines ..
167  EXTERNAL blacs_gridinfo, dcombssq, dgebr2d,
168  $ dgebs2d, dgamx2d, dgsum2d, infog2l,
169  $ pdtreecomb, zlassq
170 * ..
171 * .. External Functions ..
172  LOGICAL lsame
173  INTEGER idamax, numroc
174  EXTERNAL lsame, idamax, numroc
175 * ..
176 * .. Intrinsic Functions ..
177  INTRINSIC abs, max, min, mod, sqrt
178 * ..
179 * .. Executable Statements ..
180 *
181 * Get grid parameters.
182 *
183  ictxt = desca( ctxt_ )
184  CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
185 *
186  CALL infog2l( ia, ja, desca, nprow, npcol, myrow, mycol, ii, jj,
187  $ iarow, iacol )
188  iroff = mod( ia-1, desca( mb_ ) )
189  icoff = mod( ja-1, desca( nb_ ) )
190  mp = numroc( m+iroff, desca( mb_ ), myrow, iarow, nprow )
191  nq = numroc( n+icoff, desca( nb_ ), mycol, iacol, npcol )
192  IF( myrow.EQ.iarow )
193  $ mp = mp - iroff
194  IF( mycol.EQ.iacol )
195  $ nq = nq - icoff
196  lda = desca( lld_ )
197 *
198  IF( min( m, n ).EQ.0 ) THEN
199 *
200  VALUE = zero
201 *
202 ************************************************************************
203 * max norm
204 *
205  ELSE IF( lsame( norm, 'M' ) ) THEN
206 *
207 * Find max(abs(A(i,j))).
208 *
209  VALUE = zero
210  IF( nq.GT.0 .AND. mp.GT.0 ) THEN
211  ioffa = (jj-1)*lda
212  DO 20 j = jj, jj+nq-1
213  DO 10 i = ii, mp+ii-1
214  VALUE = max( VALUE, abs( a( ioffa+i ) ) )
215  10 CONTINUE
216  ioffa = ioffa + lda
217  20 CONTINUE
218  END IF
219  CALL dgamx2d( ictxt, 'All', ' ', 1, 1, VALUE, 1, i, j, -1,
220  $ 0, 0 )
221 *
222 ************************************************************************
223 * one norm
224 *
225  ELSE IF( lsame( norm, 'O' ) .OR. norm.EQ.'1' ) THEN
226 *
227 * Find norm1( sub( A ) ).
228 *
229  IF( nq.GT.0 ) THEN
230  ioffa = ( jj - 1 ) * lda
231  DO 40 j = jj, jj+nq-1
232  sum = zero
233  IF( mp.GT.0 ) THEN
234  DO 30 i = ii, mp+ii-1
235  sum = sum + abs( a( ioffa+i ) )
236  30 CONTINUE
237  END IF
238  ioffa = ioffa + lda
239  work( j-jj+1 ) = sum
240  40 CONTINUE
241  END IF
242 *
243 * Find sum of global matrix columns and store on row 0 of
244 * process grid
245 *
246  CALL dgsum2d( ictxt, 'Columnwise', ' ', 1, nq, work, 1,
247  $ 0, mycol )
248 *
249 * Find maximum sum of columns for 1-norm
250 *
251  IF( myrow.EQ.0 ) THEN
252  IF( nq.GT.0 ) THEN
253  VALUE = work( idamax( nq, work, 1 ) )
254  ELSE
255  VALUE = zero
256  END IF
257  CALL dgamx2d( ictxt, 'Rowwise', ' ', 1, 1, VALUE, 1, i, j,
258  $ -1, 0, 0 )
259  END IF
260 *
261 ************************************************************************
262 * inf norm
263 *
264  ELSE IF( lsame( norm, 'I' ) ) THEN
265 *
266 * Find normI( sub( A ) ).
267 *
268  IF( mp.GT.0 ) THEN
269  ioffa = ii + ( jj - 1 ) * lda
270  DO 60 i = ii, ii+mp-1
271  sum = zero
272  IF( nq.GT.0 ) THEN
273  DO 50 j = ioffa, ioffa + nq*lda - 1, lda
274  sum = sum + abs( a( j ) )
275  50 CONTINUE
276  END IF
277  work( i-ii+1 ) = sum
278  ioffa = ioffa + 1
279  60 CONTINUE
280  END IF
281 *
282 * Find sum of global matrix rows and store on column 0 of
283 * process grid
284 *
285  CALL dgsum2d( ictxt, 'Rowwise', ' ', mp, 1, work, max( 1, mp ),
286  $ myrow, 0 )
287 *
288 * Find maximum sum of rows for supnorm
289 *
290  IF( mycol.EQ.0 ) THEN
291  IF( mp.GT.0 ) THEN
292  VALUE = work( idamax( mp, work, 1 ) )
293  ELSE
294  VALUE = zero
295  END IF
296  CALL dgamx2d( ictxt, 'Columnwise', ' ', 1, 1, VALUE, 1, i,
297  $ j, -1, 0, 0 )
298  END IF
299 *
300 ************************************************************************
301 * Frobenius norm
302 * SSQ(1) is scale
303 * SSQ(2) is sum-of-squares
304 *
305  ELSE IF( ( lsame( norm, 'F' ) ) .OR. ( lsame( norm, 'E' ) ) ) THEN
306 *
307 * Find normF( sub( A ) ).
308 *
309  ssq(1) = zero
310  ssq(2) = one
311  ioffa = ii + ( jj - 1 ) * lda
312  IF( nq.GT.0 ) THEN
313  DO 70 j = ioffa, ioffa + nq*lda - 1, lda
314  colssq(1) = zero
315  colssq(2) = one
316  CALL zlassq( mp, a( j ), 1, colssq(1), colssq(2) )
317  CALL dcombssq( ssq, colssq )
318  70 CONTINUE
319  END IF
320 *
321 * Perform the global scaled sum
322 *
323  CALL pdtreecomb( ictxt, 'All', 2, ssq, 0, 0, dcombssq )
324  VALUE = ssq( 1 ) * sqrt( ssq( 2 ) )
325 *
326  END IF
327 *
328  IF( myrow.EQ.0 .AND. mycol.EQ.0 ) THEN
329  CALL dgebs2d( ictxt, 'All', ' ', 1, 1, VALUE, 1 )
330  ELSE
331  CALL dgebr2d( ictxt, 'All', ' ', 1, 1, VALUE, 1, 0, 0 )
332  END IF
333 *
334  pzlange = VALUE
335 *
336  RETURN
337 *
338 * End of PZLANGE
339 *
340  END
max
#define max(A, B)
Definition: pcgemr.c:180
pdtreecomb
subroutine pdtreecomb(ICTXT, SCOPE, N, MINE, RDEST0, CDEST0, SUBPTR)
Definition: pdtreecomb.f:3
pzlange
double precision function pzlange(NORM, M, N, A, IA, JA, DESCA, WORK)
Definition: pzlange.f:3
infog2l
subroutine infog2l(GRINDX, GCINDX, DESC, NPROW, NPCOL, MYROW, MYCOL, LRINDX, LCINDX, RSRC, CSRC)
Definition: infog2l.f:3
lsame
logical function lsame(CA, CB)
Definition: tools.f:1724
numroc
integer function numroc(N, NB, IPROC, ISRCPROC, NPROCS)
Definition: numroc.f:2
dcombssq
subroutine dcombssq(V1, V2)
Definition: pdtreecomb.f:259
min
#define min(A, B)
Definition: pcgemr.c:181