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 context using
219  * cs_sles_set_error_handler(), in which case it will be called whenever
220  * convergence fails.
221  *
222  * Remark: in the advent of re-trying with different parameters,
223  * it is recommended that either the parameters be sufficiently similar that
224  * performance logging will not be affected, so the info reported to
225  * the user is not biased. Complex strategies involving different
226  * solver types should not be based on the use of an error handler, but
227  * built-into the solver function, with appropriate performance logging
228  * (though they may use an error handler in case of final failure).
229  *
230  * parameters:
231  * context <-> pointer to solver context
232  * state <-- convergence state
233  * name <-- name of linear system
234  * a <-- matrix
235  * rotation_mode <-- Halo update option for rotational periodicity
236  * rhs <-- Right hand side
237  * vx <-- System solution
238  *----------------------------------------------------------------------------*/
239 
240 typedef void
241 (cs_sles_error_handler_t) (void *context,
243  const char *name,
244  const cs_matrix_t *a,
245  cs_halo_rotation_t rotation_mode,
246  const cs_real_t *rhs,
247  cs_real_t *vx);
248 
249 /*----------------------------------------------------------------------------
250  * Function pointer for the default definition of a sparse
251  * linear equation solver
252  *
253  * The function may be associated using cs_sles_set_default_define(), so
254  * that it may provide a definition that will be used when
255  * cs_sles_setup() or cs_sles_solve() is used for a system for which
256  * no matching call to cs_sles_define() has been done.
257  *
258  * The function should call cs_sles_define() with arguments f_id
259  * and name, and appropriately chosen function pointers.
260  *
261  * A pointer to the matrix of the system to be solved is also provided,
262  * so that the corresponding information may be used to better choose
263  * defaults.
264  *
265  * parameters:
266  * f_id <-- associated field id, or < 0
267  * name <-- associated name if f_id < 0, or NULL
268  * a <-- Matrix
269  *----------------------------------------------------------------------------*/
270 
271 typedef void
272 (cs_sles_define_t) (int f_id,
273  const char *name,
274  const cs_matrix_t *a);
275 
276 /*----------------------------------------------------------------------------
277  * Function pointer for the default definition of a sparse
278  * linear equation solver's verbosity
279  *
280  * The function may be associated using cs_sles_set_default_verbosity(), so
281  * that it may provide a definition that will be used when
282  * cs_sles_default_verbosity() is called.
283  *
284  * parameters:
285  * f_id <-- associated field id, or < 0
286  * name <-- associated name if f_id < 0, or NULL
287  *
288  * returns:
289  * default verbosity value
290  *----------------------------------------------------------------------------*/
291 
292 typedef int
294  const char *name);
295 
296 /*============================================================================
297  * Global variables
298  *============================================================================*/
299 
300 /*=============================================================================
301  * Public function prototypes for Fortran API
302  *============================================================================*/
303 
304 /*=============================================================================
305  * Public function prototypes
306  *============================================================================*/
307 
308 /*----------------------------------------------------------------------------
309  * Initialize sparse linear equation solver API.
310  *----------------------------------------------------------------------------*/
311 
312 void
313 cs_sles_initialize(void);
314 
315 /*----------------------------------------------------------------------------
316  * Finalize sparse linear equation solver API.
317  *----------------------------------------------------------------------------*/
318 
319 void
320 cs_sles_finalize(void);
321 
322 /*----------------------------------------------------------------------------
323  * Log sparse linear equation solver info
324  *
325  * parameters:
326  * log_type <-- log type (CS_LOG_SETUP or CS_LOG_PERFORMANCE)
327  *----------------------------------------------------------------------------*/
328 
329 void
330 cs_sles_log(cs_log_t log_type);
331 
332 /*----------------------------------------------------------------------------
333  * Return pointer to linear system object, based on matching field id or
334  * system name.
335  *
336  * If this system did not previously exist, NULL is returned.
337  *
338  * parameters:
339  * f_id <-- associated field id, or < 0
340  * name <-- associated name if f_id < 0, or NULL
341  *
342  * returns:
343  * pointer to associated linear system object, or NULL
344  *----------------------------------------------------------------------------*/
345 
346 cs_sles_t *
347 cs_sles_find(int f_id,
348  const char *name);
349 
350 /*----------------------------------------------------------------------------
351  * Return pointer to linear system object, based on matching field id or
352  * system name.
353  *
354  * If this system did not previously exist, it is created and added to
355  * to the list of "known" systems. In this case, it will be usable
356  * only if cs_sles_define() is called for the same field id and name
357  * (in which case calling the present function is redundant), or if
358  * cs_sles_set_sefault_define() has been previously used to define
359  * the default solver policy.
360  *
361  * parameters:
362  * f_id <-- associated field id, or < 0
363  * name <-- associated name if f_id < 0, or NULL
364  *
365  * returns:
366  * pointer to associated linear system object.
367  *----------------------------------------------------------------------------*/
368 
369 cs_sles_t *
370 cs_sles_find_or_add(int f_id,
371  const char *name);
372 
373 /*----------------------------------------------------------------------------
374  * Temporarily replace field id with name for matching calls
375  * to cs_sles_setup(), cs_sles_solve(), cs_sles_free(), and other operations
376  * involving access through a field id.
377  *
378  * Deprecated.
379  * This function is provided to allow some peculiar
380  * calling sequences, in which codits is called with a nonzero
381  * ivar value, but specific solver options must still be set.
382  * In the future, a more consistent mechanism (using a zero ivar
383  * value or designing a cleaner method to handle those exceptional cases)
384  * is preferred. As such, only a stack depth of 1 is allowed.
385  *
386  * parameters:
387  * f_id associated field id, or < 0
388  * name associated name if f_id < 0, or NULL
389  */
390 /*----------------------------------------------------------------------------*/
391 
392 void
393 cs_sles_push(int f_id,
394  const char *name);
395 
396 /*----------------------------------------------------------------------------
397  * Restore behavior temporarily modified by cs_slesh_push().
398  *
399  * Deprecated.
400  * This function matches cs_sles_push(), which is deprecated.
401  *
402  * parameters:
403  * f_id associated field id, or < 0
404  *----------------------------------------------------------------------------*/
405 
406 void
407 cs_sles_pop(int f_id);
408 
409 /*----------------------------------------------------------------------------
410  * Define sparse linear equation solver for a given field or equation name.
411  *
412  * If this system did not previously exist, it is added to the list of
413  * "known" systems. If it existed previously,
414  *
415  * The context pointer is used to point to a structure adapted to the function
416  * pointers given here, and combined with those functions, allows using
417  * both built-in, external, or user-defined solvers.
418  *
419  * It is recommended the context type name provided here directly relate
420  * to the associated structure type (for example, "cs_sles_it_t" or
421  * "cs_multigrid_t").
422  *
423  * parameters:
424  * f_id <-- associated field id, or < 0
425  * name <-- associated name if f_id < 0, or NULL
426  * context <-> pointer to solver context management structure
427  * (cs_sles subsystem becomes owner)
428  * setup_func <-- pointer to system setup function
429  * solve_func <-- pointer to system solution function
430  * (also calls setup_func if not done yet
431  * or free_func called since last solve)
432  * free_func <-- pointer function freeing system setup
433  * log_func <-- pointer to system info logging function
434  * (optional, but recommended)
435  * copy_func <-- pointer to settings copy function (optional)
436  * destroy_func <-- pointer to function destroying solver context
437  * (called with cs_sles_finalize() or with
438  * a new call to this function for the same
439  * system)
440  *
441  * returns:
442  * pointer to associated linear system object
443  *----------------------------------------------------------------------------*/
444 
445 cs_sles_t *
446 cs_sles_define(int f_id,
447  const char *name,
448  void *context,
449  const char *type_name,
450  cs_sles_setup_t *setup_func,
451  cs_sles_solve_t *solve_func,
452  cs_sles_free_t *free_func,
453  cs_sles_log_t *log_func,
454  cs_sles_copy_t *copy_func,
455  cs_sles_destroy_t *destroy_func);
456 
457 /*----------------------------------------------------------------------------
458  * Set the verbosity for a given linear equation solver.
459  *
460  * This verbosity will be used by cs_sles_setup() and cs_sles_solve().
461  *
462  * By default, the verbosity is set to 0, or the value returned by the
463  * function set with cs_sles_set_default_define().
464  *
465  * parameters
466  * sles <-> pointer to solver object
467  * verbosity <-- verbosity level
468  *----------------------------------------------------------------------------*/
469 
470 void
472  int verbosity);
473 
474 /*----------------------------------------------------------------------------
475  * Return type name of solver context.
476  *
477  * The returned string is intended to help determine which type is associated
478  * with the void * pointer returned by cs_sles_get_context() for a given
479  * solver definition, so as to be able to call additional specific functions
480  * beyond the generic functions assigned to a cs_sles_t object.
481  *
482  * If no type_name string was associated to the solver upon its definition by
483  * cs_sles_define(), or it has not been defined yet, the string returned
484  * is "<undefined>". It is recommended the type name provided
485  * cs_sles_define() directly relate to the associated structure type
486  * (for example, "cs_sles_it_t" or "cs_multigrid_t").
487  *
488  * parameters
489  * sles pointer to solver object
490  *
491  * returns:
492  * pointer to linear system solver specific type name
493  *----------------------------------------------------------------------------*/
494 
495 const char *
497 
498 /*----------------------------------------------------------------------------
499  * Return pointer to solver context structure pointer.
500  *
501  * The context structure depends on the type of solver used, which may in
502  * turn be determined by the string returned by cs_sles_get_type().
503  * If may be used by appropriate functions specific to that type.
504  *
505  * parameters
506  * sles pointer to solver object
507  *
508  * returns:
509  * pointer to solver-specific linear system info and context
510 *----------------------------------------------------------------------------*/
511 
512 void *
514 
515 /*----------------------------------------------------------------------------
516  * Setup sparse linear equation solver.
517  *
518  * If no options were previously provided for the matching system,
519  * default options will be used.
520  *
521  * parameters:
522  * sles <-> pointer to solver object
523  * a <-- matrix
524  *
525  * returns:
526  * convergence state
527  *----------------------------------------------------------------------------*/
528 
529 void
531  const cs_matrix_t *a);
532 
533 /*----------------------------------------------------------------------------
534  * General sparse linear system resolution.
535  *
536  * If no options were previously provided for the matching system,
537  * default options will be used.
538  *
539  * Note that if cs_sles_setup() was previously called for this
540  * system, and cs_sles_free() has not been called since, the matrix
541  * provided should be the same. The optional separation between the
542  * two stages is intended to allow amortizing the cost of setup
543  * over multiple solutions.
544  *
545  * The system is considered to have converged when
546  * residue/r_norm <= precision, residue being the L2 norm of a.vx-rhs.
547  *
548  * parameters:
549  * sles <-> pointer to solver object
550  * a <-- matrix
551  * rotation_mode <-- halo update option for rotational periodicity
552  * precision <-- solver precision
553  * r_norm <-- residue normalization
554  * n_iter --> number of "equivalent" iterations
555  * residue --> residue
556  * rhs <-- right hand side
557  * vx <-> system solution
558  * aux_size <-- size of aux_vectors (in bytes)
559  * aux_vectors --- optional working area (internal allocation if NULL)
560  *
561  * returns:
562  * convergence state
563  *----------------------------------------------------------------------------*/
564 
567  const cs_matrix_t *a,
568  cs_halo_rotation_t rotation_mode,
569  double precision,
570  double r_norm,
571  int *n_iter,
572  double *residue,
573  const cs_real_t *rhs,
574  cs_real_t *vx,
575  size_t aux_size,
576  void *aux_vectors);
577 
578 /*----------------------------------------------------------------------------
579  * Free sparse linear equation solver setup.
580  *
581  * This function frees resolution-related data, such as multigrid hierarchy,
582  * preconditioning, and any other temporary arrays or objects required for
583  * resolution, but maintains context information such as that used for
584  * logging (especially performance data).
585  *
586  * parameters:
587  * sles <-> pointer to solver object
588  *----------------------------------------------------------------------------*/
589 
590 void
591 cs_sles_free(cs_sles_t *sles);
592 
593 /*----------------------------------------------------------------------------
594  * Copy the definition of a sparse linear equation solver to another.
595  *
596  * The intended use of this function is to allow associating different
597  * solvers to related systems, so as to differentiate logging, while using
598  * the same setttings by default.
599  *
600  * If the source solver does not provide a cs_sles_copy_t function,
601  * No modification is done to the solver. If the copy function is available,
602  * the context is copied, as are the matching function pointers.
603  *
604  * If previous settings have been defined and used, they are saved as
605  * per cs_sles_define().
606  *
607  * parameters:
608  * dest <-> pointer to destination solver object
609  * src <-- pointer to source solver object
610  *
611  * returns:
612  * 0 in case of success, 1 in case of failure
613  *----------------------------------------------------------------------------*/
614 
615 int
616 cs_sles_copy(cs_sles_t *dest,
617  const cs_sles_t *src);
618 
619 /*----------------------------------------------------------------------------
620  * Associate a convergence error handler to a given sparse linear
621  * equation solver.
622  *
623  * The error will be called whenever convergence fails. To dissassociate
624  * the error handler, this function may be called with \p handler = NULL.
625  *
626  * The association will only be successful if the matching solver
627  * has already been defined.
628  *
629  * parameters:
630  * sles <-> pointer to solver object
631  * error_handler_func <-- pointer to convergence error handler function
632  *----------------------------------------------------------------------------*/
633 
634 void
636  cs_sles_error_handler_t *error_handler_func);
637 
638 /*----------------------------------------------------------------------------
639  * Set default sparse linear solver definition function.
640  *
641  * The provided function will be used to provide a definition when
642  * cs_sles_setup() or cs_sles_solve() is used for a system for which no
643  * matching call to cs_sles_define() has been done.
644  *
645  * parameters:
646  * define_func <-- pointer to default definition function
647  *----------------------------------------------------------------------------*/
648 
649 void
651 
652 /*----------------------------------------------------------------------------
653  * Set default verbosity definition function.
654  *
655  * The provided function will be used to define the verbosity when
656  * cs_sles_default_verbosity() is called.
657  *
658  * parameters:
659  * verbosity_func <-- pointer to default verbosity function
660  *----------------------------------------------------------------------------*/
661 
662 void
664 
665 /*----------------------------------------------------------------------------
666  * Test if a linear system needs solving or if the residue is already
667  * within convergence criteria.
668  *
669  * parameters:
670  * solver_name <-- name of the solver calling the test
671  * system_name <-- name of the linear system tested
672  * verbosity <-- verbosity level
673  * precision <-- solver precision
674  * r_norm <-- residue normalization
675  * residue <-- residue
676  *
677  * returns:
678  * 1 if solving is required, 0 if the residue is already zero within
679  * tolerance criteria (precision of residue normalization)
680  *----------------------------------------------------------------------------*/
681 
682 int
683 cs_sles_needs_solving(const char *solver_name,
684  const char *system_name,
685  int verbosity,
686  double precision,
687  double r_norm,
688  double residue);
689 
690 /*----------------------------------------------------------------------------
691  * Output default post-processing data for failed system convergence.
692  *
693  * parameters:
694  * name <-- variable name
695  * mesh_id <-- id of error output mesh, or 0 if none
696  * rotation_mode <-- halo update option for rotational periodicity
697  * a <-- linear equation matrix
698  * rhs <-- right hand side
699  * vx <-> current system solution
700  *----------------------------------------------------------------------------*/
701 
702 void
703 cs_sles_post_error_output_def(const char *name,
704  int mesh_id,
705  cs_halo_rotation_t rotation_mode,
706  const cs_matrix_t *a,
707  const cs_real_t *rhs,
708  cs_real_t *vx);
709 
710 /*----------------------------------------------------------------------------
711  * Output post-processing variable for failed system convergence.
712  *
713  * parameters:
714  * name <-- variable name
715  * diag_block_size <-- block size for diagonal
716  * mesh_id <-- id of error output mesh, or 0 if none
717  * var <-> variable values
718  *----------------------------------------------------------------------------*/
719 
720 void
721 cs_sles_post_error_output_var(const char *name,
722  int mesh_id,
723  int diag_block_size,
724  cs_real_t *var);
725 
726 /*----------------------------------------------------------------------------*
727  * Return base name associated to a field id, name couple.
728  *
729  * This is simply a utility function which will return its name argument
730  * if f_id < 0, and the associated field's name or label otherwise.
731  *
732  * parameters:
733  * f_id <-- associated field id, or < 0
734  * name <-- associated name if f_id < 0, or NULL
735  *
736  * returns:
737  * pointer to associated linear system object, or NULL
738  *----------------------------------------------------------------------------*/
739 
740 const char *
741 cs_sles_base_name(int f_id,
742  const char *name);
743 
744 /*----------------------------------------------------------------------------*/
745 
747 
748 #endif /* __CS_SLES_H__ */
cs_sles_t * cs_sles_find_or_add(int f_id, const char *name)
Return pointer to linear system object, based on matching field id or system name.
Definition: cs_sles.c:1026
void cs_sles_log(cs_log_t log_type)
Log sparse linear equation solver info.
Definition: cs_sles.c:822
const char * cs_sles_get_type(cs_sles_t *sles)
Return type name of solver context.
Definition: cs_sles.c:1233
cs_halo_rotation_t
Definition: cs_halo.h:59
void cs_sles_setup(cs_sles_t *sles, const cs_matrix_t *a)
Setup sparse linear equation solver.
Definition: cs_sles.c:1275
void cs_sles_initialize(void)
Initialize sparse linear equation solver API.
Definition: cs_sles.c:777
void() cs_sles_define_t(int f_id, const char *name, const cs_matrix_t *a)
Function pointer for the default definition of a sparse linear equation solver.
Definition: cs_sles.h:272
Definition: cs_sles.h:60
#define BEGIN_C_DECLS
Definition: cs_defs.h:419
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:1544
void() cs_sles_free_t(void *context)
Function pointer for freeing of a linear system&#39;s context data.
Definition: cs_sles.h:158
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:961
void cs_sles_set_default_verbosity(cs_sles_verbosity_t *verbosity_func)
Set default verbosity definition function.
Definition: cs_sles.c:1525
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:1064
int cs_sles_needs_solving(const char *solver_name, const char *system_name, int verbosity, double precision, double r_norm, double residue)
void cs_sles_free(cs_sles_t *sles)
Free sparse linear equation solver setup.
Definition: cs_sles.c:1400
struct _cs_matrix_t cs_matrix_t
Definition: cs_matrix.h:86
void * cs_sles_get_context(cs_sles_t *sles)
Return pointer to solver context structure pointer.
Definition: cs_sles.c:1253
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:1488
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
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:788
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:1635
double precision, save a
Definition: cs_fuel_incl.f90:146
int cs_sles_copy(cs_sles_t *dest, const cs_sles_t *src)
Copy the definition of a sparse linear equation solver to another.
Definition: cs_sles.c:1431
void() cs_sles_destroy_t(void **context)
Definition: cs_sles.h:208
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:293
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
cs_log_t
Definition: cs_log.h:47
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
#define END_C_DECLS
Definition: cs_defs.h:420
double cs_real_t
Definition: cs_defs.h:296
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:1722
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:1508
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_set_verbosity(cs_sles_t *sles, int verbosity)
Set the verbosity for a given linear equation solver.
Definition: cs_sles.c:1205
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:1152
void cs_sles_pop(int f_id)
Restore behavior temporarily modified by cs_sles_push.
Definition: cs_sles.c:1100
void() cs_sles_error_handler_t(void *context, cs_sles_convergence_state_t state, const char *name, const cs_matrix_t *a, cs_halo_rotation_t rotation_mode, const cs_real_t *rhs, cs_real_t *vx)
Function pointer for handling of non-convergence when solving a linear system.
Definition: cs_sles.h:241
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:1323