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