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-2016 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 using
219  * cs_sles_set_error_handler(), in which case it will be called whenever
220  * convergence fails.
221  *
222  * parameters:
223  * sles <-> pointer to solver object
224  * state <-- convergence status
225  * a <-- matrix
226  * rotation_mode <-- Halo update option for rotational periodicity
227  * rhs <-- Right hand side
228  * vx <-- System solution
229  *
230  * returns:
231  * true if solve should be re-executed, false otherwise
232  *----------------------------------------------------------------------------*/
233 
234 typedef bool
237  const cs_matrix_t *a,
238  cs_halo_rotation_t rotation_mode,
239  const cs_real_t *rhs,
240  cs_real_t *vx);
241 
242 /*----------------------------------------------------------------------------
243  * Function pointer for the default definition of a sparse
244  * linear equation solver
245  *
246  * The function may be associated using cs_sles_set_default_define(), so
247  * that it may provide a definition that will be used when
248  * cs_sles_setup() or cs_sles_solve() is used for a system for which
249  * no matching call to cs_sles_define() has been done.
250  *
251  * The function should call cs_sles_define() with arguments f_id
252  * and name, and appropriately chosen function pointers.
253  *
254  * A pointer to the matrix of the system to be solved is also provided,
255  * so that the corresponding information may be used to better choose
256  * defaults.
257  *
258  * parameters:
259  * f_id <-- associated field id, or < 0
260  * name <-- associated name if f_id < 0, or NULL
261  * a <-- Matrix
262  *----------------------------------------------------------------------------*/
263 
264 typedef void
265 (cs_sles_define_t) (int f_id,
266  const char *name,
267  const cs_matrix_t *a);
268 
269 /*----------------------------------------------------------------------------
270  * Function pointer for the default definition of a sparse
271  * linear equation solver's verbosity
272  *
273  * The function may be associated using cs_sles_set_default_verbosity(), so
274  * that it may provide a definition that will be used when
275  * cs_sles_default_verbosity() is called.
276  *
277  * parameters:
278  * f_id <-- associated field id, or < 0
279  * name <-- associated name if f_id < 0, or NULL
280  *
281  * returns:
282  * default verbosity value
283  *----------------------------------------------------------------------------*/
284 
285 typedef int
287  const char *name);
288 
289 /*============================================================================
290  * Global variables
291  *============================================================================*/
292 
293 /*=============================================================================
294  * Public function prototypes for Fortran API
295  *============================================================================*/
296 
297 /*=============================================================================
298  * Public function prototypes
299  *============================================================================*/
300 
301 /*----------------------------------------------------------------------------
302  * \brief Initialize sparse linear equation solver API.
303  *----------------------------------------------------------------------------*/
304 
305 void
306 cs_sles_initialize(void);
307 
308 /*----------------------------------------------------------------------------
309  * \brief Finalize sparse linear equation solver API.
310  *----------------------------------------------------------------------------*/
311 
312 void
313 cs_sles_finalize(void);
314 
315 /*----------------------------------------------------------------------------*/
321 /*----------------------------------------------------------------------------*/
322 
323 void
324 cs_sles_log(cs_log_t log_type);
325 
326 /*----------------------------------------------------------------------------*/
338 /*----------------------------------------------------------------------------*/
339 
340 cs_sles_t *
341 cs_sles_find(int f_id,
342  const char *name);
343 
344 /*----------------------------------------------------------------------------*/
361 /*----------------------------------------------------------------------------*/
362 
363 cs_sles_t *
364 cs_sles_find_or_add(int f_id,
365  const char *name);
366 
367 /*----------------------------------------------------------------------------*/
383 /*----------------------------------------------------------------------------*/
384 
385 void
386 cs_sles_push(int f_id,
387  const char *name);
388 
389 /*----------------------------------------------------------------------------*/
397 /*----------------------------------------------------------------------------*/
398 
399 void
400 cs_sles_pop(int f_id);
401 
402 /*----------------------------------------------------------------------------*/
437 /*----------------------------------------------------------------------------*/
438 
439 cs_sles_t *
440 cs_sles_define(int f_id,
441  const char *name,
442  void *context,
443  const char *type_name,
444  cs_sles_setup_t *setup_func,
445  cs_sles_solve_t *solve_func,
446  cs_sles_free_t *free_func,
447  cs_sles_log_t *log_func,
448  cs_sles_copy_t *copy_func,
449  cs_sles_destroy_t *destroy_func);
450 
451 /*----------------------------------------------------------------------------*/
463 /*----------------------------------------------------------------------------*/
464 
465 void
467  int verbosity);
468 
469 /*----------------------------------------------------------------------------*/
479 /*----------------------------------------------------------------------------*/
480 
481 int
483 
484 /*----------------------------------------------------------------------------*/
503 /*----------------------------------------------------------------------------*/
504 
505 const char *
507 
508 /*----------------------------------------------------------------------------*/
520 /*----------------------------------------------------------------------------*/
521 
522 void *
524 
525 /*----------------------------------------------------------------------------*/
533 /*----------------------------------------------------------------------------*/
534 
535 int
536 cs_sles_get_f_id(const cs_sles_t *sles);
537 
538 /*----------------------------------------------------------------------------*/
549 /*----------------------------------------------------------------------------*/
550 
551 const char *
552 cs_sles_get_name(const cs_sles_t *sles);
553 
554 /*----------------------------------------------------------------------------*/
568 /*----------------------------------------------------------------------------*/
569 
570 void
572  const cs_matrix_t *a);
573 
574 /*----------------------------------------------------------------------------*/
605 /*----------------------------------------------------------------------------*/
606 
609  const cs_matrix_t *a,
610  cs_halo_rotation_t rotation_mode,
611  double precision,
612  double r_norm,
613  int *n_iter,
614  double *residue,
615  const cs_real_t *rhs,
616  cs_real_t *vx,
617  size_t aux_size,
618  void *aux_vectors);
619 
620 /*----------------------------------------------------------------------------*/
631 /*----------------------------------------------------------------------------*/
632 
633 void
634 cs_sles_free(cs_sles_t *sles);
635 
636 /*----------------------------------------------------------------------------*/
656 /*----------------------------------------------------------------------------*/
657 
658 int
659 cs_sles_copy(cs_sles_t *dest,
660  const cs_sles_t *src);
661 
662 /*----------------------------------------------------------------------------*/
677 /*----------------------------------------------------------------------------*/
678 
679 void
681  cs_sles_error_handler_t *error_handler_func);
682 
683 /*----------------------------------------------------------------------------*/
693 /*----------------------------------------------------------------------------*/
694 
697 
698 /*----------------------------------------------------------------------------*/
708 /*----------------------------------------------------------------------------*/
709 
710 void
712 
713 /*----------------------------------------------------------------------------*/
722 /*----------------------------------------------------------------------------*/
723 
724 void
726 
727 /*----------------------------------------------------------------------------*/
738 /*----------------------------------------------------------------------------*/
739 
740 void
741 cs_sles_post_error_output_def(const char *name,
742  int mesh_id,
743  cs_halo_rotation_t rotation_mode,
744  const cs_matrix_t *a,
745  const cs_real_t *rhs,
746  cs_real_t *vx);
747 
748 /*----------------------------------------------------------------------------*/
757 /*----------------------------------------------------------------------------*/
758 
759 void
760 cs_sles_post_error_output_var(const char *name,
761  int mesh_id,
762  int diag_block_size,
763  cs_real_t *var);
764 
765 /*----------------------------------------------------------------------------*/
777 /*----------------------------------------------------------------------------*/
778 
779 const char *
780 cs_sles_base_name(int f_id,
781  const char *name);
782 
783 /*----------------------------------------------------------------------------*/
792 /*----------------------------------------------------------------------------*/
793 
794 const char *
795 cs_sles_name(int f_id,
796  const char *name);
797 
798 /*----------------------------------------------------------------------------*/
799 
801 
802 #endif /* __CS_SLES_H__ */
cs_sles_define_t * cs_sles_get_default_define(void)
Return pointer to default sparse linear solver definition function.
Definition: cs_sles.c:1585
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:1033
void cs_sles_log(cs_log_t log_type)
Log sparse linear equation solver info.
Definition: cs_sles.c:829
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:1256
cs_halo_rotation_t
Definition: cs_halo.h:59
void( cs_sles_free_t)(void *context)
Function pointer for freeing of a linear system&#39;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:1333
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&#39;s context.
Definition: cs_sles.h:87
void cs_sles_initialize(void)
Initialize sparse linear equation solver API.
Definition: cs_sles.c:775
Definition: cs_sles.h:60
#define BEGIN_C_DECLS
Definition: cs_defs.h:448
int( cs_sles_verbosity_t)(int f_id, const char *name)
Function pointer for the default definition of a sparse linear equation solver&#39;s verbosity.
Definition: cs_sles.h:286
const char * cs_sles_get_name(const cs_sles_t *sles)
Return name associated with a given sparse linear equation solver.
Definition: cs_sles.c:1311
const char * cs_sles_name(int f_id, const char *name)
Return name associated to a field id, name couple.
Definition: cs_sles.c:1842
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:1639
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:968
void cs_sles_set_default_verbosity(cs_sles_verbosity_t *verbosity_func)
Set default verbosity definition function.
Definition: cs_sles.c:1620
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:1069
double cs_real_t
Floating-point value.
Definition: cs_defs.h:296
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:1477
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:1276
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:1565
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:793
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:1730
int cs_sles_get_f_id(const cs_sles_t *sles)
Return field id associated with a given sparse linear equation solver.
Definition: cs_sles.c:1292
double precision, save a
Definition: cs_fuel_incl.f90:146
int cs_sles_get_verbosity(cs_sles_t *sles)
Get the verbosity for a given linear equation solver.
Definition: cs_sles.c:1229
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:1508
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:265
cs_log_t
Definition: cs_log.h:47
#define END_C_DECLS
Definition: cs_defs.h:449
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:1817
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:1603
void cs_sles_set_verbosity(cs_sles_t *sles, int verbosity)
Set the verbosity for a given linear equation solver.
Definition: cs_sles.c:1210
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:1157
void cs_sles_pop(int f_id)
Restore behavior temporarily modified by cs_sles_push.
Definition: cs_sles.c:1105
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:1385
bool( cs_sles_error_handler_t)(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)
Function pointer for handling of non-convergence when solving a linear system.
Definition: cs_sles.h:235