8.3
general 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-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_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
56typedef enum {
57
63
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 * residual/r_norm <= precision, residual 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 <-- residual normalization
118 * n_iter --> number of "equivalent" iterations
119 * residual --> residual
120 * rhs <-- right hand side
121 * vx_ini <-- initial solution (vx if nonzero, nullptr if zero)
122 * vx --> system solution
123 * aux_size <-- number of elements in aux_vectors
124 * aux_vectors <-- optional working area (internal allocation if NULL)
125 *
126 * returns:
127 * convergence status
128 *----------------------------------------------------------------------------*/
129
131(cs_sles_solve_t) (void *context,
132 const char *name,
133 const cs_matrix_t *a,
134 int verbosity,
135 double precision,
136 double r_norm,
137 int *n_iter,
138 double *residual,
139 const cs_real_t *rhs,
140 cs_real_t *vx_ini,
141 cs_real_t *vx,
142 size_t aux_size,
143 void *aux_vectors);
144
145/*----------------------------------------------------------------------------
146 * Function pointer for freeing of a linear system's context data.
147 *
148 * Note that this function should free resolution-related data, such as
149 * multigrid hierarchy, preconditioning, and any other temporary arrays or
150 * objects required for resolution, but should not free the whole context,
151 * as info used for logging (especially performance data) should be
152 * maintained.
153 *
154 * parameters:
155 * context <-> pointer to solver context
156 *----------------------------------------------------------------------------*/
157
158typedef void
159(cs_sles_free_t) (void *context);
160
161/*----------------------------------------------------------------------------
162 * Function pointer for logging of linear solver setup,
163 * history and performance data.
164 *
165 * This function will be called for each solver when cs_sles_finalize()
166 * is called.
167 *
168 * parameters:
169 * context <-- pointer to solver context
170 * log_type <-- log type
171 *----------------------------------------------------------------------------*/
172
173typedef void
174(cs_sles_log_t) (const void *context,
175 cs_log_t log_type);
176
177/*----------------------------------------------------------------------------
178 * Function pointer for creation of a solver context based on the copy
179 * of another.
180 *
181 * The new context copies the settings of the copied context, but not
182 * its setup data and logged info, such as performance data.
183 *
184 * This type of function is optional, but enables associating different
185 * solvers to related systems (to differentiate logging) while using
186 * the same settings by default.
187 *
188 * parameters:
189 * context <-- source context
190 *
191 * returns:
192 * pointer to newly created context
193 *----------------------------------------------------------------------------*/
194
195typedef void *
196(cs_sles_copy_t) (const void *context);
197
198/*----------------------------------------------------------------------------
199 * Function pointer for destruction of a linear system solver context.
200 *
201 * This function should free all context data, and will be called for each
202 * system when cs_sles_finalize() is called.
203 *
204 * parameters:
205 * context <-> pointer to solver context
206 *----------------------------------------------------------------------------*/
207
208typedef void
209(cs_sles_destroy_t) (void **context);
210
211/*----------------------------------------------------------------------------
212 * Function pointer for handling of non-convegence when solving
213 * a linear system.
214 *
215 * Such a function is optional, and may be used for a variety of purposes,
216 * such as logging, postprocessing, re-trying with different parameters,
217 * aborting the run, or any combination thereof.
218 *
219 * An error handler may be associated with a given solver using
220 * cs_sles_set_error_handler(), in which case it will be called whenever
221 * convergence fails.
222 *
223 * parameters:
224 * sles <-> pointer to solver object
225 * state <-- convergence status
226 * a <-- matrix
227 * rhs <-- Right hand side
228 * vx <-- System solution
229 *
230 * returns:
231 * true if solve should be re-executed, false otherwise
232 *----------------------------------------------------------------------------*/
233
234typedef bool
237 const cs_matrix_t *a,
238 const cs_real_t *rhs,
239 cs_real_t *vx);
240
241/*----------------------------------------------------------------------------
242 * Function pointer for the default definition of a sparse
243 * linear equation solver
244 *
245 * The function may be associated using cs_sles_set_default_define(), so
246 * that it may provide a definition that will be used when
247 * cs_sles_setup() or cs_sles_solve() is used for a system for which
248 * no matching call to cs_sles_define() has been done.
249 *
250 * The function should call cs_sles_define() with arguments f_id
251 * and name, and appropriately chosen function pointers.
252 *
253 * A pointer to the matrix of the system to be solved is also provided,
254 * so that the corresponding information may be used to better choose
255 * defaults.
256 *
257 * parameters:
258 * f_id <-- associated field id, or < 0
259 * name <-- associated name if f_id < 0, or NULL
260 * a <-- Matrix
261 *----------------------------------------------------------------------------*/
262
263typedef void
264(cs_sles_define_t) (int f_id,
265 const char *name,
266 const cs_matrix_t *a);
267
268/*----------------------------------------------------------------------------
269 * Function pointer for the default definition of a sparse
270 * linear equation solver's verbosity
271 *
272 * The function may be associated using cs_sles_set_default_verbosity(), so
273 * that it may provide a definition that will be used when
274 * cs_sles_default_verbosity() is called.
275 *
276 * parameters:
277 * f_id <-- associated field id, or < 0
278 * name <-- associated name if f_id < 0, or NULL
279 *
280 * returns:
281 * default verbosity value
282 *----------------------------------------------------------------------------*/
283
284typedef int
285(cs_sles_verbosity_t) (int f_id,
286 const char *name);
287
288/*============================================================================
289 * Global variables
290 *============================================================================*/
291
292/*=============================================================================
293 * Public function prototypes for Fortran API
294 *============================================================================*/
295
296/*=============================================================================
297 * Public function prototypes
298 *============================================================================*/
299
300/*----------------------------------------------------------------------------*/
304/*----------------------------------------------------------------------------*/
305
306void
307cs_sles_set_epzero(double new_value);
308
309/*----------------------------------------------------------------------------*/
316/*----------------------------------------------------------------------------*/
317
318double
320
321/*----------------------------------------------------------------------------
322 * \brief Initialize sparse linear equation solver API.
323 *----------------------------------------------------------------------------*/
324
325void
327
328/*----------------------------------------------------------------------------
329 * \brief Finalize sparse linear equation solver API.
330 *----------------------------------------------------------------------------*/
331
332void
333cs_sles_finalize(void);
334
335/*----------------------------------------------------------------------------*/
341/*----------------------------------------------------------------------------*/
342
343void
344cs_sles_log(cs_log_t log_type);
345
346/*----------------------------------------------------------------------------*/
347/*
348 * \brief Return pointer to linear system object, based on matching field id or
349 * system name.
350 *
351 * If this system did not previously exist, nullptr is returned.
352 *
353 * \param[in] f_id associated field id, or < 0
354 * \param[in] name associated name if f_id < 0, or NULL
355 *
356 * \return pointer to associated linear system object, or NULL
357 */
358/*----------------------------------------------------------------------------*/
359
360cs_sles_t *
361cs_sles_find(int f_id,
362 const char *name);
363
364/*----------------------------------------------------------------------------*/
365/*
366 * \brief Return pointer to linear system object, based on matching field id or
367 * system name.
368 *
369 * If this system did not previously exist, it is created and added to
370 * to the list of "known" systems. In this case, it will be usable
371 * only if cs_sles_define() is called for the same field id and name
372 * (in which case calling the present function is redundant), or if
373 * cs_sles_set_sefault_define() has been previously used to define
374 * the default solver policy.
375 *
376 * \param[in] f_id associated field id, or < 0
377 * \param[in] name associated name if f_id < 0, or NULL
378 *
379 * \return pointer to associated linear system object, or NULL
380 */
381/*----------------------------------------------------------------------------*/
382
383cs_sles_t *
384cs_sles_find_or_add(int f_id,
385 const char *name);
386
387/*----------------------------------------------------------------------------*/
388/*
389 * \brief Temporarily replace field id with name for matching calls
390 * to \ref cs_sles_setup, \ref cs_sles_solve, \ref cs_sles_free, and other
391 * operations involving access through a field id.
392 *
393 * This function is provided to allow some peculiar calling sequences,
394 * in which \ref cs_equation_iterative_solve_scalar is called with a given
395 * field id, but specific solver options must still be set.
396 * In the future, a cleaner method to handle those exceptional cases
397 * would be preferred. As such, only a stack depth of 1 is allowed.
398 *
399 * \param[in] f_id associated field id, or < 0
400 * \param[in] name associated name if f_id < 0, or NULL
401 */
402/*----------------------------------------------------------------------------*/
403
404void
405cs_sles_push(int f_id,
406 const char *name);
407
408/*----------------------------------------------------------------------------*/
409/*
410 * \brief Restore behavior temporarily modified by \ref cs_sles_push.
411 *
412 * \deprecated This function matches \ref cs_sles_push, which is deprecated.
413 *
414 * \param[in] f_id associated field id, or < 0
415 */
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/*----------------------------------------------------------------------------*/
504/*
505 * \brief Activate postprocessing output for a given linear equation solver.
506 *
507 * This allows the output of the residual at the end of each solution
508 * series, using a single postprocessing writer.
509 * By default, no output is activated.
510 *
511 * \param[in, out] sles pointer to solver object
512 * \param[in] writer_id id of the writer
513 */
514/*----------------------------------------------------------------------------*/
515
516void
518 int writer_id);
519
520/*----------------------------------------------------------------------------*/
521/*
522 * \brief Return the id of the associated writer if postprocessing output
523 * is active for a given linear equation solver.
524 *
525 * \param[in] sles pointer to solver object
526 *
527 * \return id od associated writer, or 0
528 */
529/*----------------------------------------------------------------------------*/
530
531int
533
534/*----------------------------------------------------------------------------*/
535/*
536 * \brief Return type name of solver context.
537 *
538 * The returned string is intended to help determine which type is associated
539 * with the void * pointer returned by \ref cs_sles_get_context for a given
540 * solver definition, so as to be able to call additional specific functions
541 * beyond the generic functions assigned to a cs_sles_t object.
542 *
543 * If no type_name string was associated to the solver upon its definition by
544 * \ref cs_sles_define, or it has not been defined yet, the string returned
545 * is "<undefined>". It is recommended the type name provided
546 * \ref cs_sles_define directly relate to the associated structure type
547 * (for example, "cs_sles_it_t" or "cs_multigrid_t").
548 *
549 * \param[in] sles pointer to solver object
550 *
551 * \return pointer to linear system solver specific type name
552 */
553/*----------------------------------------------------------------------------*/
554
555const char *
557
558/*----------------------------------------------------------------------------*/
559/*
560 * \brief Return pointer to solver context structure pointer.
561 *
562 * The context structure depends on the type of solver used, which may in
563 * turn be determined by the string returned by cs_sles_get_type().
564 * If may be used by appropriate functions specific to that type.
565 *
566 * \param[in] sles pointer to solver object
567 *
568 * \return pointer to solver-specific linear system info and context
569 */
570/*----------------------------------------------------------------------------*/
571
572void *
574
575/*----------------------------------------------------------------------------*/
576/*
577 * \brief Return field id associated with a given sparse linear equation solver.
578 *
579 * \param[in] sles pointer to solver object
580 *
581 * \return associated field id (or -1 if defined by name)
582 */
583/*----------------------------------------------------------------------------*/
584
585int
586cs_sles_get_f_id(const cs_sles_t *sles);
587
588/*----------------------------------------------------------------------------*/
589/*
590 * \brief Return name associated with a given sparse linear equation solver.
591 *
592 * This is simply a utility function which will return its name argument
593 * if f_id < 0, and the associated field's name or label otherwise.
594 *
595 * \param[in] sles pointer to solver object
596 *
597 * \return pointer to associated linear system object name
598 */
599/*----------------------------------------------------------------------------*/
600
601const char *
602cs_sles_get_name(const cs_sles_t *sles);
603
604/*----------------------------------------------------------------------------*/
605/*
606 * \brief Query if immediate_return ("no-op") is allowed when initial
607 * guess is zero (solve by increments) and the RHS is already zero within the
608 * normalized tolerance criteria.
609 *
610 * \param[in] sles pointer to solver object
611 *
612 * \return true if immediate return is allowed, false if at least one
613 * iteration is required
614 */
615/*----------------------------------------------------------------------------*/
616
617bool
619
620/*----------------------------------------------------------------------------*/
621/*
622 * \brief Indicate if immediate_return ("no-op") is allowed when initial
623 * guess is zero (solve by increments) and the RHS is already zero within the
624 * normalized tolerance criteria.
625 *
626 * \param[in, out] sles pointer to solver object
627 * \param[in] allow_no_op true if immediate return is allowed,
628 * false if at least one iteration is required
629 */
630/*----------------------------------------------------------------------------*/
631
632void
634 bool allow_no_op);
635
636/*----------------------------------------------------------------------------*/
637/*
638 * \brief Setup sparse linear equation solver.
639 *
640 * Use of this function is optional: if a \ref cs_sles_solve is called
641 * for the same system before this function is called, the latter will be
642 * called automatically.
643 *
644 * If no options were previously provided for the matching system,
645 * default options will be used.
646 *
647 * \param[in, out] sles pointer to solver object
648 * \param[in] a matrix
649 */
650/*----------------------------------------------------------------------------*/
651
652void
654 const cs_matrix_t *a);
655
656/*----------------------------------------------------------------------------*/
657/*
658 * \brief General sparse linear system resolution.
659 *
660 * If no options were previously provided for the matching system,
661 * default options will be used.
662 *
663 * Note that if \ref cs_sles_setup was previously called for this
664 * system, and \ref cs_sles_free has not been called since, the matrix
665 * provided should be the same. The optional separation between the
666 * two stages is intended to allow amortizing the cost of setup
667 * over multiple solutions.
668 *
669 * If the initial solution iz zero, setting vx_ini_0 to true may
670 * allow extra optimizations with some solvers, and avoids
671 * requiring an upstream initialization.
672 *
673 * The system is considered to have converged when
674 * residual/r_norm <= precision, residual being the L2 norm of a.vx-rhs.
675 *
676 * \param[in, out] sles pointer to solver object
677 * \param[in] a matrix
678 * \param[in] precision solver precision
679 * \param[in] r_norm residual normalization
680 * \param[in] vx_ini_0 if true, initialize solution to 0
681 * \param[out] n_iter number of "equivalent" iterations
682 * \param[out] residual residual
683 * \param[in] rhs right hand side
684 * \param[in, out] vx system solution
685 * \param[in] aux_size size of aux_vectors (in bytes)
686 * \param aux_vectors optional working area
687 * (internal allocation if NULL)
688 *
689 * \return convergence state
690 */
691/*----------------------------------------------------------------------------*/
692
695 const cs_matrix_t *a,
696 double precision,
697 double r_norm,
698 bool vx_ini_0,
699 int *n_iter,
700 double *residual,
701 const cs_real_t *rhs,
702 cs_real_t *vx,
703 size_t aux_size,
704 void *aux_vectors);
705
706/*----------------------------------------------------------------------------*/
707/*
708 * \brief Free sparse linear equation solver setup.
709 *
710 * This function frees resolution-related data, such as multigrid hierarchy,
711 * preconditioning, and any other temporary arrays or objects required for
712 * resolution, but maintains context information such as that used for
713 * logging (especially performance data).
714 *
715 * \param[in, out] sles pointer to solver object
716 */
717/*----------------------------------------------------------------------------*/
718
719void
721
722/*----------------------------------------------------------------------------*/
723/*
724 * \brief Copy the definition of a sparse linear equation solver to another.
725 *
726 * The intended use of this function is to allow associating different
727 * solvers to related systems, so as to differentiate logging, while using
728 * the same settings by default.
729 *
730 * If the source solver does not provide a \ref cs_sles_copy_t function,
731 * No modification is done to the solver. If the copy function is available,
732 * the context is copied, as are the matching function pointers.
733 *
734 * If previous settings have been defined and used, they are saved as
735 * per \ref cs_sles_define.
736 *
737 * \param[in, out] dest pointer to destination solver object
738 * \param[in] src pointer to source solver object
739 *
740 * \return 0 in case of success, 1 in case of failure
741 */
742/*----------------------------------------------------------------------------*/
743
744int
746 const cs_sles_t *src);
747
748/*----------------------------------------------------------------------------*/
749/*
750 * \brief Associate a convergence error handler to a given sparse linear
751 * equation solver.
752 *
753 * The error will be called whenever convergence fails. To dissassociate
754 * the error handler, this function may be called with \p handler = NULL.
755 *
756 * The association will only be successful if the matching solver
757 * has already been defined.
758 *
759 * \param[in, out] sles pointer to solver object
760 * \param[in] error_handler_func pointer to convergence error
761 * handler function
762 */
763/*----------------------------------------------------------------------------*/
764
765void
767 cs_sles_error_handler_t *error_handler_func);
768
769/*----------------------------------------------------------------------------*/
770/*
771 * \brief Return pointer to default sparse linear solver definition function.
772 *
773 * The associated function will be used to provide a definition when
774 * \ref cs_sles_setup or \ref cs_sles_solve is used for a system for which no
775 * matching call to \ref cs_sles_define has been done.
776 *
777 * \return define_func pointer to default definition function
778 */
779/*----------------------------------------------------------------------------*/
780
783
784/*----------------------------------------------------------------------------*/
785/*
786 * \brief Set default sparse linear solver definition function.
787 *
788 * The provided function will be used to provide a definition when
789 * \ref cs_sles_setup or \ref cs_sles_solve is used for a system for which no
790 * matching call to \ref cs_sles_define has been done.
791 *
792 * \param[in] define_func pointer to default definition function
793 */
794/*----------------------------------------------------------------------------*/
795
796void
798
799/*----------------------------------------------------------------------------*/
800/*
801 * \brief Set default verbosity definition function.
802 *
803 * The provided function will be used to define the verbosity when
804 * \ref cs_sles_find_or_add is called.
805 *
806 * \param[in] verbosity_func pointer to default verbosity function
807 */
808/*----------------------------------------------------------------------------*/
809
810void
812
813/*----------------------------------------------------------------------------*/
814/*
815 * \brief Output default post-processing data for failed system convergence.
816 *
817 * \param[in] name variable name
818 * \param[in] mesh_id id of error output mesh, or 0 if none
819 * \param[in] a linear equation matrix
820 * \param[in] rhs right hand side
821 * \param[in, out] vx current system solution
822 */
823/*----------------------------------------------------------------------------*/
824
825void
826cs_sles_post_error_output_def(const char *name,
827 int mesh_id,
828 const cs_matrix_t *a,
829 const cs_real_t *rhs,
830 cs_real_t *vx);
831
832/*----------------------------------------------------------------------------*/
833/*
834 * \brief Output post-processing variable related to system convergence.
835 *
836 * \param[in] name variable name
837 * \param[in] mesh_id id of error output mesh, or 0 if none
838 * \param[in] location_id mesh location id (cells or vertices)
839 * \param[in] writer_id id of specified associated writer, or
840 * \ref CS_POST_WRITER_ALL_ASSOCIATED for all
841 * \param[in] diag_block_size block size for diagonal
842 * \param[in, out] var variable values
843 */
844/*----------------------------------------------------------------------------*/
845
846void
847cs_sles_post_output_var(const char *name,
848 int mesh_id,
849 int location_id,
850 int writer_id,
851 int diag_block_size,
852 cs_real_t var[]);
853
854/*----------------------------------------------------------------------------*/
855/*
856 * \brief Return base name associated to a field id, name couple.
857 *
858 * This is simply a utility function which will return its name argument
859 * if f_id < 0, and the associated field's name or label otherwise.
860 *
861 * \param[in] f_id associated field id, or < 0
862 * \param[in] name associated name if f_id < 0, or NULL
863 *
864 * \return pointer to associated linear system object, or NULL
865 */
866/*----------------------------------------------------------------------------*/
867
868const char *
869cs_sles_base_name(int f_id,
870 const char *name);
871
872/*----------------------------------------------------------------------------*/
873/*
874 * \brief Return name associated to a field id, name couple.
875 *
876 * \param[in] f_id associated field id, or < 0
877 * \param[in] name associated name if f_id < 0, or NULL
878 *
879 * \return pointer to name associated to the field id, name couple
880 */
881/*----------------------------------------------------------------------------*/
882
883const char *
884cs_sles_name(int f_id,
885 const char *name);
886
887/*----------------------------------------------------------------------------*/
888
890
891#endif /* __CS_SLES_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
cs_log_t
Definition: cs_log.h:48
struct _cs_matrix_t cs_matrix_t
Definition: cs_matrix.h:110
void() cs_sles_destroy_t(void **context)
Definition: cs_sles.h:209
cs_sles_convergence_state_t cs_sles_solve(cs_sles_t *sles, const cs_matrix_t *a, double precision, double r_norm, bool vx_ini_0, int *n_iter, double *residual, const cs_real_t *rhs, cs_real_t *vx, size_t aux_size, void *aux_vectors)
General sparse linear system resolution.
Definition: cs_sles.cpp:1736
void * cs_sles_get_context(cs_sles_t *sles)
Return pointer to solver context structure pointer.
Definition: cs_sles.cpp:1568
void cs_sles_set_default_define(cs_sles_define_t *define_func)
Set default sparse linear solver definition function.
Definition: cs_sles.cpp:2010
cs_sles_convergence_state_t
Definition: cs_sles.h:56
@ CS_SLES_MAX_ITERATION
Definition: cs_sles.h:60
@ CS_SLES_BREAKDOWN
Definition: cs_sles.h:59
@ CS_SLES_DIVERGED
Definition: cs_sles.h:58
@ CS_SLES_CONVERGED
Definition: cs_sles.h:62
@ CS_SLES_ITERATING
Definition: cs_sles.h:61
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:174
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.cpp:1915
void cs_sles_setup(cs_sles_t *sles, const cs_matrix_t *a)
Setup sparse linear equation solver.
Definition: cs_sles.cpp:1663
void cs_sles_log(cs_log_t log_type)
Log sparse linear equation solver info.
Definition: cs_sles.cpp:1032
double cs_sles_get_epzero(void)
Get the current threshold value used in the detection of immediate exit.
Definition: cs_sles.cpp:960
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.cpp:2158
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:196
const char * cs_sles_get_name(const cs_sles_t *sles)
Return name associated with a given sparse linear equation solver.
Definition: cs_sles.cpp:1603
int cs_sles_get_verbosity(cs_sles_t *sles)
Get the verbosity for a given linear equation solver.
Definition: cs_sles.cpp:1469
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:285
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_pop(int f_id)
Restore behavior temporarily modified by cs_sles_push.
Definition: cs_sles.cpp:1344
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.cpp:1396
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.cpp:1308
void cs_sles_free(cs_sles_t *sles)
Free sparse linear equation solver setup.
Definition: cs_sles.cpp:1884
const char * cs_sles_name(int f_id, const char *name)
Return name associated to a field id, name couple.
Definition: cs_sles.cpp:2299
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.cpp:1516
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 *residual, const cs_real_t *rhs, cs_real_t *vx_ini, cs_real_t *vx, size_t aux_size, void *aux_vectors)
Function pointer for resolution of a linear system.
Definition: cs_sles.h:131
cs_sles_define_t * cs_sles_get_default_define(void)
Return pointer to default sparse linear solver definition function.
Definition: cs_sles.cpp:1992
void cs_sles_set_verbosity(cs_sles_t *sles, int verbosity)
Set the verbosity for a given linear equation solver.
Definition: cs_sles.cpp:1450
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.cpp:2274
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.cpp:1622
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.cpp:1273
void cs_sles_set_epzero(double new_value)
Set the threshold value used in the detection of immediate exit.
Definition: cs_sles.cpp:945
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.cpp:1208
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.cpp:1584
void cs_sles_finalize(void)
Finalize sparse linear equation solver API.
Definition: cs_sles.cpp:992
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.cpp:2045
const char * cs_sles_get_type(cs_sles_t *sles)
Return type name of solver context.
Definition: cs_sles.cpp:1548
void() cs_sles_free_t(void *context)
Function pointer for freeing of a linear system's context data.
Definition: cs_sles.h:159
void cs_sles_set_default_verbosity(cs_sles_verbosity_t *verbosity_func)
Set default verbosity definition function.
Definition: cs_sles.cpp:2027
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.cpp:1972
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:235
struct _cs_sles_t cs_sles_t
Definition: cs_sles.h:68
void cs_sles_initialize(void)
Initialize sparse linear equation solver API.
Definition: cs_sles.cpp:972
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.cpp:1488
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:264
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.cpp:1640