programmer's 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-2018 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 
59  CS_SLES_PCG, /* Preconditionned conjugate gradient */
60  CS_SLES_IPCG, /* Inexact preconditionned conjugate gradient */
61  CS_SLES_JACOBI, /* Jacobi */
62  CS_SLES_BICGSTAB, /* Bi-conjugate gradient stabilized */
63  CS_SLES_BICGSTAB2, /* Bi-conjugate gradient stabilized - 2*/
64  CS_SLES_GMRES, /* Generalized minimal residual */
65  CS_SLES_P_GAUSS_SEIDEL, /* Process-local Gauss-Seidel */
66  CS_SLES_P_SYM_GAUSS_SEIDEL, /* Process-local symmetric Gauss-Seidel */
67  CS_SLES_PCR3, /* 3-layer conjugate residual */
68  CS_SLES_N_IT_TYPES /* Number of resolution algorithms */
69 
71 
72 /* Iterative linear solver context (opaque) */
73 
74 typedef struct _cs_sles_it_t cs_sles_it_t;
75 
76 /*============================================================================
77  * Global variables
78  *============================================================================*/
79 
80 /* Short names for matrix types */
81 
82 extern const char *cs_sles_it_type_name[];
83 
84 /*=============================================================================
85  * Public function prototypes
86  *============================================================================*/
87 
88 /*----------------------------------------------------------------------------
89  * Define and associate an iterative sparse linear system solver
90  * for a given field or equation name.
91  *
92  * If this system did not previously exist, it is added to the list of
93  * "known" systems. Otherwise, its definition is replaced by the one
94  * defined here.
95  *
96  * This is a utility function: if finer control is needed, see
97  * cs_sles_define() and cs_sles_it_create().
98  *
99  * Note that this function returns a pointer directly to the iterative solver
100  * management structure. This may be used to set further options,
101  * for example using cs_sles_set_plot_options(). If needed, cs_sles_find()
102  * may be used to obtain a pointer to the matching cs_sles_t container.
103  *
104  * parameters:
105  * f_id <-- associated field id, or < 0
106  * name <-- associated name if f_id < 0, or NULL
107  * solver_type <-- type of solver (PCG, Jacobi, ...)
108  * poly_degree <-- preconditioning polynomial degree
109  * (0: diagonal; -1: non-preconditioned)
110  * n_max_iter <-- maximum number of iterations
111  *
112  * returns:
113  * pointer to newly created iterative solver info object.
114  *----------------------------------------------------------------------------*/
115 
116 cs_sles_it_t *
117 cs_sles_it_define(int f_id,
118  const char *name,
119  cs_sles_it_type_t solver_type,
120  int poly_degree,
121  int n_max_iter);
122 
123 /*----------------------------------------------------------------------------
124  * Create iterative sparse linear system solver info and context.
125  *
126  * parameters:
127  * solver_type <-- type of solver (PCG, Jacobi, ...)
128  * poly_degree <-- preconditioning polynomial degree
129  * (0: diagonal; -1: non-preconditioned)
130  * n_max_iter <-- maximum number of iterations
131  * update_stats <-- automatic solver statistics indicator
132  *
133  * returns:
134  * pointer to newly created solver info object.
135  *----------------------------------------------------------------------------*/
136 
137 cs_sles_it_t *
139  int poly_degree,
140  int n_max_iter,
141  bool update_stats);
142 
143 /*----------------------------------------------------------------------------
144  * Destroy iterative sparse linear system solver info and context.
145  *
146  * parameters:
147  * context <-> pointer to iterative sparse linear solver info
148  * (actual type: cs_sles_it_t **)
149  *----------------------------------------------------------------------------*/
150 
151 void
152 cs_sles_it_destroy(void **context);
153 
154 /*----------------------------------------------------------------------------
155  * Create iterative sparse linear system solver info and context
156  * based on existing info and context.
157  *
158  * parameters:
159  * context <-- pointer to reference info and context
160  * (actual type: cs_sles_it_t *)
161  *
162  * returns:
163  * pointer to newly created solver info object
164  * (actual type: cs_sles_it_t *)
165  *----------------------------------------------------------------------------*/
166 
167 void *
168 cs_sles_it_copy(const void *context);
169 
170 /*----------------------------------------------------------------------------
171  * Setup iterative sparse linear equation solver.
172  *
173  * parameters:
174  * context <-> pointer to iterative sparse linear solver info
175  * (actual type: cs_sles_it_t *)
176  * name <-- pointer to system name
177  * a <-- associated matrix
178  * verbosity <-- verbosity level
179  *----------------------------------------------------------------------------*/
180 
181 void
182 cs_sles_it_setup(void *context,
183  const char *name,
184  const cs_matrix_t *a,
185  int verbosity);
186 
187 /*----------------------------------------------------------------------------
188  * Call iterative sparse linear equation solver.
189  *
190  * parameters:
191  * context <-> pointer to iterative sparse linear solver info
192  * (actual type: cs_sles_it_t *)
193  * name <-- pointer to system name
194  * a <-- matrix
195  * verbosity <-- verbosity level
196  * rotation_mode <-- halo update option for rotational periodicity
197  * precision <-- solver precision
198  * r_norm <-- residue normalization
199  * n_iter --> number of iterations
200  * residue --> residue
201  * rhs <-- right hand side
202  * vx <-> system solution
203  * aux_size <-- number of elements in aux_vectors (in bytes)
204  * aux_vectors --- optional working area (internal allocation if NULL)
205  *
206  * returns:
207  * convergence state
208  *----------------------------------------------------------------------------*/
209 
211 cs_sles_it_solve(void *context,
212  const char *name,
213  const cs_matrix_t *a,
214  int verbosity,
215  cs_halo_rotation_t rotation_mode,
216  double precision,
217  double r_norm,
218  int *n_iter,
219  double *residue,
220  const cs_real_t *rhs,
221  cs_real_t *vx,
222  size_t aux_size,
223  void *aux_vectors);
224 
225 /*----------------------------------------------------------------------------
226  * Free iterative sparse linear equation solver setup context.
227  *
228  * This function frees resolution-related data, such as
229  * buffers and preconditioning but does not free the whole context,
230  * as info used for logging (especially performance data) is maintained.
231 
232  * parameters:
233  * context <-> pointer to iterative sparse linear solver info
234  * (actual type: cs_sles_it_t *)
235  *----------------------------------------------------------------------------*/
236 
237 void
238 cs_sles_it_free(void *context);
239 
240 /*----------------------------------------------------------------------------
241  * Log sparse linear equation solver info.
242  *
243  * parameters:
244  * context <-> pointer to iterative sparse linear solver info
245  * (actual type: cs_sles_it_t *)
246  * log_type <-- log type
247  *----------------------------------------------------------------------------*/
248 
249 void
250 cs_sles_it_log(const void *context,
251  cs_log_t log_type);
252 
253 /*----------------------------------------------------------------------------
254  * Return iterative solver type.
255  *
256  * parameters:
257  * context <-- pointer to iterative solver info and context
258  *
259  * returns:
260  * selected solver type
261  *----------------------------------------------------------------------------*/
262 
264 cs_sles_it_get_type(const cs_sles_it_t *context);
265 
266 /*----------------------------------------------------------------------------
267  * Return the initial residue for the previous solve operation with a solver.
268  *
269  * This is useful for convergence tests when this solver is used as
270  * a preconditioning smoother.
271  *
272  * This operation is only valid between calls to cs_sles_it_setup()
273  * (or cs_sles_it_solve()) and cs_sles_it_free().
274  * It returns -1 otherwise.
275  *
276  * parameters:
277  * context <-- pointer to iterative solver info and context
278  *
279  * returns:
280  * initial residue from last call to \ref cs_sles_solve with this solver
281  *----------------------------------------------------------------------------*/
282 
283 double
285 
286 /*----------------------------------------------------------------------------
287  * Return a preconditioner context for an iterative sparse linear
288  * equation solver.
289  *
290  * This allows modifying parameters of a non default (Jacobi or polynomial)
291  * preconditioner.
292  *
293  * parameters:
294  * context <-- pointer to iterative solver info and context
295  *
296  * returns:
297  * pointer to preconditoner context
298  *----------------------------------------------------------------------------*/
299 
300 cs_sles_pc_t *
302 
303 /*----------------------------------------------------------------------------
304  * Assign a preconditioner to an iterative sparse linear equation
305  * solver, transfering its ownership to to solver context.
306  *
307  * This allows assigning a non default (Jacobi or polynomial) preconditioner.
308  *
309  * The input pointer is set to NULL to make it clear the caller does not
310  * own the preconditioner anymore, though the context can be accessed using
311  * cs_sles_it_get_cp().
312  *
313  * parameters:
314  * context <-> pointer to iterative solver info and context
315  * pc <-> pointer to preconditoner context
316  *----------------------------------------------------------------------------*/
317 
318 void
320  cs_sles_pc_t **pc);
321 
322 /*----------------------------------------------------------------------------
323  * Copy options from one iterative sparse linear system solver info
324  * and context to another.
325  *
326  * Optional plotting contexts are shared between the source and destination
327  * contexts.
328  *
329  * Preconditioner settings are to be handled separately.
330  *
331  * parameters:
332  * src <-- pointer to source info and context
333  * dest <-> pointer to destination info and context
334  *----------------------------------------------------------------------------*/
335 
336 void
338  cs_sles_it_t *dest);
339 
340 /*----------------------------------------------------------------------------
341  * Associate a similar info and context object with which some setup
342  * data may be shared.
343  *
344  * This is especially useful for sharing preconditioning data between
345  * similar solver contexts (for example ascending and descending multigrid
346  * smoothers based on the same matrix).
347  *
348  * For preconditioning data to be effectively shared, cs_sles_it_setup()
349  * (or cs_sles_it_solve()) must be called on "shareable" before being
350  * called on "context" (without cs_sles_it_free() being called in between,
351  * of course).
352  *
353  * It is the caller's responsibility to ensure the context is not used
354  * for a cs_sles_it_setup() or cs_sles_it_solve() operation after the
355  * shareable object has been destroyed (normally by cs_sles_it_destroy()).
356  *
357  * parameters:
358  * context <-> pointer to iterative sparse linear system solver info
359  * shareable <-- pointer to iterative solver info and context
360  * whose context may be shared
361  *----------------------------------------------------------------------------*/
362 
363 void
365  const cs_sles_it_t *shareable);
366 
367 #if defined(HAVE_MPI)
368 
369 /*----------------------------------------------------------------------------
370  * Set MPI communicator for dot products.
371  *
372  * parameters:
373  * context <-> pointer to iterative sparse linear system solver info
374  * comm <-- MPI communicator
375  *----------------------------------------------------------------------------*/
376 
377 void
378 cs_sles_it_set_mpi_reduce_comm(cs_sles_it_t *context,
379  MPI_Comm comm);
380 
381 #endif /* defined(HAVE_MPI) */
382 
383 /*----------------------------------------------------------------------------
384  * Assign ordering to iterative solver.
385  *
386  * The solver context takes ownership of the order array (i.e. it will
387  * handle its later deallocation).
388  *
389  * This is useful only for Block Gauss-Seidel.
390  *
391  * parameters:
392  * context <-> pointer to iterative solver info and context
393  * order <-> pointer to ordering array
394  *----------------------------------------------------------------------------*/
395 
396 void
398  cs_lnum_t **order);
399 
400 /*----------------------------------------------------------------------------*/
414 /*----------------------------------------------------------------------------*/
415 
416 void
418  cs_sles_convergence_state_t threshold);
419 
420 /*----------------------------------------------------------------------------
421  * Query mean number of rows under which Conjugate Gradient algorithm
422  * uses the single-reduction variant.
423  *
424  * The single-reduction variant requires only one parallel sum per
425  * iteration (instead of 2), at the cost of additional vector operations,
426  * so it tends to be more expensive when the number of matrix rows per
427  * MPI rank is high, then becomes cheaper when the MPI latency cost becomes
428  * more significant.
429  *
430  * This option is ignored for non-parallel runs, so 0 is returned.
431  *
432  * return:
433  * mean number of rows per active rank under which the
434  * single-reduction variant will be used
435  *----------------------------------------------------------------------------*/
436 
437 cs_lnum_t
439 
440 /*----------------------------------------------------------------------------
441  * Set mean number of rows under which Conjugate Gradient algorithm
442  * should use the single-reduction variant.
443  *
444  * The single-reduction variant requires only one parallel sum per
445  * iteration (instead of 2), at the cost of additional vector operations,
446  * so it tends to be more expensive when the number of matrix rows per
447  * MPI rank is high, then becomes cheaper when the MPI latency cost becomes
448  * more significant.
449  *
450  * This option is ignored for non-parallel runs.
451  *
452  * parameters:
453  * threshold <-- mean number of rows per active rank under which the
454  * single-reduction variant will be used
455  *----------------------------------------------------------------------------*/
456 
457 void
459 
460 /*----------------------------------------------------------------------------
461  * Log the current global settings relative to parallelism.
462  *----------------------------------------------------------------------------*/
463 
464 void
466 
467 /*----------------------------------------------------------------------------
468  * Error handler for iterative sparse linear equation solver.
469  *
470  * In case of divergence or breakdown, this error handler outputs
471  * postprocessing data to assist debugging, then aborts the run.
472  * It does nothing in case the maximum iteration count is reached.
473  *
474  * parameters:
475  * sles <-> pointer to solver object
476  * state <-- convergence state
477  * a <-- matrix
478  * rotation_mode <-- halo update option for rotational periodicity
479  * rhs <-- right hand side
480  * vx <-> system solution
481  *
482  * returns:
483  * false (do not attempt new solve)
484  *----------------------------------------------------------------------------*/
485 
486 bool
489  const cs_matrix_t *a,
490  cs_halo_rotation_t rotation_mode,
491  const cs_real_t *rhs,
492  cs_real_t *vx);
493 
494 /*----------------------------------------------------------------------------
495  * Set plotting options for an iterative sparse linear equation solver.
496  *
497  * parameters:
498  * context <-> pointer to iterative solver info and context
499  * base_name <-- base plot name to activate, NULL otherwise
500  * use_iteration <-- if true, use iteration as time stamp
501  * otherwise, use wall clock time
502  *----------------------------------------------------------------------------*/
503 
504 void
506  const char *base_name,
507  bool use_iteration);
508 
509 /*----------------------------------------------------------------------------
510  * Assign existing time plot to iterative sparse linear equation solver.
511  *
512  * This is useful mainly when a time plot has a longer lifecycle than
513  * the linear solver context, such as inside a multigrid solver.
514  *
515  * parameters:
516  * context <-> pointer to iterative solver info and context
517  * time_plot <-- pointer to time plot structure
518  * time_stamp <-- associated time stamp
519  *----------------------------------------------------------------------------*/
520 
521 void
523  cs_time_plot_t *time_plot,
524  int time_stamp);
525 
526 /*----------------------------------------------------------------------------*/
527 
529 
530 #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:4557
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:5080
cs_halo_rotation_t
Definition: cs_halo.h:60
cs_sles_it_type_t
Iterative solver types.
Definition: cs_sles_it.h:57
Definition: cs_sles_it.h:64
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:4995
#define BEGIN_C_DECLS
Definition: cs_defs.h:453
Definition: cs_sles_it.h:60
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:5227
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:5048
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:5252
struct _cs_sles_it_t cs_sles_it_t
Definition: cs_sles_it.h:74
double cs_real_t
Floating-point value.
Definition: cs_defs.h:297
struct _cs_matrix_t cs_matrix_t
Definition: cs_matrix.h:90
void cs_sles_it_destroy(void **context)
Destroy iterative sparse linear system solver info and context.
Definition: cs_sles_it.c:4377
struct _cs_time_plot_t cs_time_plot_t
Definition: cs_time_plot.h:48
void cs_sles_it_free(void *context)
Free iterative sparse linear equation solver setup context.
Definition: cs_sles_it.c:4934
void cs_sles_it_log(const void *context, cs_log_t log_type)
Log sparse linear equation solver info.
Definition: cs_sles_it.c:4453
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
Definition: cs_sles_it.h:68
void cs_sles_it_assign_plot(cs_sles_it_t *context, cs_time_plot_t *time_plot, int time_stamp)
Assign existing time plot to iterative sparse linear equation solver.
Definition: cs_sles_it.c:5418
Definition: cs_sles_it.h:63
Definition: cs_sles_it.h:59
double precision, save a
Definition: cs_fuel_incl.f90:146
Definition: cs_sles_it.h:67
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:5019
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:4609
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:5326
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:4970
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:5294
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:5374
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:5280
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:293
Definition: cs_sles_it.h:66
Definition: cs_sles_it.h:62
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:5125
#define END_C_DECLS
Definition: cs_defs.h:454
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:5188
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:4287
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:4232
Definition: cs_sles_it.h:65
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:4415
Definition: cs_sles_it.h:61