programmer's documentation
cs_sles.h
Go to the documentation of this file.
1 #ifndef __CS_SLES_H__
2 #define __CS_SLES_H__
3 
4 /*============================================================================
5  * Sparse Linear Equation Solver driver
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_log.h"
36 #include "cs_halo_perio.h"
37 #include "cs_matrix.h"
38 
39 /*----------------------------------------------------------------------------*/
40 
42 
43 /*============================================================================
44  * Macro definitions
45  *============================================================================*/
46 
47 /*============================================================================
48  * Type definitions
49  *============================================================================*/
50 
51 /*----------------------------------------------------------------------------
52  * Convergence status
53  *----------------------------------------------------------------------------*/
54 
55 typedef enum {
56 
62 
64 
65 /* General linear solver context (opaque) */
66 
67 typedef struct _cs_sles_t cs_sles_t;
68 
69 /*----------------------------------------------------------------------------
70  * Function pointer for pre-resolution setup of a linear system solvers's
71  * context.
72  *
73  * This setup may include building a multigrid hierarchy, or a preconditioner.
74  *
75  * Use of this type of function is optional: the context is expected to
76  * maintain state, so that if a cs_sles_solve_t function is called before a
77  * cs_sles_setup_t function, the latter will be called automatically.
78  *
79  * parameters:
80  * context <-> pointer to solver context
81  * name <-- pointer to name of linear system
82  * a <-- matrix
83  * verbosity <-- associated verbosity
84  *----------------------------------------------------------------------------*/
85 
86 typedef void
87 (cs_sles_setup_t) (void *context,
88  const char *name,
89  const cs_matrix_t *a,
90  int verbosity);
91 
92 /*----------------------------------------------------------------------------
93  * Function pointer for resolution of a linear system.
94  *
95  * If the associated cs_sles_setup_t function has not been called before
96  * this function, it will be called automatically.
97  *
98  * The solution context setup by this call (or that of the matching setup
99  * function) will be maintained until the matching cs_sles_free_t function
100  * is called.
101  *
102  * The matrix is not expected to change between successive calls, although
103  * the right hand side may. If the matrix changes, the associated
104  * cs_sles_setup_t or cs_sles_free_t function must be called between
105  * solves.
106  *
107  * The system is considered to have converged when
108  * residue/r_norm <= precision, residue being the L2 norm of a.vx-rhs.
109  *
110  * parameters:
111  * context <-> pointer to solver context
112  * name <-- pointer to name of linear system
113  * a <-- matrix
114  * verbosity <-- associated verbosity
115  * rotation_mode <-- halo update option for rotational periodicity
116  * precision <-- solver precision
117  * r_norm <-- residue normalization
118  * n_iter --> number of "equivalent" iterations
119  * residue --> residue
120  * rhs <-- right hand side
121  * vx <-- system solution
122  * aux_size <-- number of elements in aux_vectors
123  * aux_vectors <-- optional working area (internal allocation if NULL)
124  *
125  * returns:
126  * convergence status
127  *----------------------------------------------------------------------------*/
128 
130 (cs_sles_solve_t) (void *context,
131  const char *name,
132  const cs_matrix_t *a,
133  int verbosity,
134  cs_halo_rotation_t rotation_mode,
135  double precision,
136  double r_norm,
137  int *n_iter,
138  double *residue,
139  const cs_real_t *rhs,
140  cs_real_t *vx,
141  size_t aux_size,
142  void *aux_vectors);
143 
144 /*----------------------------------------------------------------------------
145  * Function pointer for freeing of a linear system's context data.
146  *
147  * Note that this function should free resolution-related data, such as
148  * multigrid hierarchy, preconditioning, and any other temporary arrays or
149  * objects required for resolution, but should not free the whole context,
150  * as info used for logging (especially performance data) should be
151  * maintained.
152  *
153  * parameters:
154  * context <-> pointer to solver context
155  *----------------------------------------------------------------------------*/
156 
157 typedef void
158 (cs_sles_free_t) (void *context);
159 
160 /*----------------------------------------------------------------------------
161  * Function pointer for logging of linear solver setup,
162  * history and performance data.
163  *
164  * This function will be called for each solver when cs_sles_finalize()
165  * is called.
166  *
167  * parameters:
168  * context <-- pointer to solver context
169  * log_type <-- log type
170  *----------------------------------------------------------------------------*/
171 
172 typedef void
173 (cs_sles_log_t) (const void *context,
174  cs_log_t log_type);
175 
176 /*----------------------------------------------------------------------------
177  * Function pointer for creation of a solver context based on the copy
178  * of another.
179  *
180  * The new context copies the settings of the copied context, but not
181  * its setup data and logged info, such as performance data.
182  *
183  * This type of function is optional, but enables associating different
184  * solvers to related systems (to differentiate logging) while using
185  * the same settings by default.
186  *
187  * parameters:
188  * context <-- source context
189  *
190  * returns:
191  * pointer to newly created context
192  *----------------------------------------------------------------------------*/
193 
194 typedef void *
195 (cs_sles_copy_t) (const void *context);
196 
197 /*----------------------------------------------------------------------------
198  * Function pointer for destruction of a linear system solver context.
199  *
200  * This function should free all context data, and will be called for each
201  * system when cs_sles_finalize() is called.
202  *
203  * parameters:
204  * context <-> pointer to solver context
205  *----------------------------------------------------------------------------*/
206 
207 typedef void
208 (cs_sles_destroy_t) (void **context);
209 
210 /*----------------------------------------------------------------------------
211  * Function pointer for handling of non-convegence when solving
212  * a linear system.
213  *
214  * Such a function is optional, and may be used for a variety of purposes,
215  * such as logging, postprocessing, re-trying with different parameters,
216  * aborting the run, or any combination thereof.
217  *
218  * An error handler may be associated with a given solver context using
219  * cs_sles_set_error_handler(), in which case it will be called whenever
220  * convergence fails.
221  *
222  * Remark: in the advent of re-trying with different parameters,
223  * it is recommended that either the parameters be sufficiently similar that
224  * performance logging will not be affected, so the info reported to
225  * the user is not biased. Complex strategies involving different
226  * solver types should not be based on the use of an error handler, but
227  * built-into the solver function, with appropriate performance logging
228  * (though they may use an error handler in case of final failure).
229  *
230  * parameters:
231  * context <-> pointer to solver context
232  * state <-- convergence state
233  * name <-- name of linear system
234  * a <-- matrix
235  * rotation_mode <-- Halo update option for rotational periodicity
236  * rhs <-- Right hand side
237  * vx <-- System solution
238  *----------------------------------------------------------------------------*/
239 
240 typedef void
241 (cs_sles_error_handler_t) (void *context,
243  const char *name,
244  const cs_matrix_t *a,
245  cs_halo_rotation_t rotation_mode,
246  const cs_real_t *rhs,
247  cs_real_t *vx);
248 
249 /*----------------------------------------------------------------------------
250  * Function pointer for the default definition of a sparse
251  * linear equation solver
252  *
253  * The function may be associated using cs_sles_set_default_define(), so
254  * that it may provide a definition that will be used when
255  * cs_sles_setup() or cs_sles_solve() is used for a system for which
256  * no matching call to cs_sles_define() has been done.
257  *
258  * The function should call cs_sles_define() with arguments f_id
259  * and name, and appropriately chosen function pointers.
260  *
261  * A pointer to the matrix of the system to be solved is also provided,
262  * so that the corresponding information may be used to better choose
263  * defaults.
264  *
265  * parameters:
266  * f_id <-- associated field id, or < 0
267  * name <-- associated name if f_id < 0, or NULL
268  * a <-- Matrix
269  *----------------------------------------------------------------------------*/
270 
271 typedef void
272 (cs_sles_define_t) (int f_id,
273  const char *name,
274  const cs_matrix_t *a);
275 
276 /*----------------------------------------------------------------------------
277  * Function pointer for the default definition of a sparse
278  * linear equation solver's verbosity
279  *
280  * The function may be associated using cs_sles_set_default_verbosity(), so
281  * that it may provide a definition that will be used when
282  * cs_sles_default_verbosity() is called.
283  *
284  * parameters:
285  * f_id <-- associated field id, or < 0
286  * name <-- associated name if f_id < 0, or NULL
287  *
288  * returns:
289  * default verbosity value
290  *----------------------------------------------------------------------------*/
291 
292 typedef int
294  const char *name);
295 
296 /*============================================================================
297  * Global variables
298  *============================================================================*/
299 
300 /*=============================================================================
301  * Public function prototypes for Fortran API
302  *============================================================================*/
303 
304 /*=============================================================================
305  * Public function prototypes
306  *============================================================================*/
307 
308 /*----------------------------------------------------------------------------
309  * \brief Initialize sparse linear equation solver API.
310  *----------------------------------------------------------------------------*/
311 
312 void
313 cs_sles_initialize(void);
314 
315 /*----------------------------------------------------------------------------
316  * \brief Finalize sparse linear equation solver API.
317  *----------------------------------------------------------------------------*/
318 
319 void
320 cs_sles_finalize(void);
321 
322 /*----------------------------------------------------------------------------*/
328 /*----------------------------------------------------------------------------*/
329 
330 void
331 cs_sles_log(cs_log_t log_type);
332 
333 /*----------------------------------------------------------------------------*/
345 /*----------------------------------------------------------------------------*/
346 
347 cs_sles_t *
348 cs_sles_find(int f_id,
349  const char *name);
350 
351 /*----------------------------------------------------------------------------*/
368 /*----------------------------------------------------------------------------*/
369 
370 cs_sles_t *
371 cs_sles_find_or_add(int f_id,
372  const char *name);
373 
374 /*----------------------------------------------------------------------------*/
390 /*----------------------------------------------------------------------------*/
391 
392 void
393 cs_sles_push(int f_id,
394  const char *name);
395 
396 /*----------------------------------------------------------------------------*/
404 /*----------------------------------------------------------------------------*/
405 
406 void
407 cs_sles_pop(int f_id);
408 
409 /*----------------------------------------------------------------------------*/
444 /*----------------------------------------------------------------------------*/
445 
446 cs_sles_t *
447 cs_sles_define(int f_id,
448  const char *name,
449  void *context,
450  const char *type_name,
451  cs_sles_setup_t *setup_func,
452  cs_sles_solve_t *solve_func,
453  cs_sles_free_t *free_func,
454  cs_sles_log_t *log_func,
455  cs_sles_copy_t *copy_func,
456  cs_sles_destroy_t *destroy_func);
457 
458 /*----------------------------------------------------------------------------*/
470 /*----------------------------------------------------------------------------*/
471 
472 void
474  int verbosity);
475 
476 /*----------------------------------------------------------------------------*/
495 /*----------------------------------------------------------------------------*/
496 
497 const char *
499 
500 /*----------------------------------------------------------------------------*/
512 /*----------------------------------------------------------------------------*/
513 
514 void *
516 
517 /*----------------------------------------------------------------------------*/
531 /*----------------------------------------------------------------------------*/
532 
533 void
535  const cs_matrix_t *a);
536 
537 /*----------------------------------------------------------------------------*/
568 /*----------------------------------------------------------------------------*/
569 
572  const cs_matrix_t *a,
573  cs_halo_rotation_t rotation_mode,
574  double precision,
575  double r_norm,
576  int *n_iter,
577  double *residue,
578  const cs_real_t *rhs,
579  cs_real_t *vx,
580  size_t aux_size,
581  void *aux_vectors);
582 
583 /*----------------------------------------------------------------------------*/
594 /*----------------------------------------------------------------------------*/
595 
596 void
597 cs_sles_free(cs_sles_t *sles);
598 
599 /*----------------------------------------------------------------------------*/
619 /*----------------------------------------------------------------------------*/
620 
621 int
622 cs_sles_copy(cs_sles_t *dest,
623  const cs_sles_t *src);
624 
625 /*----------------------------------------------------------------------------*/
640 /*----------------------------------------------------------------------------*/
641 
642 void
644  cs_sles_error_handler_t *error_handler_func);
645 
646 /*----------------------------------------------------------------------------*/
656 /*----------------------------------------------------------------------------*/
657 
658 void
660 
661 /*----------------------------------------------------------------------------*/
670 /*----------------------------------------------------------------------------*/
671 
672 void
674 
675 /*----------------------------------------------------------------------------*/
686 /*----------------------------------------------------------------------------*/
687 
688 void
689 cs_sles_post_error_output_def(const char *name,
690  int mesh_id,
691  cs_halo_rotation_t rotation_mode,
692  const cs_matrix_t *a,
693  const cs_real_t *rhs,
694  cs_real_t *vx);
695 
696 /*----------------------------------------------------------------------------*/
705 /*----------------------------------------------------------------------------*/
706 
707 void
708 cs_sles_post_error_output_var(const char *name,
709  int mesh_id,
710  int diag_block_size,
711  cs_real_t *var);
712 
713 /*----------------------------------------------------------------------------*/
725 /*----------------------------------------------------------------------------*/
726 
727 const char *
728 cs_sles_base_name(int f_id,
729  const char *name);
730 
731 /*----------------------------------------------------------------------------*/
732 
734 
735 #endif /* __CS_SLES_H__ */
cs_sles_t * cs_sles_find_or_add(int f_id, const char *name)
Return pointer to linear system object, based on matching field id or system name.
Definition: cs_sles.c:1040
void cs_sles_log(cs_log_t log_type)
Log sparse linear equation solver info.
Definition: cs_sles.c:836
void *( cs_sles_copy_t)(const void *context)
Function pointer for creation of a solver context based on the copy of another.
Definition: cs_sles.h:195
const char * cs_sles_get_type(cs_sles_t *sles)
Return type name of solver context.
Definition: cs_sles.c:1245
cs_halo_rotation_t
Definition: cs_halo.h:59
void( cs_sles_error_handler_t)(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)
Function pointer for handling of non-convergence when solving a linear system.
Definition: cs_sles.h:241
void( cs_sles_free_t)(void *context)
Function pointer for freeing of a linear system's context data.
Definition: cs_sles.h:158
void cs_sles_setup(cs_sles_t *sles, const cs_matrix_t *a)
Setup sparse linear equation solver.
Definition: cs_sles.c:1287
void( cs_sles_setup_t)(void *context, const char *name, const cs_matrix_t *a, int verbosity)
Function pointer for pre-resolution setup of a linear system's context.
Definition: cs_sles.h:87
void cs_sles_initialize(void)
Initialize sparse linear equation solver API.
Definition: cs_sles.c:782
Definition: cs_sles.h:60
#define BEGIN_C_DECLS
Definition: cs_defs.h:429
int( cs_sles_verbosity_t)(int f_id, const char *name)
Function pointer for the default definition of a sparse linear equation solver's verbosity.
Definition: cs_sles.h:293
void cs_sles_post_error_output_def(const char *name, int mesh_id, cs_halo_rotation_t rotation_mode, const cs_matrix_t *a, const cs_real_t *rhs, cs_real_t *vx)
Output default post-processing data for failed system convergence.
Definition: cs_sles.c:1569
cs_sles_t * cs_sles_find(int f_id, const char *name)
Return pointer to linear system object, based on matching field id or system name.
Definition: cs_sles.c:975
void cs_sles_set_default_verbosity(cs_sles_verbosity_t *verbosity_func)
Set default verbosity definition function.
Definition: cs_sles.c:1550
void cs_sles_push(int f_id, const char *name)
Temporarily replace field id with name for matching calls to cs_sles_setup, cs_sles_solve, cs_sles_free, and other operations involving access through a field id.
Definition: cs_sles.c:1076
cs_sles_convergence_state_t( cs_sles_solve_t)(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)
Function pointer for resolution of a linear system.
Definition: cs_sles.h:130
void cs_sles_free(cs_sles_t *sles)
Free sparse linear equation solver setup.
Definition: cs_sles.c:1425
struct _cs_matrix_t cs_matrix_t
Definition: cs_matrix.h:89
void * cs_sles_get_context(cs_sles_t *sles)
Return pointer to solver context structure pointer.
Definition: cs_sles.c:1265
void cs_sles_set_error_handler(cs_sles_t *sles, cs_sles_error_handler_t *error_handler_func)
Associate a convergence error handler to a given sparse linear equation solver.
Definition: cs_sles.c:1513
cs_sles_convergence_state_t
Convergence status indicator.
Definition: cs_sles.h:55
struct _cs_sles_t cs_sles_t
Definition: cs_sles.h:67
Definition: cs_sles.h:59
void cs_sles_finalize(void)
Finalize sparse linear equation solver API.
Definition: cs_sles.c:800
void cs_sles_post_error_output_var(const char *name, int mesh_id, int diag_block_size, cs_real_t *var)
Output post-processing variable for failed system convergence.
Definition: cs_sles.c:1660
double precision, save a
Definition: cs_fuel_incl.f90:146
int cs_sles_copy(cs_sles_t *dest, const cs_sles_t *src)
Copy the definition of a sparse linear equation solver to another.
Definition: cs_sles.c:1456
void( cs_sles_define_t)(int f_id, const char *name, const cs_matrix_t *a)
Function pointer for the default definition of a sparse linear equation solver.
Definition: cs_sles.h:272
cs_log_t
Definition: cs_log.h:47
#define END_C_DECLS
Definition: cs_defs.h:430
double cs_real_t
Definition: cs_defs.h:296
void( cs_sles_destroy_t)(void **context)
Definition: cs_sles.h:208
Definition: cs_sles.h:61
const char * cs_sles_base_name(int f_id, const char *name)
Return base name associated to a field id, name couple.
Definition: cs_sles.c:1747
Definition: cs_sles.h:58
void cs_sles_set_default_define(cs_sles_define_t *define_func)
Set default sparse linear solver definition function.
Definition: cs_sles.c:1533
void cs_sles_set_verbosity(cs_sles_t *sles, int verbosity)
Set the verbosity for a given linear equation solver.
Definition: cs_sles.c:1217
cs_sles_t * cs_sles_define(int f_id, const char *name, void *context, const char *type_name, cs_sles_setup_t *setup_func, cs_sles_solve_t *solve_func, cs_sles_free_t *free_func, cs_sles_log_t *log_func, cs_sles_copy_t *copy_func, cs_sles_destroy_t *destroy_func)
Define sparse linear equation solver for a given field or equation name.
Definition: cs_sles.c:1164
void cs_sles_pop(int f_id)
Restore behavior temporarily modified by cs_sles_push.
Definition: cs_sles.c:1112
void( cs_sles_log_t)(const void *context, cs_log_t log_type)
Function pointer for logging of linear solver history and performance data.
Definition: cs_sles.h:173
Definition: cs_sles.h:57
cs_sles_convergence_state_t cs_sles_solve(cs_sles_t *sles, const cs_matrix_t *a, 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)
General sparse linear system resolution.
Definition: cs_sles.c:1339