ScaLAPACK 2.1  2.1
ScaLAPACK: Scalable Linear Algebra PACKage
pdgetrs.f
Go to the documentation of this file.
1  SUBROUTINE pdgetrs( TRANS, N, NRHS, A, IA, JA, DESCA, IPIV, B,
2  $ IB, JB, DESCB, INFO )
3 *
4 * -- ScaLAPACK routine (version 1.7) --
5 * University of Tennessee, Knoxville, Oak Ridge National Laboratory,
6 * and University of California, Berkeley.
7 * May 1, 1997
8 *
9 * .. Scalar Arguments ..
10  CHARACTER TRANS
11  INTEGER IA, IB, INFO, JA, JB, N, NRHS
12 * ..
13 * .. Array Arguments ..
14  INTEGER DESCA( * ), DESCB( * ), IPIV( * )
15  DOUBLE PRECISION A( * ), B( * )
16 * ..
17 *
18 * Purpose
19 * =======
20 *
21 * PDGETRS solves a system of distributed linear equations
22 *
23 * op( sub( A ) ) * X = sub( B )
24 *
25 * with a general N-by-N distributed matrix sub( A ) using the LU
26 * factorization computed by PDGETRF.
27 * sub( A ) denotes A(IA:IA+N-1,JA:JA+N-1), op( A ) = A or A**T and
28 * sub( B ) denotes B(IB:IB+N-1,JB:JB+NRHS-1).
29 *
30 * Notes
31 * =====
32 *
33 * Each global data object is described by an associated description
34 * vector. This vector stores the information required to establish
35 * the mapping between an object element and its corresponding process
36 * and memory location.
37 *
38 * Let A be a generic term for any 2D block cyclicly distributed array.
39 * Such a global array has an associated description vector DESCA.
40 * In the following comments, the character _ should be read as
41 * "of the global array".
42 *
43 * NOTATION STORED IN EXPLANATION
44 * --------------- -------------- --------------------------------------
45 * DTYPE_A(global) DESCA( DTYPE_ )The descriptor type. In this case,
46 * DTYPE_A = 1.
47 * CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
48 * the BLACS process grid A is distribu-
49 * ted over. The context itself is glo-
50 * bal, but the handle (the integer
51 * value) may vary.
52 * M_A (global) DESCA( M_ ) The number of rows in the global
53 * array A.
54 * N_A (global) DESCA( N_ ) The number of columns in the global
55 * array A.
56 * MB_A (global) DESCA( MB_ ) The blocking factor used to distribute
57 * the rows of the array.
58 * NB_A (global) DESCA( NB_ ) The blocking factor used to distribute
59 * the columns of the array.
60 * RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
61 * row of the array A is distributed.
62 * CSRC_A (global) DESCA( CSRC_ ) The process column over which the
63 * first column of the array A is
64 * distributed.
65 * LLD_A (local) DESCA( LLD_ ) The leading dimension of the local
66 * array. LLD_A >= MAX(1,LOCr(M_A)).
67 *
68 * Let K be the number of rows or columns of a distributed matrix,
69 * and assume that its process grid has dimension p x q.
70 * LOCr( K ) denotes the number of elements of K that a process
71 * would receive if K were distributed over the p processes of its
72 * process column.
73 * Similarly, LOCc( K ) denotes the number of elements of K that a
74 * process would receive if K were distributed over the q processes of
75 * its process row.
76 * The values of LOCr() and LOCc() may be determined via a call to the
77 * ScaLAPACK tool function, NUMROC:
78 * LOCr( M ) = NUMROC( M, MB_A, MYROW, RSRC_A, NPROW ),
79 * LOCc( N ) = NUMROC( N, NB_A, MYCOL, CSRC_A, NPCOL ).
80 * An upper bound for these quantities may be computed by:
81 * LOCr( M ) <= ceil( ceil(M/MB_A)/NPROW )*MB_A
82 * LOCc( N ) <= ceil( ceil(N/NB_A)/NPCOL )*NB_A
83 *
84 * This routine requires square block data decomposition ( MB_A=NB_A ).
85 *
86 * Arguments
87 * =========
88 *
89 * TRANS (global input) CHARACTER
90 * Specifies the form of the system of equations:
91 * = 'N': sub( A ) * X = sub( B ) (No transpose)
92 * = 'T': sub( A )**T * X = sub( B ) (Transpose)
93 * = 'C': sub( A )**T * X = sub( B ) (Transpose)
94 *
95 * N (global input) INTEGER
96 * The number of rows and columns to be operated on, i.e. the
97 * order of the distributed submatrix sub( A ). N >= 0.
98 *
99 * NRHS (global input) INTEGER
100 * The number of right hand sides, i.e., the number of columns
101 * of the distributed submatrix sub( B ). NRHS >= 0.
102 *
103 * A (local input) DOUBLE PRECISION pointer into the local
104 * memory to an array of dimension (LLD_A, LOCc(JA+N-1)).
105 * On entry, this array contains the local pieces of the factors
106 * L and U from the factorization sub( A ) = P*L*U; the unit
107 * diagonal elements of L are not stored.
108 *
109 * IA (global input) INTEGER
110 * The row index in the global array A indicating the first
111 * row of sub( A ).
112 *
113 * JA (global input) INTEGER
114 * The column index in the global array A indicating the
115 * first column of sub( A ).
116 *
117 * DESCA (global and local input) INTEGER array of dimension DLEN_.
118 * The array descriptor for the distributed matrix A.
119 *
120 * IPIV (local input) INTEGER array, dimension ( LOCr(M_A)+MB_A )
121 * This array contains the pivoting information.
122 * IPIV(i) -> The global row local row i was swapped with.
123 * This array is tied to the distributed matrix A.
124 *
125 * B (local input/local output) DOUBLE PRECISION pointer into the
126 * local memory to an array of dimension
127 * (LLD_B,LOCc(JB+NRHS-1)). On entry, the right hand sides
128 * sub( B ). On exit, sub( B ) is overwritten by the solution
129 * distributed matrix X.
130 *
131 * IB (global input) INTEGER
132 * The row index in the global array B indicating the first
133 * row of sub( B ).
134 *
135 * JB (global input) INTEGER
136 * The column index in the global array B indicating the
137 * first column of sub( B ).
138 *
139 * DESCB (global and local input) INTEGER array of dimension DLEN_.
140 * The array descriptor for the distributed matrix B.
141 *
142 * INFO (global output) INTEGER
143 * = 0: successful exit
144 * < 0: If the i-th argument is an array and the j-entry had
145 * an illegal value, then INFO = -(i*100+j), if the i-th
146 * argument is a scalar and had an illegal value, then
147 * INFO = -i.
148 *
149 * =====================================================================
150 *
151 * .. Parameters ..
152  INTEGER BLOCK_CYCLIC_2D, CSRC_, CTXT_, DLEN_, DTYPE_,
153  $ lld_, mb_, m_, nb_, n_, rsrc_
154  parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
155  $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
156  $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
157  DOUBLE PRECISION ONE
158  parameter( one = 1.0d+0 )
159 * ..
160 * .. Local Scalars ..
161  LOGICAL NOTRAN
162  INTEGER IAROW, IBROW, ICOFFA, ICTXT, IROFFA, IROFFB,
163  $ mycol, myrow, npcol, nprow
164 * ..
165 * .. Local Arrays ..
166  INTEGER DESCIP( DLEN_ ), IDUM1( 1 ), IDUM2( 1 )
167 * ..
168 * .. External Subroutines ..
169  EXTERNAL blacs_gridinfo, chk1mat, descset, pchk2mat,
170  $ pdlapiv, pdtrsm, pxerbla
171 * ..
172 * .. External Functions ..
173  LOGICAL LSAME
174  INTEGER INDXG2P, NUMROC
175  EXTERNAL indxg2p, lsame, numroc
176 * ..
177 * .. Intrinsic Functions ..
178  INTRINSIC ichar, mod
179 * ..
180 * .. Executable Statements ..
181 *
182 * Get grid parameters
183 *
184  ictxt = desca( ctxt_ )
185  CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
186 *
187 * Test the input parameters
188 *
189  info = 0
190  IF( nprow.EQ.-1 ) THEN
191  info = -(700+ctxt_)
192  ELSE
193  notran = lsame( trans, 'N' )
194  CALL chk1mat( n, 2, n, 2, ia, ja, desca, 7, info )
195  CALL chk1mat( n, 2, nrhs, 3, ib, jb, descb, 12, info )
196  IF( info.EQ.0 ) THEN
197  iarow = indxg2p( ia, desca( mb_ ), myrow, desca( rsrc_ ),
198  $ nprow )
199  ibrow = indxg2p( ib, descb( mb_ ), myrow, descb( rsrc_ ),
200  $ nprow )
201  iroffa = mod( ia-1, desca( mb_ ) )
202  icoffa = mod( ja-1, desca( nb_ ) )
203  iroffb = mod( ib-1, descb( mb_ ) )
204  IF( .NOT.notran .AND. .NOT.lsame( trans, 'T' ) .AND. .NOT.
205  $ lsame( trans, 'C' ) ) THEN
206  info = -1
207  ELSE IF( iroffa.NE.0 ) THEN
208  info = -5
209  ELSE IF( icoffa.NE.0 ) THEN
210  info = -6
211  ELSE IF( desca( mb_ ).NE.desca( nb_ ) ) THEN
212  info = -(700+nb_)
213  ELSE IF( iroffb.NE.0 .OR. ibrow.NE.iarow ) THEN
214  info = -10
215  ELSE IF( descb( mb_ ).NE.desca( nb_ ) ) THEN
216  info = -(1200+nb_)
217  ELSE IF( ictxt.NE.descb( ctxt_ ) ) THEN
218  info = -(1200+ctxt_)
219  END IF
220  END IF
221  IF( notran ) THEN
222  idum1( 1 ) = ichar( 'N' )
223  ELSE IF( lsame( trans, 'T' ) ) THEN
224  idum1( 1 ) = ichar( 'T' )
225  ELSE
226  idum1( 1 ) = ichar( 'C' )
227  END IF
228  idum2( 1 ) = 1
229  CALL pchk2mat( n, 2, n, 2, ia, ja, desca, 7, n, 2, nrhs, 3,
230  $ ib, jb, descb, 12, 1, idum1, idum2, info )
231  END IF
232 *
233  IF( info.NE.0 ) THEN
234  CALL pxerbla( ictxt, 'PDGETRS', -info )
235  RETURN
236  END IF
237 *
238 * Quick return if possible
239 *
240  IF( n.EQ.0 .OR. nrhs.EQ.0 )
241  $ RETURN
242 *
243  CALL descset( descip, desca( m_ ) + desca( mb_ )*nprow, 1,
244  $ desca( mb_ ), 1, desca( rsrc_ ), mycol, ictxt,
245  $ desca( mb_ ) + numroc( desca( m_ ), desca( mb_ ),
246  $ myrow, desca( rsrc_ ), nprow ) )
247 *
248  IF( notran ) THEN
249 *
250 * Solve sub( A ) * X = sub( B ).
251 *
252 * Apply row interchanges to the right hand sides.
253 *
254  CALL pdlapiv( 'Forward', 'Row', 'Col', n, nrhs, b, ib, jb,
255  $ descb, ipiv, ia, 1, descip, idum1 )
256 *
257 * Solve L*X = sub( B ), overwriting sub( B ) with X.
258 *
259  CALL pdtrsm( 'Left', 'Lower', 'No transpose', 'Unit', n, nrhs,
260  $ one, a, ia, ja, desca, b, ib, jb, descb )
261 *
262 * Solve U*X = sub( B ), overwriting sub( B ) with X.
263 *
264  CALL pdtrsm( 'Left', 'Upper', 'No transpose', 'Non-unit', n,
265  $ nrhs, one, a, ia, ja, desca, b, ib, jb, descb )
266  ELSE
267 *
268 * Solve sub( A )' * X = sub( B ).
269 *
270 * Solve U'*X = sub( B ), overwriting sub( B ) with X.
271 *
272  CALL pdtrsm( 'Left', 'Upper', 'Transpose', 'Non-unit', n, nrhs,
273  $ one, a, ia, ja, desca, b, ib, jb, descb )
274 *
275 * Solve L'*X = sub( B ), overwriting sub( B ) with X.
276 *
277  CALL pdtrsm( 'Left', 'Lower', 'Transpose', 'Unit', n, nrhs,
278  $ one, a, ia, ja, desca, b, ib, jb, descb )
279 *
280 * Apply row interchanges to the solution vectors.
281 *
282  CALL pdlapiv( 'Backward', 'Row', 'Col', n, nrhs, b, ib, jb,
283  $ descb, ipiv, ia, 1, descip, idum1 )
284 *
285  END IF
286 *
287  RETURN
288 *
289 * End of PDGETRS
290 *
291  END
pdgetrs
subroutine pdgetrs(TRANS, N, NRHS, A, IA, JA, DESCA, IPIV, B, IB, JB, DESCB, INFO)
Definition: pdgetrs.f:3
pdlapiv
subroutine pdlapiv(DIREC, ROWCOL, PIVROC, M, N, A, IA, JA, DESCA, IPIV, IP, JP, DESCIP, IWORK)
Definition: pdlapiv.f:3
pchk2mat
subroutine pchk2mat(MA, MAPOS0, NA, NAPOS0, IA, JA, DESCA, DESCAPOS0, MB, MBPOS0, NB, NBPOS0, IB, JB, DESCB, DESCBPOS0, NEXTRA, EX, EXPOS, INFO)
Definition: pchkxmat.f:175
descset
subroutine descset(DESC, M, N, MB, NB, IRSRC, ICSRC, ICTXT, LLD)
Definition: descset.f:3
chk1mat
subroutine chk1mat(MA, MAPOS0, NA, NAPOS0, IA, JA, DESCA, DESCAPOS0, INFO)
Definition: chk1mat.f:3
pxerbla
subroutine pxerbla(ICTXT, SRNAME, INFO)
Definition: pxerbla.f:2