SpM Handbook 1.2.4
Loading...
Searching...
No Matches
spm_rhs.c
Go to the documentation of this file.
1/**
2 *
3 * @file spm_rhs.c
4 *
5 * SParse Matrix package RHS main routines.
6 *
7 * @copyright 2016-2024 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
8 * Univ. Bordeaux. All rights reserved.
9 *
10 * @version 1.2.4
11 * @author Pierre Ramet
12 * @author Mathieu Faverge
13 * @author Tony Delarue
14 * @date 2024-05-29
15 *
16 * @addtogroup spm
17 * @{
18 **/
19#include "common.h"
20
21/**
22 *******************************************************************************
23 *
24 * @brief Print a set of vector associated to an spm matrix.
25 *
26 *******************************************************************************
27 *
28 * @param[in] spm
29 * The sparse matrix.
30 *
31 * @param[in] nrhs
32 * The number of columns of x.
33 *
34 * @param[in] x
35 * The set of vectors associated to the spm of size ldx-by-nrhs.
36 *
37 * @param[in] ldx
38 * The local leading dimension of the set of vectors (ldx >= spm->n).
39 *
40 * @param[in] stream
41 * File to print the spm matrix. stdout, if stream == NULL.
42 *
43 *******************************************************************************/
44void
46 int nrhs,
47 const void *x,
48 spm_int_t ldx,
49 FILE *stream )
50{
51 if (stream == NULL) {
52 stream = stdout;
53 }
54
55 switch(spm->flttype)
56 {
57 case SpmPattern:
58 /* Not handled for now */
59 break;
60 case SpmFloat:
61 s_spmPrintRHS( stream, spm, nrhs, x, ldx );
62 break;
63 case SpmComplex32:
64 c_spmPrintRHS( stream, spm, nrhs, x, ldx );
65 break;
66 case SpmComplex64:
67 z_spmPrintRHS( stream, spm, nrhs, x, ldx );
68 break;
69 case SpmDouble:
70 default:
71 d_spmPrintRHS( stream, spm, nrhs, x, ldx );
72 }
73
74 return;
75}
76
77/**
78 *******************************************************************************
79 *
80 * @brief Generate right hand side vectors associated to a given matrix.
81 *
82 * The vectors can be initialized randomly or to get a specific solution.
83 *
84 *******************************************************************************
85 *
86 * @param[in] type
87 * Defines how to compute the vector b.
88 * - SpmRhsOne: b is computed such that x = 1 [ + I ]
89 * - SpmRhsI: b is computed such that x = i [ + i * I ]
90 * - SpmRhsRndX: b is computed by matrix-vector product, such that
91 * x is a random vector in the range [-0.5, 0.5]
92 * - SpmRhsRndB: b is computed randomly and x is not computed.
93 *
94 * @param[in] nrhs
95 * Defines the number of right hand side that must be generated.
96 *
97 * @param[in] spm
98 * The sparse matrix used to generate the right hand side, and the
99 * solution of the full problem.
100 *
101 * @param[out] x
102 * On exit, if x != NULL, then the x vector(s) generated to compute b
103 * is returned. Must be of size at least ldx * nrhs.
104 *
105 * @param[in] ldx
106 * Defines the leading dimension of x when multiple right hand sides
107 * are available. ldx >= max( 1, spm->nexp ).
108 *
109 * @param[inout] b
110 * b must be an allocated matrix of size at least ldb * nrhs.
111 * On exit, b is initialized as defined by the type parameter.
112 *
113 * @param[in] ldb
114 * Defines the leading dimension of b when multiple right hand sides
115 * are available. ldb >= max( 1, spm->nexp ).
116 *
117 *******************************************************************************
118 *
119 * @retval SPM_SUCCESS if the b vector has been computed successfully,
120 * @retval SPM_ERR_BADPARAMETER otherwise.
121 *
122 *******************************************************************************/
123int
125 spm_int_t nrhs,
126 const spmatrix_t *spm,
127 void *x,
128 spm_int_t ldx,
129 void *b,
130 spm_int_t ldb )
131{
132 static int (*ptrfunc[4])( spm_rhstype_t, int,
133 const spmatrix_t *,
134 void *, int, void *, int ) =
135 {
137 };
138
139 int id = spm->flttype - SpmFloat;
140
141 if ( (x != NULL) && (ldx < spm_imax( 1, spm->nexp )) ) {
142 fprintf( stderr, "spmGenRHS: ldx must be >= max( 1, spm->nexp )\n" );
144 }
145 if ( ldb < spm_imax( 1, spm->nexp ) ) {
146 fprintf( stderr, "spmGenRHS: ldb must be >= max( 1, spm->nexp )\n" );
148 }
149
150 if ( (id < 0) || (id > 3) ) {
152 }
153 else {
154 return ptrfunc[id]( type, nrhs, spm, x, ldx, b, ldb );
155 }
156}
157
158/**
159 *******************************************************************************
160 *
161 * @brief Stores the local values of a global RHS in the local one thanks to spm
162 * distribution.
163 *
164 *******************************************************************************
165 *
166 * @param[in] nrhs
167 * Defines the number of right hand side.
168 *
169 * @param[in] spm
170 * The sparse matrix used to generate the right hand side.
171 *
172 * @param[in] bglob
173 * The ldbg -by- nrhs global dense matrix
174 *
175 * @param[in] ldbg
176 * The leading dimension of bglob. ldbg >= max( 1, spm->gNexp ).
177 *
178 * @param[inout] bloc
179 * The ldbl -by- nrhs local dense matrix
180 * Must be allocated on entry, on exit contains the local values of
181 * bglob.
182 *
183 * @param[in] ldbl
184 * The leading dimension of bloc. ldbl >= max( 1, spm->nexp ).
185 *
186 *******************************************************************************
187 *
188 * @retval SPM_SUCCESS if the b vector has been Reduceed successfully,
189 * @retval SPM_ERR_BADPARAMETER otherwise.
190 *
191 *******************************************************************************/
192int
194 const spmatrix_t *spm,
195 const void *bglob,
196 spm_int_t ldbg,
197 void *bloc,
198 spm_int_t ldbl )
199{
200 if ( (spm == NULL) || (spm->values == NULL) ) {
202 }
203
204 if ( bglob == NULL ) {
206 }
207
208 if ( bloc == NULL ) {
210 }
211
212 if ( ldbg < spm_imax( 1, spm->gNexp ) ) {
213 fprintf( stderr, "spmExtractLocalRHS: ldbg must be >= max( 1, spm->gNexp )\n" );
215 }
216
217 switch (spm->flttype)
218 {
219 case SpmFloat:
220 s_spmExtractLocalRHS( nrhs, spm, bglob, ldbg, bloc, ldbl );
221 break;
222 case SpmComplex32:
223 c_spmExtractLocalRHS( nrhs, spm, bglob, ldbg, bloc, ldbl );
224 break;
225 case SpmComplex64:
226 z_spmExtractLocalRHS( nrhs, spm, bglob, ldbg, bloc, ldbl );
227 break;
228 case SpmDouble:
229 default:
230 d_spmExtractLocalRHS( nrhs, spm, bglob, ldbg, bloc, ldbl );
231 break;
232 }
233
234 return SPM_SUCCESS;
235}
236
237/**
238 *******************************************************************************
239 *
240 * @brief Reduce an RHS thanks to spm distribution.
241 *
242 *******************************************************************************
243 *
244 * @param[in] nrhs
245 * Defines the number of right hand side.
246 *
247 * @param[in] spm
248 * The sparse matrix used to generate the right hand side.
249 *
250 * @param[inout] bglob
251 * The ldbg -by- nrhs global dense matrix
252 * On entry, the partial accumulation of a matrix vector product.
253 * On exit, each node conatins the reduction sum of the matrices.
254 *
255 * @param[in] ldbg
256 * The leading dimension of bglob. ldbg >= max( 1, spm->gNexp ).
257 *
258 * @param[inout] bloc
259 * The ldbl -by- nrhs local dense matrix
260 * On entry, must be allocated.
261 * On exit, contains the local extraction of the final bglob.
262 *
263 * @param[in] ldbl
264 * The leading dimension of bloc. ldbl >= max( 1, spm->nexp ).
265 *
266 *******************************************************************************
267 *
268 * @retval SPM_SUCCESS if the bglob vector has been Reduceed successfully,
269 * @retval SPM_ERR_BADPARAMETER otherwise.
270 *
271 *******************************************************************************/
272int
274 const spmatrix_t *spm,
275 void *bglob,
276 spm_int_t ldbg,
277 void *bloc,
278 spm_int_t ldbl )
279{
280 if ( (spm == NULL) || (spm->values == NULL) ) {
282 }
283
284 if ( bglob == NULL ) {
286 }
287
288 if ( bloc == NULL ) {
290 }
291
292 if ( ldbg < spm_imax( 1, spm->gNexp ) ) {
293 fprintf( stderr, "spmReduceRHS: ldbg must be >= max( 1, spm->gNexp )\n" );
295 }
296
297 switch (spm->flttype)
298 {
299 case SpmFloat:
300 s_spmReduceRHS( nrhs, spm, bglob, ldbg, bloc, ldbl );
301 break;
302 case SpmComplex32:
303 c_spmReduceRHS( nrhs, spm, bglob, ldbg, bloc, ldbl );
304 break;
305 case SpmComplex64:
306 z_spmReduceRHS( nrhs, spm, bglob, ldbg, bloc, ldbl );
307 break;
308 case SpmDouble:
309 default:
310 d_spmReduceRHS( nrhs, spm, bglob, ldbg, bloc, ldbl );
311 break;
312 }
313
314 return SPM_SUCCESS;
315}
316
317/**
318 *******************************************************************************
319 *
320 * @brief Gather an RHS thanks to spm distribution.
321 *
322 *******************************************************************************
323 *
324 * @param[in] nrhs
325 * Defines the number of right hand side.
326 *
327 * @param[in] spm
328 * The sparse matrix used to generate the right hand side.
329 *
330 * @param[in] bloc
331 * The ldbl -by- nrhs local dense matrix
332 * On entry, the local dense matrix to gather
333 *
334 * @param[in] ldbl
335 * The leading dimension of bloc. ldbl >= max( 1, spm->nexp ).
336 *
337 * @param[in] root
338 * Clustnum where the complete vector will be gathered.
339 * Set it to -1 if you want to gather the global RHS on all nodes.
340 *
341 * @param[inout] bglob
342 * The ldbg -by- nrhs global dense matrix
343 * On entry, the allocated matrix if on root. Non referenced otherwise.
344 * On exit, on the root node, contains the global dense matrix.
345 *
346 * @param[in] ldbg
347 * The leading dimension of bglob. ldbg >= max( 1, spm->gNexp ).
348 *
349 *******************************************************************************
350 *
351 * @retval SPM_SUCCESS if the b vector has been gathered successfully,
352 * @retval SPM_ERR_BADPARAMETER otherwise.
353 *
354 *******************************************************************************/
355int
357 const spmatrix_t *spm,
358 const void *bloc,
359 spm_int_t ldbl,
360 int root,
361 void *bglob,
362 spm_int_t ldbg )
363{
364 if ( (spm == NULL) || (spm->values == NULL) ) {
366 }
367
368 if ( bloc == NULL ) {
370 }
371
372 if ( ldbl < spm_imax( 1, spm->nexp ) ) {
373 fprintf( stderr, "spmGatherRHS: ldbl must be >= max( 1, spm->nexp )\n" );
375 }
376
377 if ( ((root == -1) || (root == spm->clustnum)) && (bglob == NULL) ) {
379 }
380
381 if ( ldbg < spm_imax( 1, spm->gNexp ) ) {
382 fprintf( stderr, "spmGatherRHS: ldbg must be >= max( 1, spm->gNexp )\n" );
384 }
385
386 switch (spm->flttype)
387 {
388 case SpmFloat:
389 s_spmGatherRHS( nrhs, spm, bloc, ldbl, root, bglob, ldbg );
390 break;
391 case SpmComplex32:
392 c_spmGatherRHS( nrhs, spm, bloc, ldbl, root, bglob, ldbg );
393 break;
394 case SpmComplex64:
395 z_spmGatherRHS( nrhs, spm, bloc, ldbl, root, bglob, ldbg );
396 break;
397 case SpmDouble:
398 default:
399 d_spmGatherRHS( nrhs, spm, bloc, ldbl, root, bglob, ldbg );
400 break;
401 }
402
403 return SPM_SUCCESS;
404}
405
406/**
407 * @}
408 */
static spm_int_t spm_imax(spm_int_t a, spm_int_t b)
Internal function to compute max(a,b)
Definition datatypes.h:119
enum spm_rhstype_e spm_rhstype_t
How to generate RHS.
@ SPM_ERR_BADPARAMETER
Definition const.h:89
@ SPM_SUCCESS
Definition const.h:82
@ SpmComplex64
Definition const.h:66
@ SpmDouble
Definition const.h:64
@ SpmFloat
Definition const.h:63
@ SpmComplex32
Definition const.h:65
@ SpmPattern
Definition const.h:62
void s_spmPrintRHS(FILE *f, const spmatrix_t *spm, int nrhs, const void *x, spm_int_t ldx)
Write into a file the vectors associated to a spm.
void z_spmPrintRHS(FILE *f, const spmatrix_t *spm, int nrhs, const void *x, spm_int_t ldx)
Write into a file the vectors associated to a spm.
void d_spmPrintRHS(FILE *f, const spmatrix_t *spm, int nrhs, const void *x, spm_int_t ldx)
Write into a file the vectors associated to a spm.
void c_spmPrintRHS(FILE *f, const spmatrix_t *spm, int nrhs, const void *x, spm_int_t ldx)
Write into a file the vectors associated to a spm.
int c_spmGenRHS(spm_rhstype_t type, int nrhs, const spmatrix_t *spm, void *x, int ldx, void *b, int ldb)
Generate nrhs right hand side vectors associated to a given matrix to test a problem with a solver.
int s_spmGenRHS(spm_rhstype_t type, int nrhs, const spmatrix_t *spm, void *x, int ldx, void *b, int ldb)
Generate nrhs right hand side vectors associated to a given matrix to test a problem with a solver.
void c_spmReduceRHS(int nrhs, const spmatrix_t *spm, spm_complex32_t *bglob, spm_int_t ldbg, spm_complex32_t *bloc, spm_int_t ldbl)
Reduce all the global coefficients of a rhs and store the local ones.
Definition c_spm_rhs.c:109
int d_spmGenRHS(spm_rhstype_t type, int nrhs, const spmatrix_t *spm, void *x, int ldx, void *b, int ldb)
Generate nrhs right hand side vectors associated to a given matrix to test a problem with a solver.
int z_spmGenRHS(spm_rhstype_t type, int nrhs, const spmatrix_t *spm, void *x, int ldx, void *b, int ldb)
Generate nrhs right hand side vectors associated to a given matrix to test a problem with a solver.
void z_spmReduceRHS(int nrhs, const spmatrix_t *spm, spm_complex64_t *bglob, spm_int_t ldbg, spm_complex64_t *bloc, spm_int_t ldbl)
Reduce all the global coefficients of a rhs and store the local ones.
Definition z_spm_rhs.c:109
void s_spmExtractLocalRHS(int nrhs, const spmatrix_t *spm, const float *bglob, spm_int_t ldbg, float *bloc, spm_int_t ldbl)
Stores the local values of a global RHS in the local one.
Definition s_spm_rhs.c:51
void z_spmGatherRHS(int nrhs, const spmatrix_t *spm, const spm_complex64_t *x, spm_int_t ldx, int root, spm_complex64_t *gx, spm_int_t ldgx)
Gather all the global C coefficients and store the good ones in local.
Definition z_spm_rhs.c:165
void d_spmReduceRHS(int nrhs, const spmatrix_t *spm, double *bglob, spm_int_t ldbg, double *bloc, spm_int_t ldbl)
Reduce all the global coefficients of a rhs and store the local ones.
Definition d_spm_rhs.c:109
void c_spmExtractLocalRHS(int nrhs, const spmatrix_t *spm, const spm_complex32_t *bglob, spm_int_t ldbg, spm_complex32_t *bloc, spm_int_t ldbl)
Stores the local values of a global RHS in the local one.
Definition c_spm_rhs.c:51
void c_spmGatherRHS(int nrhs, const spmatrix_t *spm, const spm_complex32_t *x, spm_int_t ldx, int root, spm_complex32_t *gx, spm_int_t ldgx)
Gather all the global C coefficients and store the good ones in local.
Definition c_spm_rhs.c:165
void d_spmGatherRHS(int nrhs, const spmatrix_t *spm, const double *x, spm_int_t ldx, int root, double *gx, spm_int_t ldgx)
Gather all the global C coefficients and store the good ones in local.
Definition d_spm_rhs.c:165
void d_spmExtractLocalRHS(int nrhs, const spmatrix_t *spm, const double *bglob, spm_int_t ldbg, double *bloc, spm_int_t ldbl)
Stores the local values of a global RHS in the local one.
Definition d_spm_rhs.c:51
void s_spmGatherRHS(int nrhs, const spmatrix_t *spm, const float *x, spm_int_t ldx, int root, float *gx, spm_int_t ldgx)
Gather all the global C coefficients and store the good ones in local.
Definition s_spm_rhs.c:165
void s_spmReduceRHS(int nrhs, const spmatrix_t *spm, float *bglob, spm_int_t ldbg, float *bloc, spm_int_t ldbl)
Reduce all the global coefficients of a rhs and store the local ones.
Definition s_spm_rhs.c:109
void z_spmExtractLocalRHS(int nrhs, const spmatrix_t *spm, const spm_complex64_t *bglob, spm_int_t ldbg, spm_complex64_t *bloc, spm_int_t ldbl)
Stores the local values of a global RHS in the local one.
Definition z_spm_rhs.c:51
spm_coeftype_t flttype
Definition spm.h:67
spm_int_t nexp
Definition spm.h:78
void * values
Definition spm.h:92
spm_int_t gNexp
Definition spm.h:77
int clustnum
Definition spm.h:95
int spmReduceRHS(spm_int_t nrhs, const spmatrix_t *spm, void *Bg, spm_int_t ldbg, void *Bl, spm_int_t ldbl)
Reduce an RHS thanks to spm distribution.
Definition spm_rhs.c:273
int spmGenRHS(spm_rhstype_t type, spm_int_t nrhs, const spmatrix_t *spm, void *opt_X, spm_int_t opt_ldx, void *B, spm_int_t ldb)
Generate right hand side vectors associated to a given matrix.
Definition spm_rhs.c:124
int spm_int_t
The main integer datatype used in spm arrays.
Definition datatypes.h:70
void spmPrintRHS(const spmatrix_t *spm, int nrhs, const void *x, spm_int_t ldx, FILE *stream)
Print a set of vector associated to an spm matrix.
Definition spm_rhs.c:45
int spmExtractLocalRHS(spm_int_t nrhs, const spmatrix_t *spm, const void *Bg, spm_int_t ldbg, void *Bl, spm_int_t ldbl)
Stores the local values of a global RHS in the local one thanks to spm distribution.
Definition spm_rhs.c:193
int spmGatherRHS(spm_int_t nrhs, const spmatrix_t *spm, const void *Bl, spm_int_t ldbl, int root, void *Bg, spm_int_t ldbg)
Gather an RHS thanks to spm distribution.
Definition spm_rhs.c:356
The sparse matrix data structure.
Definition spm.h:64