8.3
general 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-2024 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_matrix.h"
36#include "cs_sles.h"
37#include "cs_sles_pc.h"
38
39/*----------------------------------------------------------------------------*/
40
42
43/*============================================================================
44 * Macro definitions
45 *============================================================================*/
46
47/*============================================================================
48 * Type definitions
49 *============================================================================*/
50
51/*----------------------------------------------------------------------------
52 * Solver types
53 *----------------------------------------------------------------------------*/
54
55typedef enum {
56
83
84/* Iterative linear solver context (opaque) */
85
86typedef struct _cs_sles_it_t cs_sles_it_t;
87
88/* Forward type declarations */
89
90typedef struct _cs_sles_it_convergence_t cs_sles_it_convergence_t;
91
92/*============================================================================
93 * Global variables
94 *============================================================================*/
95
96/* Names for iterative solver types */
97
98extern const char *cs_sles_it_type_name[];
99
100/*=============================================================================
101 * User function prototypes
102 *============================================================================*/
103
104/*----------------------------------------------------------------------------
105 * Solution of A.vx = Rhs using a user-defined iterative solver
106 *
107 * On entry, vx is considered initialized.
108 *
109 * parameters:
110 * c <-- pointer to solver context info
111 * a <-- matrix
112 * diag_block_size <-- diagonal block size (unused here)
113 * convergence <-- convergence information structure
114 * rhs <-- right hand side
115 * vx_ini <-- initial system solution
116 * (vx if nonzero, nullptr if zero)
117 * vx <-> system solution
118 * aux_size <-- number of elements in aux_vectors (in bytes)
119 * aux_vectors --- optional working area (allocation otherwise)
120 *
121 * returns:
122 * convergence state
123 *----------------------------------------------------------------------------*/
124
127 const cs_matrix_t *a,
128 cs_lnum_t diag_block_size,
129 cs_sles_it_convergence_t *convergence,
130 const cs_real_t *rhs,
131 cs_real_t *vx_ini,
132 cs_real_t *vx,
133 size_t aux_size,
134 void *aux_vectors);
135
136/*=============================================================================
137 * Public function prototypes
138 *============================================================================*/
139
140/*----------------------------------------------------------------------------
141 * Define and associate an iterative sparse linear system solver
142 * for a given field or equation name.
143 *
144 * If this system did not previously exist, it is added to the list of
145 * "known" systems. Otherwise, its definition is replaced by the one
146 * defined here.
147 *
148 * This is a utility function: if finer control is needed, see
149 * cs_sles_define() and cs_sles_it_create().
150 *
151 * Note that this function returns a pointer directly to the iterative solver
152 * management structure. This may be used to set further options,
153 * for example using cs_sles_set_plot_options(). If needed, cs_sles_find()
154 * may be used to obtain a pointer to the matching cs_sles_t container.
155 *
156 * parameters:
157 * f_id <-- associated field id, or < 0
158 * name <-- associated name if f_id < 0, or NULL
159 * solver_type <-- type of solver (PCG, Jacobi, ...)
160 * poly_degree <-- preconditioning polynomial degree
161 * (0: diagonal; -1: non-preconditioned)
162 * n_max_iter <-- maximum number of iterations
163 *
164 * returns:
165 * pointer to newly created iterative solver info object.
166 *----------------------------------------------------------------------------*/
167
169cs_sles_it_define(int f_id,
170 const char *name,
171 cs_sles_it_type_t solver_type,
172 int poly_degree,
173 int n_max_iter);
174
175/*----------------------------------------------------------------------------
176 * Create iterative sparse linear system solver info and context.
177 *
178 * parameters:
179 * solver_type <-- type of solver (PCG, Jacobi, ...)
180 * poly_degree <-- preconditioning polynomial degree
181 * (0: diagonal; -1: non-preconditioned)
182 * n_max_iter <-- maximum number of iterations
183 * update_stats <-- automatic solver statistics indicator
184 *
185 * returns:
186 * pointer to newly created solver info object.
187 *----------------------------------------------------------------------------*/
188
191 int poly_degree,
192 int n_max_iter,
193 bool update_stats);
194
195/*----------------------------------------------------------------------------
196 * Destroy iterative sparse linear system solver info and context.
197 *
198 * parameters:
199 * context <-> pointer to iterative sparse linear solver info
200 * (actual type: cs_sles_it_t **)
201 *----------------------------------------------------------------------------*/
202
203void
204cs_sles_it_destroy(void **context);
205
206/*----------------------------------------------------------------------------
207 * Create iterative sparse linear system solver info and context
208 * based on existing info and context.
209 *
210 * parameters:
211 * context <-- pointer to reference info and context
212 * (actual type: cs_sles_it_t *)
213 *
214 * returns:
215 * pointer to newly created solver info object
216 * (actual type: cs_sles_it_t *)
217 *----------------------------------------------------------------------------*/
218
219void *
220cs_sles_it_copy(const void *context);
221
222/*----------------------------------------------------------------------------
223 * Setup iterative sparse linear equation solver.
224 *
225 * parameters:
226 * context <-> pointer to iterative sparse linear solver info
227 * (actual type: cs_sles_it_t *)
228 * name <-- pointer to system name
229 * a <-- associated matrix
230 * verbosity <-- verbosity level
231 *----------------------------------------------------------------------------*/
232
233void
234cs_sles_it_setup(void *context,
235 const char *name,
236 const cs_matrix_t *a,
237 int verbosity);
238
239/*----------------------------------------------------------------------------
240 * Call iterative sparse linear equation solver.
241 *
242 * parameters:
243 * context <-> pointer to iterative sparse linear solver info
244 * (actual type: cs_sles_it_t *)
245 * name <-- pointer to system name
246 * a <-- matrix
247 * verbosity <-- verbosity level
248 * precision <-- solver precision
249 * r_norm <-- residual normalization
250 * n_iter --> number of iterations
251 * residual --> residual
252 * rhs <-- right hand side
253 * vx_ini <-- initial solution (vx if nonzero, nullptr if zero)
254 * vx <-> system solution
255 * aux_size <-- number of elements in aux_vectors (in bytes)
256 * aux_vectors --- optional working area (internal allocation if NULL)
257 *
258 * returns:
259 * convergence state
260 *----------------------------------------------------------------------------*/
261
263cs_sles_it_solve(void *context,
264 const char *name,
265 const cs_matrix_t *a,
266 int verbosity,
267 double precision,
268 double r_norm,
269 int *n_iter,
270 double *residual,
271 const cs_real_t *rhs,
272 cs_real_t *vx,
273 cs_real_t *vx_ini,
274 size_t aux_size,
275 void *aux_vectors);
276
277/*----------------------------------------------------------------------------
278 * Free iterative sparse linear equation solver setup context.
279 *
280 * This function frees resolution-related data, such as
281 * buffers and preconditioning but does not free the whole context,
282 * as info used for logging (especially performance data) is maintained.
283
284 * parameters:
285 * context <-> pointer to iterative sparse linear solver info
286 * (actual type: cs_sles_it_t *)
287 *----------------------------------------------------------------------------*/
288
289void
290cs_sles_it_free(void *context);
291
292/*----------------------------------------------------------------------------
293 * Log sparse linear equation solver info.
294 *
295 * parameters:
296 * context <-> pointer to iterative sparse linear solver info
297 * (actual type: cs_sles_it_t *)
298 * log_type <-- log type
299 *----------------------------------------------------------------------------*/
300
301void
302cs_sles_it_log(const void *context,
303 cs_log_t log_type);
304
305/*----------------------------------------------------------------------------
306 * Return iterative solver type.
307 *
308 * parameters:
309 * context <-- pointer to iterative solver info and context
310 *
311 * returns:
312 * selected solver type
313 *----------------------------------------------------------------------------*/
314
316cs_sles_it_get_type(const cs_sles_it_t *context);
317
318/*----------------------------------------------------------------------------
319 * Return the initial residual for the previous solve operation with a solver.
320 *
321 * This is useful for convergence tests when this solver is used as
322 * a preconditioning smoother.
323 *
324 * This operation is only valid between calls to cs_sles_it_setup()
325 * (or cs_sles_it_solve()) and cs_sles_it_free().
326 * It returns -1 otherwise.
327 *
328 * parameters:
329 * context <-- pointer to iterative solver info and context
330 *
331 * returns:
332 * initial residual from last call to \ref cs_sles_solve with this solver
333 *----------------------------------------------------------------------------*/
334
335double
337
338/*----------------------------------------------------------------------------
339 * Return a preconditioner context for an iterative sparse linear
340 * equation solver.
341 *
342 * This allows modifying parameters of a non default (Jacobi or polynomial)
343 * preconditioner.
344 *
345 * parameters:
346 * context <-- pointer to iterative solver info and context
347 *
348 * returns:
349 * pointer to preconditoner context
350 *----------------------------------------------------------------------------*/
351
354
355/*----------------------------------------------------------------------------
356 * \brief Assign a preconditioner to an iterative sparse linear equation
357 * solver, transfering its ownership to the solver context.
358 *
359 * This allows assigning a non default (Jacobi or polynomial) preconditioner.
360 *
361 * The input pointer is set to NULL to make it clear the caller does not
362 * own the preconditioner anymore, though the context can be accessed using
363 * \ref cs_sles_it_get_pc.
364 *
365 * \param[in, out] context pointer to iterative solver info and context
366 * \param[in, out] pc pointer to preconditoner context
367 *----------------------------------------------------------------------------*/
368
369void
371 cs_sles_pc_t **pc);
372
373/*----------------------------------------------------------------------------
374 * Copy options from one iterative sparse linear system solver info
375 * and context to another.
376 *
377 * Optional plotting contexts are shared between the source and destination
378 * contexts.
379 *
380 * Preconditioner settings are to be handled separately.
381 *
382 * parameters:
383 * src <-- pointer to source info and context
384 * dest <-> pointer to destination info and context
385 *----------------------------------------------------------------------------*/
386
387void
389 cs_sles_it_t *dest);
390
391/*----------------------------------------------------------------------------
392 * Associate a similar info and context object with which some setup
393 * data may be shared.
394 *
395 * This is especially useful for sharing preconditioning data between
396 * similar solver contexts (for example ascending and descending multigrid
397 * smoothers based on the same matrix).
398 *
399 * For preconditioning data to be effectively shared, cs_sles_it_setup()
400 * (or cs_sles_it_solve()) must be called on "shareable" before being
401 * called on "context" (without cs_sles_it_free() being called in between,
402 * of course).
403 *
404 * It is the caller's responsibility to ensure the context is not used
405 * for a cs_sles_it_setup() or cs_sles_it_solve() operation after the
406 * shareable object has been destroyed (normally by cs_sles_it_destroy()).
407 *
408 * parameters:
409 * context <-> pointer to iterative sparse linear system solver info
410 * shareable <-- pointer to iterative solver info and context
411 * whose context may be shared
412 *----------------------------------------------------------------------------*/
413
414void
416 const cs_sles_it_t *shareable);
417
418#if defined(HAVE_MPI)
419
420/*----------------------------------------------------------------------------*/
432/*----------------------------------------------------------------------------*/
433
434void
435cs_sles_it_set_mpi_reduce_comm(cs_sles_it_t *context,
436 MPI_Comm comm,
437 MPI_Comm caller_comm);
438
439#endif /* defined(HAVE_MPI) */
440
441/*----------------------------------------------------------------------------
442 * Assign ordering to iterative solver.
443 *
444 * The solver context takes ownership of the order array (i.e. it will
445 * handle its later deallocation).
446 *
447 * This is useful only for Block Gauss-Seidel.
448 *
449 * parameters:
450 * context <-> pointer to iterative solver info and context
451 * order <-> pointer to ordering array
452 *----------------------------------------------------------------------------*/
453
454void
456 cs_lnum_t **order);
457
458/*----------------------------------------------------------------------------*/
465/*----------------------------------------------------------------------------*/
466
467double
469
470/*----------------------------------------------------------------------------*/
477/*----------------------------------------------------------------------------*/
478
479void
480cs_sles_it_set_breakdown_threshold(double threshold);
481
482/*----------------------------------------------------------------------------*/
499/*----------------------------------------------------------------------------*/
500
501void
504 int n_iter_max);
505
506/*----------------------------------------------------------------------------*/
514/*----------------------------------------------------------------------------*/
515
516void
518 int interval);
519
520/*----------------------------------------------------------------------------*/
527/*----------------------------------------------------------------------------*/
528
529void
531 int n_max_iter);
532
533/*----------------------------------------------------------------------------
534 * Query mean number of rows under which Conjugate Gradient algorithm
535 * uses the single-reduction variant.
536 *
537 * The single-reduction variant requires only one parallel sum per
538 * iteration (instead of 2), at the cost of additional vector operations,
539 * so it tends to be more expensive when the number of matrix rows per
540 * MPI rank is high, then becomes cheaper when the MPI latency cost becomes
541 * more significant.
542 *
543 * This option is ignored for non-parallel runs, so 0 is returned.
544 *
545 * return:
546 * mean number of rows per active rank under which the
547 * single-reduction variant will be used
548 *----------------------------------------------------------------------------*/
549
552
553/*----------------------------------------------------------------------------
554 * Set mean number of rows under which Conjugate Gradient algorithm
555 * should use the single-reduction variant.
556 *
557 * The single-reduction variant requires only one parallel sum per
558 * iteration (instead of 2), at the cost of additional vector operations,
559 * so it tends to be more expensive when the number of matrix rows per
560 * MPI rank is high, then becomes cheaper when the MPI latency cost becomes
561 * more significant.
562 *
563 * This option is ignored for non-parallel runs.
564 *
565 * parameters:
566 * threshold <-- mean number of rows per active rank under which the
567 * single-reduction variant will be used
568 *----------------------------------------------------------------------------*/
569
570void
572
573/*----------------------------------------------------------------------------
574 * Log the current global settings relative to parallelism.
575 *----------------------------------------------------------------------------*/
576
577void
579
580/*----------------------------------------------------------------------------
581 * Error handler for iterative sparse linear equation solver.
582 *
583 * In case of divergence or breakdown, this error handler outputs
584 * postprocessing data to assist debugging, then aborts the run.
585 * It does nothing in case the maximum iteration count is reached.
586 *
587 * parameters:
588 * sles <-> pointer to solver object
589 * state <-- convergence state
590 * a <-- matrix
591 * rhs <-- right hand side
592 * vx <-> system solution
593 *
594 * returns:
595 * false (do not attempt new solve)
596 *----------------------------------------------------------------------------*/
597
598bool
601 const cs_matrix_t *a,
602 const cs_real_t *rhs,
603 cs_real_t *vx);
604
605/*----------------------------------------------------------------------------
606 * Set plotting options for an iterative sparse linear equation solver.
607 *
608 * parameters:
609 * context <-> pointer to iterative solver info and context
610 * base_name <-- base plot name to activate, NULL otherwise
611 * use_iteration <-- if true, use iteration as time stamp
612 * otherwise, use wall clock time
613 *----------------------------------------------------------------------------*/
614
615void
617 const char *base_name,
618 bool use_iteration);
619
620/*----------------------------------------------------------------------------
621 * Convergence test.
622 *
623 * parameters:
624 * c <-- pointer to solver context info
625 * n_iter <-- Number of iterations done
626 * residual <-- Non normalized residual
627 * convergence <-> Convergence information structure
628 *
629 * returns:
630 * convergence status.
631 *----------------------------------------------------------------------------*/
632
635 unsigned n_iter,
636 double residual,
637 cs_sles_it_convergence_t *convergence);
638
639/*----------------------------------------------------------------------------*/
640
642
643#endif /* __CS_SLES_IT_H__ */
#define BEGIN_C_DECLS
Definition: cs_defs.h:542
double cs_real_t
Floating-point value.
Definition: cs_defs.h:342
#define END_C_DECLS
Definition: cs_defs.h:543
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:335
cs_log_t
Definition: cs_log.h:48
struct _cs_matrix_t cs_matrix_t
Definition: cs_matrix.h:110
cs_sles_convergence_state_t
Definition: cs_sles.h:56
struct _cs_sles_t cs_sles_t
Definition: cs_sles.h:68
double cs_sles_it_get_last_initial_residual(const cs_sles_it_t *context)
Return the initial residual for the previous solve operation with a solver.
Definition: cs_sles_it.cpp:5074
cs_sles_convergence_state_t cs_sles_it_solve(void *context, const char *name, const cs_matrix_t *a, int verbosity, double precision, double r_norm, int *n_iter, double *residual, const cs_real_t *rhs, cs_real_t *vx, cs_real_t *vx_ini, size_t aux_size, void *aux_vectors)
Call iterative sparse linear equation solver.
Definition: cs_sles_it.cpp:4767
double cs_sles_it_get_breakdown_threshold(void)
Retrieve the threshold value under which a breakdown happens in solvers like BiCGStab or BiCGStab2.
Definition: cs_sles_it.cpp:5312
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.cpp:5127
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.cpp:5098
void cs_sles_it_transfer_parameters(const cs_sles_it_t *src, cs_sles_it_t *dest)
Copy options from one iterative sparse linear system solver info and context to another.
Definition: cs_sles_it.cpp:5159
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.cpp:4276
bool cs_sles_it_error_post_and_abort(cs_sles_t *sles, cs_sles_convergence_state_t state, const cs_matrix_t *a, const cs_real_t *rhs, cs_real_t *vx)
Error handler for iterative sparse linear equation solver.
Definition: cs_sles_it.cpp:5491
void cs_sles_it_set_n_max_iter(cs_sles_it_t *context, int n_max_iter)
Define the max. number of iterations before stopping the algorithm.
Definition: cs_sles_it.cpp:5390
void cs_sles_it_set_fallback_threshold(cs_sles_it_t *context, cs_sles_convergence_state_t threshold, int n_iter_max)
Define convergence level under which the fallback to another solver may be used if applicable.
Definition: cs_sles_it.cpp:5352
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.cpp:4571
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.cpp:5205
void cs_sles_it_set_restart_interval(cs_sles_it_t *context, int interval)
Define the number of iterations to be done before restarting the solver. Useful only for GCR or GMRES...
Definition: cs_sles_it.cpp:5371
void cs_sles_it_log(const void *context, cs_log_t log_type)
Log sparse linear equation solver info.
Definition: cs_sles_it.cpp:4455
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.cpp:4223
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.cpp:4413
cs_sles_convergence_state_t cs_sles_it_convergence_test(cs_sles_it_t *c, unsigned n_iter, double residual, cs_sles_it_convergence_t *convergence)
Definition: cs_sles_it.cpp:5582
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.cpp:5538
struct _cs_sles_it_t cs_sles_it_t
Definition: cs_sles_it.h:86
void cs_sles_it_log_parallel_options(void)
Log the current global settings relative to parallelism.
Definition: cs_sles_it.cpp:5460
void cs_sles_it_free(void *context)
Free iterative sparse linear equation solver setup context.
Definition: cs_sles_it.cpp:5007
void cs_sles_it_assign_order(cs_sles_it_t *context, cs_lnum_t **order)
Assign ordering to iterative solver.
Definition: cs_sles_it.cpp:5280
const char * cs_sles_it_type_name[]
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.cpp:5418
struct _cs_sles_it_convergence_t cs_sles_it_convergence_t
Definition: cs_sles_it.h:90
cs_sles_it_type_t
Definition: cs_sles_it.h:55
@ CS_SLES_P_SYM_GAUSS_SEIDEL
Definition: cs_sles_it.h:69
@ CS_SLES_JACOBI
Definition: cs_sles_it.h:61
@ CS_SLES_PCR3
Definition: cs_sles_it.h:70
@ CS_SLES_GMRES
Definition: cs_sles_it.h:66
@ CS_SLES_GCR
Definition: cs_sles_it.h:65
@ CS_SLES_P_GAUSS_SEIDEL
Definition: cs_sles_it.h:68
@ CS_SLES_BICGSTAB
Definition: cs_sles_it.h:62
@ CS_SLES_N_IT_TYPES
Definition: cs_sles_it.h:73
@ CS_SLES_BICGSTAB2
Definition: cs_sles_it.h:64
@ CS_SLES_USER_DEFINED
Definition: cs_sles_it.h:71
@ CS_SLES_PCG
Definition: cs_sles_it.h:57
@ CS_SLES_TS_B_GAUSS_SEIDEL
Definition: cs_sles_it.h:77
@ CS_SLES_IPCG
Definition: cs_sles_it.h:60
@ CS_SLES_N_SMOOTHER_TYPES
Definition: cs_sles_it.h:79
@ CS_SLES_FCG
Definition: cs_sles_it.h:58
@ CS_SLES_TS_F_GAUSS_SEIDEL
Definition: cs_sles_it.h:76
void cs_sles_it_set_breakdown_threshold(double threshold)
Define the threshold value under which a breakdown happens in solvers like BiCGStab or BiCGStab2.
Definition: cs_sles_it.cpp:5327
cs_sles_convergence_state_t cs_user_sles_it_solver(cs_sles_it_t *c, const cs_matrix_t *a, cs_lnum_t diag_block_size, cs_sles_it_convergence_t *convergence, const cs_real_t *rhs, cs_real_t *vx_ini, cs_real_t *vx, size_t aux_size, void *aux_vectors)
void cs_sles_it_destroy(void **context)
Destroy iterative sparse linear system solver info and context.
Definition: cs_sles_it.cpp:4375
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.cpp:5446
cs_sles_it_type_t cs_sles_it_get_type(const cs_sles_it_t *context)
Return iterative solver type.
Definition: cs_sles_it.cpp:5046
struct _cs_sles_pc_t cs_sles_pc_t
Definition: cs_sles_pc.h:66