7.0
general documentation
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-2021 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 #include "cs_time_plot.h"
38 #include "cs_sles.h"
39 #include "cs_sles_pc.h"
40 
41 /*----------------------------------------------------------------------------*/
42 
44 
45 /*============================================================================
46  * Macro definitions
47  *============================================================================*/
48 
49 /*============================================================================
50  * Type definitions
51  *============================================================================*/
52 
53 /*----------------------------------------------------------------------------
54  * Solver types
55  *----------------------------------------------------------------------------*/
56 
57 typedef enum {
58 
83 
84 /* Iterative linear solver context (opaque) */
85 
86 typedef struct _cs_sles_it_t cs_sles_it_t;
87 
88 /* Forward type declarations */
89 
90 typedef struct _cs_sles_it_convergence_t cs_sles_it_convergence_t;
91 
92 /*============================================================================
93  * Global variables
94  *============================================================================*/
95 
96 /* Names for iterative solver types */
97 
98 extern const char *cs_sles_it_type_name[];
99 
100 /*=============================================================================
101  * Public function prototypes
102  *============================================================================*/
103 
104 /*----------------------------------------------------------------------------
105  * Define and associate an iterative sparse linear system solver
106  * for a given field or equation name.
107  *
108  * If this system did not previously exist, it is added to the list of
109  * "known" systems. Otherwise, its definition is replaced by the one
110  * defined here.
111  *
112  * This is a utility function: if finer control is needed, see
113  * cs_sles_define() and cs_sles_it_create().
114  *
115  * Note that this function returns a pointer directly to the iterative solver
116  * management structure. This may be used to set further options,
117  * for example using cs_sles_set_plot_options(). If needed, cs_sles_find()
118  * may be used to obtain a pointer to the matching cs_sles_t container.
119  *
120  * parameters:
121  * f_id <-- associated field id, or < 0
122  * name <-- associated name if f_id < 0, or NULL
123  * solver_type <-- type of solver (PCG, Jacobi, ...)
124  * poly_degree <-- preconditioning polynomial degree
125  * (0: diagonal; -1: non-preconditioned)
126  * n_max_iter <-- maximum number of iterations
127  *
128  * returns:
129  * pointer to newly created iterative solver info object.
130  *----------------------------------------------------------------------------*/
131 
132 cs_sles_it_t *
133 cs_sles_it_define(int f_id,
134  const char *name,
135  cs_sles_it_type_t solver_type,
136  int poly_degree,
137  int n_max_iter);
138 
139 /*----------------------------------------------------------------------------
140  * Create iterative sparse linear system solver info and context.
141  *
142  * parameters:
143  * solver_type <-- type of solver (PCG, Jacobi, ...)
144  * poly_degree <-- preconditioning polynomial degree
145  * (0: diagonal; -1: non-preconditioned)
146  * n_max_iter <-- maximum number of iterations
147  * update_stats <-- automatic solver statistics indicator
148  *
149  * returns:
150  * pointer to newly created solver info object.
151  *----------------------------------------------------------------------------*/
152 
153 cs_sles_it_t *
155  int poly_degree,
156  int n_max_iter,
157  bool update_stats);
158 
159 /*----------------------------------------------------------------------------
160  * Destroy iterative sparse linear system solver info and context.
161  *
162  * parameters:
163  * context <-> pointer to iterative sparse linear solver info
164  * (actual type: cs_sles_it_t **)
165  *----------------------------------------------------------------------------*/
166 
167 void
168 cs_sles_it_destroy(void **context);
169 
170 /*----------------------------------------------------------------------------
171  * Create iterative sparse linear system solver info and context
172  * based on existing info and context.
173  *
174  * parameters:
175  * context <-- pointer to reference info and context
176  * (actual type: cs_sles_it_t *)
177  *
178  * returns:
179  * pointer to newly created solver info object
180  * (actual type: cs_sles_it_t *)
181  *----------------------------------------------------------------------------*/
182 
183 void *
184 cs_sles_it_copy(const void *context);
185 
186 /*----------------------------------------------------------------------------
187  * Setup iterative sparse linear equation solver.
188  *
189  * parameters:
190  * context <-> pointer to iterative sparse linear solver info
191  * (actual type: cs_sles_it_t *)
192  * name <-- pointer to system name
193  * a <-- associated matrix
194  * verbosity <-- verbosity level
195  *----------------------------------------------------------------------------*/
196 
197 void
198 cs_sles_it_setup(void *context,
199  const char *name,
200  const cs_matrix_t *a,
201  int verbosity);
202 
203 /*----------------------------------------------------------------------------
204  * Call iterative sparse linear equation solver.
205  *
206  * parameters:
207  * context <-> pointer to iterative sparse linear solver info
208  * (actual type: cs_sles_it_t *)
209  * name <-- pointer to system name
210  * a <-- matrix
211  * verbosity <-- verbosity level
212  * rotation_mode <-- halo update option for rotational periodicity
213  * precision <-- solver precision
214  * r_norm <-- residue normalization
215  * n_iter --> number of iterations
216  * residue --> residue
217  * rhs <-- right hand side
218  * vx <-> system solution
219  * aux_size <-- number of elements in aux_vectors (in bytes)
220  * aux_vectors --- optional working area (internal allocation if NULL)
221  *
222  * returns:
223  * convergence state
224  *----------------------------------------------------------------------------*/
225 
227 cs_sles_it_solve(void *context,
228  const char *name,
229  const cs_matrix_t *a,
230  int verbosity,
231  cs_halo_rotation_t rotation_mode,
232  double precision,
233  double r_norm,
234  int *n_iter,
235  double *residue,
236  const cs_real_t *rhs,
237  cs_real_t *vx,
238  size_t aux_size,
239  void *aux_vectors);
240 
241 /*----------------------------------------------------------------------------
242  * Free iterative sparse linear equation solver setup context.
243  *
244  * This function frees resolution-related data, such as
245  * buffers and preconditioning but does not free the whole context,
246  * as info used for logging (especially performance data) is maintained.
247 
248  * parameters:
249  * context <-> pointer to iterative sparse linear solver info
250  * (actual type: cs_sles_it_t *)
251  *----------------------------------------------------------------------------*/
252 
253 void
254 cs_sles_it_free(void *context);
255 
256 /*----------------------------------------------------------------------------
257  * Log sparse linear equation solver info.
258  *
259  * parameters:
260  * context <-> pointer to iterative sparse linear solver info
261  * (actual type: cs_sles_it_t *)
262  * log_type <-- log type
263  *----------------------------------------------------------------------------*/
264 
265 void
266 cs_sles_it_log(const void *context,
267  cs_log_t log_type);
268 
269 /*----------------------------------------------------------------------------
270  * Return iterative solver type.
271  *
272  * parameters:
273  * context <-- pointer to iterative solver info and context
274  *
275  * returns:
276  * selected solver type
277  *----------------------------------------------------------------------------*/
278 
280 cs_sles_it_get_type(const cs_sles_it_t *context);
281 
282 /*----------------------------------------------------------------------------
283  * Return the initial residue for the previous solve operation with a solver.
284  *
285  * This is useful for convergence tests when this solver is used as
286  * a preconditioning smoother.
287  *
288  * This operation is only valid between calls to cs_sles_it_setup()
289  * (or cs_sles_it_solve()) and cs_sles_it_free().
290  * It returns -1 otherwise.
291  *
292  * parameters:
293  * context <-- pointer to iterative solver info and context
294  *
295  * returns:
296  * initial residue from last call to \ref cs_sles_solve with this solver
297  *----------------------------------------------------------------------------*/
298 
299 double
301 
302 /*----------------------------------------------------------------------------
303  * Return a preconditioner context for an iterative sparse linear
304  * equation solver.
305  *
306  * This allows modifying parameters of a non default (Jacobi or polynomial)
307  * preconditioner.
308  *
309  * parameters:
310  * context <-- pointer to iterative solver info and context
311  *
312  * returns:
313  * pointer to preconditoner context
314  *----------------------------------------------------------------------------*/
315 
316 cs_sles_pc_t *
318 
319 /*----------------------------------------------------------------------------
320  * Assign a preconditioner to an iterative sparse linear equation
321  * solver, transfering its ownership to to solver context.
322  *
323  * This allows assigning a non default (Jacobi or polynomial) preconditioner.
324  *
325  * The input pointer is set to NULL to make it clear the caller does not
326  * own the preconditioner anymore, though the context can be accessed using
327  * cs_sles_it_get_cp().
328  *
329  * parameters:
330  * context <-> pointer to iterative solver info and context
331  * pc <-> pointer to preconditoner context
332  *----------------------------------------------------------------------------*/
333 
334 void
336  cs_sles_pc_t **pc);
337 
338 /*----------------------------------------------------------------------------
339  * Copy options from one iterative sparse linear system solver info
340  * and context to another.
341  *
342  * Optional plotting contexts are shared between the source and destination
343  * contexts.
344  *
345  * Preconditioner settings are to be handled separately.
346  *
347  * parameters:
348  * src <-- pointer to source info and context
349  * dest <-> pointer to destination info and context
350  *----------------------------------------------------------------------------*/
351 
352 void
354  cs_sles_it_t *dest);
355 
356 /*----------------------------------------------------------------------------
357  * Associate a similar info and context object with which some setup
358  * data may be shared.
359  *
360  * This is especially useful for sharing preconditioning data between
361  * similar solver contexts (for example ascending and descending multigrid
362  * smoothers based on the same matrix).
363  *
364  * For preconditioning data to be effectively shared, cs_sles_it_setup()
365  * (or cs_sles_it_solve()) must be called on "shareable" before being
366  * called on "context" (without cs_sles_it_free() being called in between,
367  * of course).
368  *
369  * It is the caller's responsibility to ensure the context is not used
370  * for a cs_sles_it_setup() or cs_sles_it_solve() operation after the
371  * shareable object has been destroyed (normally by cs_sles_it_destroy()).
372  *
373  * parameters:
374  * context <-> pointer to iterative sparse linear system solver info
375  * shareable <-- pointer to iterative solver info and context
376  * whose context may be shared
377  *----------------------------------------------------------------------------*/
378 
379 void
381  const cs_sles_it_t *shareable);
382 
383 #if defined(HAVE_MPI)
384 
385 /*----------------------------------------------------------------------------*/
397 /*----------------------------------------------------------------------------*/
398 
399 void
400 cs_sles_it_set_mpi_reduce_comm(cs_sles_it_t *context,
401  MPI_Comm comm,
402  MPI_Comm caller_comm);
403 
404 #endif /* defined(HAVE_MPI) */
405 
406 /*----------------------------------------------------------------------------
407  * Assign ordering to iterative solver.
408  *
409  * The solver context takes ownership of the order array (i.e. it will
410  * handle its later deallocation).
411  *
412  * This is useful only for Block Gauss-Seidel.
413  *
414  * parameters:
415  * context <-> pointer to iterative solver info and context
416  * order <-> pointer to ordering array
417  *----------------------------------------------------------------------------*/
418 
419 void
421  cs_lnum_t **order);
422 
423 /*----------------------------------------------------------------------------*/
437 /*----------------------------------------------------------------------------*/
438 
439 void
441  cs_sles_convergence_state_t threshold);
442 
443 /*----------------------------------------------------------------------------
444  * Query mean number of rows under which Conjugate Gradient algorithm
445  * uses the single-reduction variant.
446  *
447  * The single-reduction variant requires only one parallel sum per
448  * iteration (instead of 2), at the cost of additional vector operations,
449  * so it tends to be more expensive when the number of matrix rows per
450  * MPI rank is high, then becomes cheaper when the MPI latency cost becomes
451  * more significant.
452  *
453  * This option is ignored for non-parallel runs, so 0 is returned.
454  *
455  * return:
456  * mean number of rows per active rank under which the
457  * single-reduction variant will be used
458  *----------------------------------------------------------------------------*/
459 
460 cs_lnum_t
462 
463 /*----------------------------------------------------------------------------
464  * Set mean number of rows under which Conjugate Gradient algorithm
465  * should use the single-reduction variant.
466  *
467  * The single-reduction variant requires only one parallel sum per
468  * iteration (instead of 2), at the cost of additional vector operations,
469  * so it tends to be more expensive when the number of matrix rows per
470  * MPI rank is high, then becomes cheaper when the MPI latency cost becomes
471  * more significant.
472  *
473  * This option is ignored for non-parallel runs.
474  *
475  * parameters:
476  * threshold <-- mean number of rows per active rank under which the
477  * single-reduction variant will be used
478  *----------------------------------------------------------------------------*/
479 
480 void
482 
483 /*----------------------------------------------------------------------------
484  * Log the current global settings relative to parallelism.
485  *----------------------------------------------------------------------------*/
486 
487 void
489 
490 /*----------------------------------------------------------------------------
491  * Error handler for iterative sparse linear equation solver.
492  *
493  * In case of divergence or breakdown, this error handler outputs
494  * postprocessing data to assist debugging, then aborts the run.
495  * It does nothing in case the maximum iteration count is reached.
496  *
497  * parameters:
498  * sles <-> pointer to solver object
499  * state <-- convergence state
500  * a <-- matrix
501  * rotation_mode <-- halo update option for rotational periodicity
502  * rhs <-- right hand side
503  * vx <-> system solution
504  *
505  * returns:
506  * false (do not attempt new solve)
507  *----------------------------------------------------------------------------*/
508 
509 bool
512  const cs_matrix_t *a,
513  cs_halo_rotation_t rotation_mode,
514  const cs_real_t *rhs,
515  cs_real_t *vx);
516 
517 /*----------------------------------------------------------------------------
518  * Set plotting options for an iterative sparse linear equation solver.
519  *
520  * parameters:
521  * context <-> pointer to iterative solver info and context
522  * base_name <-- base plot name to activate, NULL otherwise
523  * use_iteration <-- if true, use iteration as time stamp
524  * otherwise, use wall clock time
525  *----------------------------------------------------------------------------*/
526 
527 void
529  const char *base_name,
530  bool use_iteration);
531 
532 /*----------------------------------------------------------------------------*/
533 
535 
536 #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:3919
void cs_sles_it_transfer_parameters(const cs_sles_it_t *src, cs_sles_it_t *dest)
Copy options from one iterative sparse linear system solver info and context to another.
Definition: cs_sles_it.c:4382
cs_halo_rotation_t
Definition: cs_halo.h:60
Definition: cs_sles_it.h:79
cs_sles_it_type_t
Definition: cs_sles_it.h:57
Definition: cs_sles_it.h:71
Definition: cs_sles_it.h:76
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:4297
#define BEGIN_C_DECLS
Definition: cs_defs.h:495
Definition: cs_sles_it.h:66
void cs_sles_it_set_fallback_threshold(cs_sles_it_t *context, cs_sles_convergence_state_t threshold)
Define convergence level under which the fallback to another solver may be used if applicable...
Definition: cs_sles_it.c:4541
void cs_sles_it_transfer_pc(cs_sles_it_t *context, cs_sles_pc_t **pc)
Assign a preconditioner to an iterative sparse linear equation solver, transfering its ownership to t...
Definition: cs_sles_it.c:4350
struct _cs_sles_pc_t cs_sles_pc_t
Definition: cs_sles_pc.h:66
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:4566
struct _cs_sles_it_t cs_sles_it_t
Definition: cs_sles_it.h:86
Definition: cs_sles_it.h:60
double cs_real_t
Floating-point value.
Definition: cs_defs.h:307
struct _cs_matrix_t cs_matrix_t
Definition: cs_matrix.h:94
void cs_sles_it_destroy(void **context)
Destroy iterative sparse linear system solver info and context.
Definition: cs_sles_it.c:3739
void cs_sles_it_free(void *context)
Free iterative sparse linear equation solver setup context.
Definition: cs_sles_it.c:4236
void cs_sles_it_log(const void *context, cs_log_t log_type)
Log sparse linear equation solver info.
Definition: cs_sles_it.c:3815
Definition: cs_sles_it.h:69
cs_sles_convergence_state_t
Convergence status indicator.
Definition: cs_sles.h:56
struct _cs_sles_t cs_sles_t
Definition: cs_sles.h:68
struct _cs_sles_it_convergence_t cs_sles_it_convergence_t
Definition: cs_sles_it.h:90
Definition: cs_sles_it.h:73
Definition: cs_sles_it.h:63
double precision, save a
Definition: cs_fuel_incl.f90:146
cs_sles_pc_t * cs_sles_it_get_pc(cs_sles_it_t *context)
Return a preconditioner context for an iterative sparse linear equation solver.
Definition: cs_sles_it.c:4321
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:4068
bool cs_sles_it_error_post_and_abort(cs_sles_t *sles, cs_sles_convergence_state_t state, 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:4640
Definition: cs_sles_it.h:64
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:4272
Definition: cs_sles_it.h:67
cs_log_t
Definition: cs_log.h:48
void cs_sles_it_log_parallel_options(void)
Log the current global settings relative to parallelism.
Definition: cs_sles_it.c:4608
void cs_sles_it_set_plot_options(cs_sles_it_t *context, const char *base_name, bool use_iteration)
Set plotting options for an iterative sparse linear equation solver.
Definition: cs_sles_it.c:4688
Definition: cs_sles_it.h:59
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:4594
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:301
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:4427
#define END_C_DECLS
Definition: cs_defs.h:496
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:4502
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:3641
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:3586
Definition: cs_sles_it.h:62
Definition: cs_sles_it.h:70
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:3777
Definition: cs_sles_it.h:77