programmer's documentation
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
cs_sles_it.h
Go to the documentation of this file.
1 #ifndef __CS_SLES_IT_H__
2 #define __CS_SLES_IT_H__
3 
4 /*============================================================================
5  * Sparse Linear Equation Solvers
6  *============================================================================*/
7 
8 /*
9  This file is part of Code_Saturne, a general-purpose CFD tool.
10 
11  Copyright (C) 1998-2015 EDF S.A.
12 
13  This program is free software; you can redistribute it and/or modify it under
14  the terms of the GNU General Public License as published by the Free Software
15  Foundation; either version 2 of the License, or (at your option) any later
16  version.
17 
18  This program is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20  FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
21  details.
22 
23  You should have received a copy of the GNU General Public License along with
24  this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
25  Street, Fifth Floor, Boston, MA 02110-1301, USA.
26 */
27 
28 /*----------------------------------------------------------------------------*/
29 
30 /*----------------------------------------------------------------------------
31  * Local headers
32  *----------------------------------------------------------------------------*/
33 
34 #include "cs_base.h"
35 #include "cs_halo_perio.h"
36 #include "cs_matrix.h"
37 
38 /*----------------------------------------------------------------------------*/
39 
41 
42 /*============================================================================
43  * Macro definitions
44  *============================================================================*/
45 
46 /*============================================================================
47  * Type definitions
48  *============================================================================*/
49 
50 /*----------------------------------------------------------------------------
51  * Solver types
52  *----------------------------------------------------------------------------*/
53 
54 typedef enum {
55 
56  CS_SLES_PCG, /* Preconditionned conjugate gradient */
57  CS_SLES_JACOBI, /* Jacobi */
58  CS_SLES_BICGSTAB, /* Bi-conjugate gradient stabilized */
59  CS_SLES_BICGSTAB2, /* Bi-conjugate gradient stabilized - 2*/
60  CS_SLES_GMRES, /* Generalized minimal residual */
61  CS_SLES_B_GAUSS_SEIDEL, /* Block Gauss-Seidel */
62  CS_SLES_N_IT_TYPES /* Number of resolution algorithms */
63 
65 
66 /* Iterative linear solver context (opaque) */
67 
68 typedef struct _cs_sles_it_t cs_sles_it_t;
69 
70 /*============================================================================
71  * Global variables
72  *============================================================================*/
73 
74 /* Short names for matrix types */
75 
76 extern const char *cs_sles_it_type_name[];
77 
78 /*=============================================================================
79  * Public function prototypes
80  *============================================================================*/
81 
82 /*----------------------------------------------------------------------------
83  * Define and associate an iterative sparse linear system solver
84  * for a given field or equation name.
85  *
86  * If this system did not previously exist, it is added to the list of
87  * "known" systems. Otherwise, its definition is replaced by the one
88  * defined here.
89  *
90  * This is a utility function: if finer control is needed, see
91  * cs_sles_define() and cs_sles_it_create().
92  *
93  * Note that this function returns a pointer directly to the iterative solver
94  * management structure. This may be used to set further options,
95  * for example using cs_sles_it_set_verbosity(). If needed, cs_sles_find()
96  * may be used to obtain a pointer to the matching cs_sles_t container.
97  *
98  * parameters:
99  * f_id <-- associated field id, or < 0
100  * name <-- associated name if f_id < 0, or NULL
101  * solver_type <-- type of solver (PCG, Jacobi, ...)
102  * poly_degree <-- preconditioning polynomial degree
103  * (0: diagonal; -1: non-preconditioned)
104  * n_max_iter <-- maximum number of iterations
105  *
106  * returns:
107  * pointer to newly created iterative solver info object.
108  *----------------------------------------------------------------------------*/
109 
110 cs_sles_it_t *
111 cs_sles_it_define(int f_id,
112  const char *name,
113  cs_sles_it_type_t solver_type,
114  int poly_degree,
115  int n_max_iter);
116 
117 /*----------------------------------------------------------------------------
118  * Create iterative sparse linear system solver info and context.
119  *
120  * parameters:
121  * solver_type <-- type of solver (PCG, Jacobi, ...)
122  * poly_degree <-- preconditioning polynomial degree
123  * (0: diagonal; -1: non-preconditioned)
124  * n_max_iter <-- maximum number of iterations
125  * update_stats <-- automatic solver statistics indicator
126  *
127  * returns:
128  * pointer to newly created solver info object.
129  *----------------------------------------------------------------------------*/
130 
131 cs_sles_it_t *
133  int poly_degree,
134  int n_max_iter,
135  bool update_stats);
136 
137 /*----------------------------------------------------------------------------
138  * Destroy iterative sparse linear system solver info and context.
139  *
140  * parameters:
141  * context <-> pointer to iterative sparse linear solver info
142  * (actual type: cs_sles_it_t **)
143  *----------------------------------------------------------------------------*/
144 
145 void
146 cs_sles_it_destroy(void **context);
147 
148 /*----------------------------------------------------------------------------
149  * Create iterative sparse linear system solver info and context
150  * based on existing info and context.
151  *
152  * parameters:
153  * context <-- pointer to reference info and context
154  * (actual type: cs_sles_it_t *)
155  *
156  * returns:
157  * pointer to newly created solver info object
158  * (actual type: cs_sles_it_t *)
159  *----------------------------------------------------------------------------*/
160 
161 void *
162 cs_sles_it_copy(const void *context);
163 
164 /*----------------------------------------------------------------------------
165  * Setup iterative sparse linear equation solver.
166  *
167  * parameters:
168  * context <-> pointer to iterative sparse linear solver info
169  * (actual type: cs_sles_it_t *)
170  * name <-- pointer to system name
171  * a <-- associated matrix
172  * verbosity <-- verbosity level
173  *----------------------------------------------------------------------------*/
174 
175 void
176 cs_sles_it_setup(void *context,
177  const char *name,
178  const cs_matrix_t *a,
179  int verbosity);
180 
181 /*----------------------------------------------------------------------------
182  * Call iterative sparse linear equation solver.
183  *
184  * parameters:
185  * context <-> pointer to iterative sparse linear solver info
186  * (actual type: cs_sles_it_t *)
187  * name <-- pointer to system name
188  * a <-- matrix
189  * verbosity <-- verbosity level
190  * rotation_mode <-- halo update option for rotational periodicity
191  * precision <-- solver precision
192  * r_norm <-- residue normalization
193  * n_iter --> number of iterations
194  * residue --> residue
195  * rhs <-- right hand side
196  * vx <-> system solution
197  * aux_size <-- number of elements in aux_vectors (in bytes)
198  * aux_vectors --- optional working area (internal allocation if NULL)
199  *
200  * returns:
201  * convergence state
202  *----------------------------------------------------------------------------*/
203 
205 cs_sles_it_solve(void *context,
206  const char *name,
207  const cs_matrix_t *a,
208  int verbosity,
209  cs_halo_rotation_t rotation_mode,
210  double precision,
211  double r_norm,
212  int *n_iter,
213  double *residue,
214  const cs_real_t *rhs,
215  cs_real_t *vx,
216  size_t aux_size,
217  void *aux_vectors);
218 
219 /*----------------------------------------------------------------------------
220  * Free iterative sparse linear equation solver setup context.
221  *
222  * This function frees resolution-related data, such as
223  * buffers and preconditioning but does not free the whole context,
224  * as info used for logging (especially performance data) is maintained.
225 
226  * parameters:
227  * context <-> pointer to iterative sparse linear solver info
228  * (actual type: cs_sles_it_t *)
229  *----------------------------------------------------------------------------*/
230 
231 void
232 cs_sles_it_free(void *context);
233 
234 /*----------------------------------------------------------------------------
235  * Log sparse linear equation solver info.
236  *
237  * parameters:
238  * context <-> pointer to iterative sparse linear solver info
239  * (actual type: cs_sles_it_t *)
240  * log_type <-- log type
241  *----------------------------------------------------------------------------*/
242 
243 void
244 cs_sles_it_log(const void *context,
245  cs_log_t log_type);
246 
247 /*----------------------------------------------------------------------------
248  * Return iterative solver type.
249  *
250  * parameters:
251  * context <-- pointer to iterative solver info and context
252  *
253  * returns:
254  * selected solver type
255  *----------------------------------------------------------------------------*/
256 
258 cs_sles_it_get_type(const cs_sles_it_t *context);
259 
260 /*----------------------------------------------------------------------------
261  * Return the initial residue for the previous solve operation with a solver.
262  *
263  * This is useful for convergence tests when this solver is used as
264  * a preconditioning smoother.
265  *
266  * This operation is only valid between calls to cs_sles_it_setup()
267  * (or cs_sles_it_solve()) and cs_sles_it_free().
268  * It returns -1 otherwise.
269  *
270  * parameters:
271  * context <-- pointer to iterative solver info and context
272  *
273  * returns:
274  * initial residue from last call to \ref cs_sles_solve with this solver
275  *----------------------------------------------------------------------------*/
276 
277 double
279 
280 /*----------------------------------------------------------------------------
281  * Associate a similar info and context object with which some setup
282  * data may be shared.
283  *
284  * This is especially useful for sharing preconditioning data between
285  * similar solver contexts (for example ascending and descending multigrid
286  * smoothers based on the same matrix).
287  *
288  * For preconditioning data to be effectively shared, cs_sles_it_setup()
289  * (or cs_sles_it_solve()) must be called on "shareable" before being
290  * called on "context" (without cs_sles_it_free() being called in between,
291  * of course).
292  *
293  * It is the caller's responsibility to ensure the context is not used
294  * for a cs_sles_it_setup() or cs_sles_it_solve() operation after the
295  * shareable object has been destroyed (normally by cs_sles_it_destroy()).
296  *
297  * parameters:
298  * context <-> pointer to iterative sparse linear system solver info
299  * shareable <-- pointer to iterative solver info and context
300  * whose context may be shared
301  *----------------------------------------------------------------------------*/
302 
303 void
305  const cs_sles_it_t *shareable);
306 
307 #if defined(HAVE_MPI)
308 
309 /*----------------------------------------------------------------------------
310  * Set MPI communicator for dot products.
311  *
312  * parameters:
313  * context <-> pointer to iterative sparse linear system solver info
314  * comm <-- MPI communicator
315  *----------------------------------------------------------------------------*/
316 
317 void
318 cs_sles_it_set_mpi_reduce_comm(cs_sles_it_t *context,
319  MPI_Comm comm);
320 
321 #endif /* defined(HAVE_MPI) */
322 
323 /*----------------------------------------------------------------------------
324  * Assign ordering to iterative solver.
325  *
326  * The solver context takes ownership of the order array (i.e. it will
327  * handle its later deallocation).
328  *
329  * This is useful only for Block Gauss-Seidel.
330  *
331  * parameters:
332  * context <-> pointer to iterative solver info and context
333  * order <-> pointer to ordering array
334  *----------------------------------------------------------------------------*/
335 
336 void
338  cs_lnum_t **order);
339 
340 /*----------------------------------------------------------------------------
341  * Query mean number of rows under which Conjugate Gradient algorithm
342  * uses the single-reduction variant.
343  *
344  * The single-reduction variant requires only one parallel sum per
345  * iteration (instead of 2), at the cost of additional vector operations,
346  * so it tends to be more expensive when the number of matrix rows per
347  * MPI rank is high, then becomes cheaper when the MPI latency cost becomes
348  * more significant.
349  *
350  * This option is ignored for non-parallel runs, so 0 is returned.
351  *
352  * return:
353  * mean number of rows per active rank under which the
354  * single-reduction variant will be used
355  *----------------------------------------------------------------------------*/
356 
357 cs_lnum_t
359 
360 /*----------------------------------------------------------------------------
361  * Set mean number of rows under which Conjugate Gradient algorithm
362  * should use the single-reduction variant.
363  *
364  * The single-reduction variant requires only one parallel sum per
365  * iteration (instead of 2), at the cost of additional vector operations,
366  * so it tends to be more expensive when the number of matrix rows per
367  * MPI rank is high, then becomes cheaper when the MPI latency cost becomes
368  * more significant.
369  *
370  * This option is ignored for non-parallel runs.
371  *
372  * parameters:
373  * threshold <-- mean number of rows per active rank under which the
374  * single-reduction variant will be used
375  *----------------------------------------------------------------------------*/
376 
377 void
379 
380 /*----------------------------------------------------------------------------
381  * Log the current global settings relative to parallelism.
382  *----------------------------------------------------------------------------*/
383 
384 void
386 
387 /*----------------------------------------------------------------------------
388  * Error handler for iterative sparse linear equation solver.
389  *
390  * In case of divergence or breakdown, this error handler outputs
391  * postprocessing data to assist debugging, then aborts the run.
392  * It does nothing in case the maximum iteration count is reached.
393  *
394  * parameters:
395  * context <-> pointer to iterative sparse linear system solver info
396  * (actual type: cs_sles_it_t *)
397  * state <-- convergence state
398  * name <-- pointer to name of linear system
399  * a <-- matrix
400  * rotation_mode <-- halo update option for rotational periodicity
401  * rhs <-- right hand side
402  * vx <-> system solution
403  */
404 /*----------------------------------------------------------------------------*/
405 
406 void
407 cs_sles_it_error_post_and_abort(void *context,
409  const char *name,
410  const cs_matrix_t *a,
411  cs_halo_rotation_t rotation_mode,
412  const cs_real_t *rhs,
413  cs_real_t *vx);
414 
415 /*----------------------------------------------------------------------------*/
416 
418 
419 #endif /* __CS_SLES_IT_H__ */
void cs_sles_it_setup(void *context, const char *name, const cs_matrix_t *a, int verbosity)
Setup iterative sparse linear equation solver.
Definition: cs_sles_it.c:3342
cs_halo_rotation_t
Definition: cs_halo.h:59
cs_sles_it_type_t
Iterative solver types.
Definition: cs_sles_it.h:54
Definition: cs_sles_it.h:60
void cs_sles_it_error_post_and_abort(void *context, cs_sles_convergence_state_t state, const char *name, const cs_matrix_t *a, cs_halo_rotation_t rotation_mode, const cs_real_t *rhs, cs_real_t *vx)
Error handler for iterative sparse linear equation solver.
Definition: cs_sles_it.c:3906
double cs_sles_it_get_last_initial_residue(const cs_sles_it_t *context)
Return the initial residue for the previous solve operation with a solver.
Definition: cs_sles_it.c:3699
#define BEGIN_C_DECLS
Definition: cs_defs.h:419
cs_lnum_t cs_sles_it_get_pcg_single_reduction(void)
Query mean number of rows under which Conjugate Gradient algorithm uses the single-reduction variant...
Definition: cs_sles_it.c:3832
struct _cs_sles_it_t cs_sles_it_t
Definition: cs_sles_it.h:68
struct _cs_matrix_t cs_matrix_t
Definition: cs_matrix.h:86
void cs_sles_it_destroy(void **context)
Destroy iterative sparse linear system solver info and context.
Definition: cs_sles_it.c:3219
void cs_sles_it_free(void *context)
Free iterative sparse linear equation solver setup context.
Definition: cs_sles_it.c:3644
void cs_sles_it_log(const void *context, cs_log_t log_type)
Log sparse linear equation solver info.
Definition: cs_sles_it.c:3276
cs_sles_convergence_state_t
Convergence status indicator.
Definition: cs_sles.h:55
Definition: cs_sles_it.h:62
Definition: cs_sles_it.h:59
Definition: cs_sles_it.h:56
const char * cs_sles_it_type_name[]
cs_sles_convergence_state_t cs_sles_it_solve(void *context, const char *name, const cs_matrix_t *a, int verbosity, cs_halo_rotation_t rotation_mode, double precision, double r_norm, int *n_iter, double *residue, const cs_real_t *rhs, cs_real_t *vx, size_t aux_size, void *aux_vectors)
Call iterative sparse linear equation solver.
Definition: cs_sles_it.c:3383
cs_sles_it_type_t cs_sles_it_get_type(const cs_sles_it_t *context)
Return iterative solver type.
Definition: cs_sles_it.c:3674
cs_log_t
Definition: cs_log.h:47
void cs_sles_it_log_parallel_options(void)
Log the current global settings relative to parallelism.
Definition: cs_sles_it.c:3874
void cs_sles_it_set_pcg_single_reduction(cs_lnum_t threshold)
Set mean number of rows under which Conjugate Gradient algorithm should use the single-reduction vari...
Definition: cs_sles_it.c:3860
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:292
Definition: cs_sles_it.h:58
void cs_sles_it_set_shareable(cs_sles_it_t *context, const cs_sles_it_t *shareable)
Associate a similar info and context object with which some setup data may be shared.
Definition: cs_sles_it.c:3733
#define END_C_DECLS
Definition: cs_defs.h:420
double cs_real_t
Definition: cs_defs.h:296
Definition: cs_sles_it.h:61
void cs_sles_it_assign_order(cs_sles_it_t *context, cs_lnum_t **order)
Assign ordering to iterative solver.
Definition: cs_sles_it.c:3791
cs_sles_it_t * cs_sles_it_create(cs_sles_it_type_t solver_type, int poly_degree, int n_max_iter, bool update_stats)
Create iterative sparse linear system solver info and context.
Definition: cs_sles_it.c:3168
cs_sles_it_t * cs_sles_it_define(int f_id, const char *name, cs_sles_it_type_t solver_type, int poly_degree, int n_max_iter)
Define and associate an iterative sparse linear system solver for a given field or equation name...
Definition: cs_sles_it.c:3123
void * cs_sles_it_copy(const void *context)
Create iterative sparse linear system solver info and context based on existing info and context...
Definition: cs_sles_it.c:3247
Definition: cs_sles_it.h:57