8.0
general documentation
Loading...
Searching...
No Matches
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-2023 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#include "cs_matrix.h"
39
40/*----------------------------------------------------------------------------*/
41
43
44/*============================================================================
45 * Macro definitions
46 *============================================================================*/
47
48/*============================================================================
49 * Type definitions
50 *============================================================================*/
51
52/*----------------------------------------------------------------------------
53 * Convergence status
54 *----------------------------------------------------------------------------*/
55
65
66/* General linear solver context (opaque) */
67
68typedef struct _cs_sles_t cs_sles_t;
69
70/*----------------------------------------------------------------------------
71 * Function pointer for pre-resolution setup of a linear system solvers's
72 * context.
73 *
74 * This setup may include building a multigrid hierarchy, or a preconditioner.
75 *
76 * Use of this type of function is optional: the context is expected to
77 * maintain state, so that if a cs_sles_solve_t function is called before a
78 * cs_sles_setup_t function, the latter will be called automatically.
79 *
80 * parameters:
81 * context <-> pointer to solver context
82 * name <-- pointer to name of linear system
83 * a <-- matrix
84 * verbosity <-- associated verbosity
85 *----------------------------------------------------------------------------*/
86
87typedef void
88(cs_sles_setup_t) (void *context,
89 const char *name,
90 const cs_matrix_t *a,
91 int verbosity);
92
93/*----------------------------------------------------------------------------
94 * Function pointer for resolution of a linear system.
95 *
96 * If the associated cs_sles_setup_t function has not been called before
97 * this function, it will be called automatically.
98 *
99 * The solution context setup by this call (or that of the matching setup
100 * function) will be maintained until the matching cs_sles_free_t function
101 * is called.
102 *
103 * The matrix is not expected to change between successive calls, although
104 * the right hand side may. If the matrix changes, the associated
105 * cs_sles_setup_t or cs_sles_free_t function must be called between
106 * solves.
107 *
108 * The system is considered to have converged when
109 * residue/r_norm <= precision, residue being the L2 norm of a.vx-rhs.
110 *
111 * parameters:
112 * context <-> pointer to solver context
113 * name <-- pointer to name of linear system
114 * a <-- matrix
115 * verbosity <-- associated verbosity
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 double precision,
135 double r_norm,
136 int *n_iter,
137 double *residue,
138 const cs_real_t *rhs,
139 cs_real_t *vx,
140 size_t aux_size,
141 void *aux_vectors);
142
143/*----------------------------------------------------------------------------
144 * Function pointer for freeing of a linear system's context data.
145 *
146 * Note that this function should free resolution-related data, such as
147 * multigrid hierarchy, preconditioning, and any other temporary arrays or
148 * objects required for resolution, but should not free the whole context,
149 * as info used for logging (especially performance data) should be
150 * maintained.
151 *
152 * parameters:
153 * context <-> pointer to solver context
154 *----------------------------------------------------------------------------*/
155
156typedef void
157(cs_sles_free_t) (void *context);
158
159/*----------------------------------------------------------------------------
160 * Function pointer for logging of linear solver setup,
161 * history and performance data.
162 *
163 * This function will be called for each solver when cs_sles_finalize()
164 * is called.
165 *
166 * parameters:
167 * context <-- pointer to solver context
168 * log_type <-- log type
169 *----------------------------------------------------------------------------*/
170
171typedef void
172(cs_sles_log_t) (const void *context,
173 cs_log_t log_type);
174
175/*----------------------------------------------------------------------------
176 * Function pointer for creation of a solver context based on the copy
177 * of another.
178 *
179 * The new context copies the settings of the copied context, but not
180 * its setup data and logged info, such as performance data.
181 *
182 * This type of function is optional, but enables associating different
183 * solvers to related systems (to differentiate logging) while using
184 * the same settings by default.
185 *
186 * parameters:
187 * context <-- source context
188 *
189 * returns:
190 * pointer to newly created context
191 *----------------------------------------------------------------------------*/
192
193typedef void *
194(cs_sles_copy_t) (const void *context);
195
196/*----------------------------------------------------------------------------
197 * Function pointer for destruction of a linear system solver context.
198 *
199 * This function should free all context data, and will be called for each
200 * system when cs_sles_finalize() is called.
201 *
202 * parameters:
203 * context <-> pointer to solver context
204 *----------------------------------------------------------------------------*/
205
206typedef void
207(cs_sles_destroy_t) (void **context);
208
209/*----------------------------------------------------------------------------
210 * Function pointer for handling of non-convegence when solving
211 * a linear system.
212 *
213 * Such a function is optional, and may be used for a variety of purposes,
214 * such as logging, postprocessing, re-trying with different parameters,
215 * aborting the run, or any combination thereof.
216 *
217 * An error handler may be associated with a given solver using
218 * cs_sles_set_error_handler(), in which case it will be called whenever
219 * convergence fails.
220 *
221 * parameters:
222 * sles <-> pointer to solver object
223 * state <-- convergence status
224 * a <-- matrix
225 * rhs <-- Right hand side
226 * vx <-- System solution
227 *
228 * returns:
229 * true if solve should be re-executed, false otherwise
230 *----------------------------------------------------------------------------*/
231
232typedef bool
235 const cs_matrix_t *a,
236 const cs_real_t *rhs,
237 cs_real_t *vx);
238
239/*----------------------------------------------------------------------------
240 * Function pointer for the default definition of a sparse
241 * linear equation solver
242 *
243 * The function may be associated using cs_sles_set_default_define(), so
244 * that it may provide a definition that will be used when
245 * cs_sles_setup() or cs_sles_solve() is used for a system for which
246 * no matching call to cs_sles_define() has been done.
247 *
248 * The function should call cs_sles_define() with arguments f_id
249 * and name, and appropriately chosen function pointers.
250 *
251 * A pointer to the matrix of the system to be solved is also provided,
252 * so that the corresponding information may be used to better choose
253 * defaults.
254 *
255 * parameters:
256 * f_id <-- associated field id, or < 0
257 * name <-- associated name if f_id < 0, or NULL
258 * a <-- Matrix
259 *----------------------------------------------------------------------------*/
260
261typedef void
262(cs_sles_define_t) (int f_id,
263 const char *name,
264 const cs_matrix_t *a);
265
266/*----------------------------------------------------------------------------
267 * Function pointer for the default definition of a sparse
268 * linear equation solver's verbosity
269 *
270 * The function may be associated using cs_sles_set_default_verbosity(), so
271 * that it may provide a definition that will be used when
272 * cs_sles_default_verbosity() is called.
273 *
274 * parameters:
275 * f_id <-- associated field id, or < 0
276 * name <-- associated name if f_id < 0, or NULL
277 *
278 * returns:
279 * default verbosity value
280 *----------------------------------------------------------------------------*/
281
282typedef int
283(cs_sles_verbosity_t) (int f_id,
284 const char *name);
285
286/*============================================================================
287 * Global variables
288 *============================================================================*/
289
290/*=============================================================================
291 * Public function prototypes for Fortran API
292 *============================================================================*/
293
294/*=============================================================================
295 * Public function prototypes
296 *============================================================================*/
297
298/*----------------------------------------------------------------------------*/
302/*----------------------------------------------------------------------------*/
303
304void
305cs_sles_set_epzero(double new_value);
306
307/*----------------------------------------------------------------------------*/
314/*----------------------------------------------------------------------------*/
315
316double
318
319/*----------------------------------------------------------------------------
320 * \brief Initialize sparse linear equation solver API.
321 *----------------------------------------------------------------------------*/
322
323void
325
326/*----------------------------------------------------------------------------
327 * \brief Finalize sparse linear equation solver API.
328 *----------------------------------------------------------------------------*/
329
330void
331cs_sles_finalize(void);
332
333/*----------------------------------------------------------------------------*/
339/*----------------------------------------------------------------------------*/
340
341void
342cs_sles_log(cs_log_t log_type);
343
344/*----------------------------------------------------------------------------*/
356/*----------------------------------------------------------------------------*/
357
358cs_sles_t *
359cs_sles_find(int f_id,
360 const char *name);
361
362/*----------------------------------------------------------------------------*/
379/*----------------------------------------------------------------------------*/
380
381cs_sles_t *
382cs_sles_find_or_add(int f_id,
383 const char *name);
384
385/*----------------------------------------------------------------------------*/
402/*----------------------------------------------------------------------------*/
403
404void
405cs_sles_push(int f_id,
406 const char *name);
407
408/*----------------------------------------------------------------------------*/
416/*----------------------------------------------------------------------------*/
417
418void
419cs_sles_pop(int f_id);
420
421/*----------------------------------------------------------------------------*/
456/*----------------------------------------------------------------------------*/
457
458cs_sles_t *
459cs_sles_define(int f_id,
460 const char *name,
461 void *context,
462 const char *type_name,
463 cs_sles_setup_t *setup_func,
464 cs_sles_solve_t *solve_func,
465 cs_sles_free_t *free_func,
466 cs_sles_log_t *log_func,
467 cs_sles_copy_t *copy_func,
468 cs_sles_destroy_t *destroy_func);
469
470/*----------------------------------------------------------------------------*/
482/*----------------------------------------------------------------------------*/
483
484void
486 int verbosity);
487
488/*----------------------------------------------------------------------------*/
498/*----------------------------------------------------------------------------*/
499
500int
502
503/*----------------------------------------------------------------------------*/
514/*----------------------------------------------------------------------------*/
515
516void
518 int writer_id);
519
520/*----------------------------------------------------------------------------*/
529/*----------------------------------------------------------------------------*/
530
531int
533
534/*----------------------------------------------------------------------------*/
553/*----------------------------------------------------------------------------*/
554
555const char *
557
558/*----------------------------------------------------------------------------*/
570/*----------------------------------------------------------------------------*/
571
572void *
574
575/*----------------------------------------------------------------------------*/
583/*----------------------------------------------------------------------------*/
584
585int
586cs_sles_get_f_id(const cs_sles_t *sles);
587
588/*----------------------------------------------------------------------------*/
599/*----------------------------------------------------------------------------*/
600
601const char *
602cs_sles_get_name(const cs_sles_t *sles);
603
604/*----------------------------------------------------------------------------*/
615/*----------------------------------------------------------------------------*/
616
617bool
619
620/*----------------------------------------------------------------------------*/
630/*----------------------------------------------------------------------------*/
631
632void
634 bool allow_no_op);
635
636/*----------------------------------------------------------------------------*/
650/*----------------------------------------------------------------------------*/
651
652void
654 const cs_matrix_t *a);
655
656/*----------------------------------------------------------------------------*/
686/*----------------------------------------------------------------------------*/
687
690 const cs_matrix_t *a,
691 double precision,
692 double r_norm,
693 int *n_iter,
694 double *residue,
695 const cs_real_t *rhs,
696 cs_real_t *vx,
697 size_t aux_size,
698 void *aux_vectors);
699
700/*----------------------------------------------------------------------------*/
711/*----------------------------------------------------------------------------*/
712
713void
715
716/*----------------------------------------------------------------------------*/
736/*----------------------------------------------------------------------------*/
737
738int
740 const cs_sles_t *src);
741
742/*----------------------------------------------------------------------------*/
757/*----------------------------------------------------------------------------*/
758
759void
761 cs_sles_error_handler_t *error_handler_func);
762
763/*----------------------------------------------------------------------------*/
773/*----------------------------------------------------------------------------*/
774
777
778/*----------------------------------------------------------------------------*/
788/*----------------------------------------------------------------------------*/
789
790void
792
793/*----------------------------------------------------------------------------*/
802/*----------------------------------------------------------------------------*/
803
804void
806
807/*----------------------------------------------------------------------------*/
817/*----------------------------------------------------------------------------*/
818
819void
820cs_sles_post_error_output_def(const char *name,
821 int mesh_id,
822 const cs_matrix_t *a,
823 const cs_real_t *rhs,
824 cs_real_t *vx);
825
826/*----------------------------------------------------------------------------*/
838/*----------------------------------------------------------------------------*/
839
840void
841cs_sles_post_output_var(const char *name,
842 int mesh_id,
843 int location_id,
844 int writer_id,
845 int diag_block_size,
846 cs_real_t var[]);
847
848/*----------------------------------------------------------------------------*/
860/*----------------------------------------------------------------------------*/
861
862const char *
863cs_sles_base_name(int f_id,
864 const char *name);
865
866/*----------------------------------------------------------------------------*/
875/*----------------------------------------------------------------------------*/
876
877const char *
878cs_sles_name(int f_id,
879 const char *name);
880
881/*----------------------------------------------------------------------------*/
882
884
885#endif /* __CS_SLES_H__ */
#define BEGIN_C_DECLS
Definition cs_defs.h:509
double cs_real_t
Floating-point value.
Definition cs_defs.h:319
#define bool
Definition cs_defs.h:196
#define END_C_DECLS
Definition cs_defs.h:510
cs_log_t
Definition cs_log.h:48
struct _cs_matrix_t cs_matrix_t
Definition cs_matrix.h:110
void * cs_sles_get_context(cs_sles_t *sles)
Return pointer to solver context structure pointer.
Definition cs_sles.c:1535
void cs_sles_set_default_define(cs_sles_define_t *define_func)
Set default sparse linear solver definition function.
Definition cs_sles.c:1958
cs_sles_convergence_state_t
Definition cs_sles.h:56
@ CS_SLES_ITERATING
Definition cs_sles.h:61
@ CS_SLES_MAX_ITERATION
Definition cs_sles.h:60
@ CS_SLES_DIVERGED
Definition cs_sles.h:58
@ CS_SLES_CONVERGED
Definition cs_sles.h:62
@ CS_SLES_BREAKDOWN
Definition cs_sles.h:59
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:1863
void cs_sles_setup(cs_sles_t *sles, const cs_matrix_t *a)
Setup sparse linear equation solver.
Definition cs_sles.c:1630
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:262
void cs_sles_log(cs_log_t log_type)
Log sparse linear equation solver info.
Definition cs_sles.c:997
double cs_sles_get_epzero(void)
Get the current threshold value used in the detection of immediate exit.
Definition cs_sles.c:925
void cs_sles_post_output_var(const char *name, int mesh_id, int location_id, int writer_id, int diag_block_size, cs_real_t var[])
Output post-processing variable related to system convergence.
Definition cs_sles.c:2103
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:1570
int cs_sles_get_verbosity(cs_sles_t *sles)
Get the verbosity for a given linear equation solver.
Definition cs_sles.c:1436
void cs_sles_pop(int f_id)
Restore behavior temporarily modified by cs_sles_push.
Definition cs_sles.c:1311
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:1363
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,...
Definition cs_sles.c:1275
void cs_sles_free(cs_sles_t *sles)
Free sparse linear equation solver setup.
Definition cs_sles.c:1832
const char * cs_sles_name(int f_id, const char *name)
Return name associated to a field id, name couple.
Definition cs_sles.c:2244
int cs_sles_get_post_output(cs_sles_t *sles)
Return the id of the associated writer if postprocessing output is active for a given linear equation...
Definition cs_sles.c:1483
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:194
cs_sles_define_t * cs_sles_get_default_define(void)
Return pointer to default sparse linear solver definition function.
Definition cs_sles.c:1940
void cs_sles_set_verbosity(cs_sles_t *sles, int verbosity)
Set the verbosity for a given linear equation solver.
Definition cs_sles.c:1417
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:2219
bool cs_sles_get_allow_no_op(const cs_sles_t *sles)
Query if immediate_return ("no-op") is allowed when initial guess is zero (solve by increments) and t...
Definition cs_sles.c:1589
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:1238
void cs_sles_set_epzero(double new_value)
Set the threshold value used in the detection of immediate exit.
Definition cs_sles.c:910
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:1173
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:1551
void cs_sles_finalize(void)
Finalize sparse linear equation solver API.
Definition cs_sles.c:957
void cs_sles_post_error_output_def(const char *name, int mesh_id, 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:1993
const char * cs_sles_get_type(cs_sles_t *sles)
Return type name of solver context.
Definition cs_sles.c:1515
bool cs_sles_error_handler_t(cs_sles_t *sles, cs_sles_convergence_state_t state, const cs_matrix_t *a, 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:233
cs_sles_convergence_state_t cs_sles_solve(cs_sles_t *sles, const cs_matrix_t *a, 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:1698
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:283
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:172
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:88
void cs_sles_set_default_verbosity(cs_sles_verbosity_t *verbosity_func)
Set default verbosity definition function.
Definition cs_sles.c:1975
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:1920
void cs_sles_destroy_t(void **context)
Definition cs_sles.h:207
struct _cs_sles_t cs_sles_t
Definition cs_sles.h:68
cs_sles_convergence_state_t cs_sles_solve_t(void *context, const char *name, const cs_matrix_t *a, int verbosity, 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_initialize(void)
Initialize sparse linear equation solver API.
Definition cs_sles.c:937
void cs_sles_set_post_output(cs_sles_t *sles, int writer_id)
Activate postprocessing output for a given linear equation solver.
Definition cs_sles.c:1455
void cs_sles_free_t(void *context)
Function pointer for freeing of a linear system's context data.
Definition cs_sles.h:157
void cs_sles_set_allow_no_op(cs_sles_t *sles, bool allow_no_op)
Indicate if immediate_return ("no-op") is allowed when initial guess is zero (solve by increments) an...
Definition cs_sles.c:1607