7.2
general documentation
cs_sles_petsc.h
Go to the documentation of this file.
1 #ifndef __CS_SLES_PETSC_H__
2 #define __CS_SLES_PETSC_H__
3 
4 /*============================================================================
5  * Sparse Linear Equation Solvers using PETSc
6  *============================================================================*/
7 
8 /*
9  This file is part of code_saturne, a general-purpose CFD tool.
10 
11  Copyright (C) 1998-2022 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 
38 /*----------------------------------------------------------------------------*/
39 
41 
42 /*============================================================================
43  * Macro definitions
44  *============================================================================*/
45 
46 /*============================================================================
47  * Type definitions
48  *============================================================================*/
49 
50 /*----------------------------------------------------------------------------
51  * Function pointer for user settings of a PETSc KSP solver setup.
52  *
53  * This function is called at the end of the setup stage for a KSP solver.
54  *
55  * Note that using the advanced KSPSetPostSolve and KSPSetPreSolve functions,
56  * this also allows setting further function pointers for pre and post-solve
57  * operations (see the PETSc documentation).
58  *
59  * Note: if the context pointer is non-NULL, it must point to valid data
60  * when the selection function is called so that value or structure should
61  * not be temporary (i.e. local);
62  *
63  * parameters:
64  * context <-> pointer to optional (untyped) value or structure
65  * ksp <-> pointer to PETSc KSP context
66  *----------------------------------------------------------------------------*/
67 
68 typedef void
69 (cs_sles_petsc_setup_hook_t) (void *context,
70  void *ksp);
71 
72 /* Iterative linear solver context (opaque) */
73 
74 typedef struct _cs_sles_petsc_t cs_sles_petsc_t;
75 
76 /*============================================================================
77  * Global variables
78  *============================================================================*/
79 
80 /*=============================================================================
81  * User function prototypes
82  *============================================================================*/
83 
84 /*----------------------------------------------------------------------------
85  * Function pointer for user settings of a PETSc KSP solver setup.
86  *
87  * This function is called at the end of the setup stage for a KSP solver.
88  *
89  * Note that using the advanced KSPSetPostSolve and KSPSetPreSolve functions,
90  * this also allows setting further function pointers for pre and post-solve
91  * operations (see the PETSc documentation).
92  *
93  * Note: if the context pointer is non-NULL, it must point to valid data
94  * when the selection function is called so that value or structure should
95  * not be temporary (i.e. local);
96  *
97  * parameters:
98  * context <-> pointer to optional (untyped) value or structure
99  * ksp <-> pointer to PETSc KSP context
100  *----------------------------------------------------------------------------*/
101 
102 void
103 cs_user_sles_petsc_hook(void *context,
104  void *ksp);
105 
106 /*=============================================================================
107  * Public function prototypes
108  *============================================================================*/
109 
110 /*----------------------------------------------------------------------------
111  * Initialize PETSc if needed (calls cs_matrix_petsc_ensure_init).
112  *----------------------------------------------------------------------------*/
113 
114 void
115 cs_sles_petsc_init(void);
116 
117 /*----------------------------------------------------------------------------*/
147 /*----------------------------------------------------------------------------*/
148 
150 cs_sles_petsc_define(int f_id,
151  const char *name,
152  const char *matrix_type,
153  cs_sles_petsc_setup_hook_t *setup_hook,
154  void *context);
155 
156 /*----------------------------------------------------------------------------*/
171 /*----------------------------------------------------------------------------*/
172 
174 cs_sles_petsc_create(const char *matrix_type,
175  cs_sles_petsc_setup_hook_t *setup_hook,
176  void *context);
177 
178 /*----------------------------------------------------------------------------*/
189 /*----------------------------------------------------------------------------*/
190 
191 void *
192 cs_sles_petsc_copy(const void *context);
193 
194 /*----------------------------------------------------------------------------*/
201 /*----------------------------------------------------------------------------*/
202 
203 void
204 cs_sles_petsc_destroy(void **context);
205 
206 /*----------------------------------------------------------------------------*/
216 /*----------------------------------------------------------------------------*/
217 
218 void
219 cs_sles_petsc_setup(void *context,
220  const char *name,
221  const cs_matrix_t *a,
222  int verbosity);
223 
224 /*----------------------------------------------------------------------------*/
245 /*----------------------------------------------------------------------------*/
246 
248 cs_sles_petsc_solve(void *context,
249  const char *name,
250  const cs_matrix_t *a,
251  int verbosity,
252  double precision,
253  double r_norm,
254  int *n_iter,
255  double *residue,
256  const cs_real_t *rhs,
257  cs_real_t *vx,
258  size_t aux_size,
259  void *aux_vectors);
260 
261 /*----------------------------------------------------------------------------*/
272 /*----------------------------------------------------------------------------*/
273 
274 void
275 cs_sles_petsc_free(void *context);
276 
277 /*----------------------------------------------------------------------------*/
293 /*----------------------------------------------------------------------------*/
294 
295 bool
298  const cs_matrix_t *a,
299  const cs_real_t *rhs,
300  cs_real_t *vx);
301 
302 /*----------------------------------------------------------------------------*/
310 /*----------------------------------------------------------------------------*/
311 
312 void
313 cs_sles_petsc_log(const void *context,
314  cs_log_t log_type);
315 
316 /*----------------------------------------------------------------------------
317  * \brief Output the settings of a KSP structure
318  *
319  * \param[in] ksp Krylov SubSpace structure
320  *----------------------------------------------------------------------------*/
321 
322 void
323 cs_sles_petsc_log_setup(void *ksp);
324 
325 /*----------------------------------------------------------------------------*/
334 /*----------------------------------------------------------------------------*/
335 
336 void
337 cs_sles_petsc_set_cvg_criteria(const void *context,
338  double rtol,
339  int max_iter);
340 
341 /*----------------------------------------------------------------------------*/
351 /*----------------------------------------------------------------------------*/
352 
353 const char *
354 cs_sles_petsc_get_mat_type(void *context);
355 
356 /*----------------------------------------------------------------------------*/
362 /*----------------------------------------------------------------------------*/
363 
364 void
366 
367 /*----------------------------------------------------------------------------*/
368 
370 
371 #endif /* __CS_SLES_PETSC_H__ */
cs_sles_convergence_state_t cs_sles_petsc_solve(void *context, const char *name, const cs_matrix_t *a, int verbosity, double precision, double r_norm, int *n_iter, double *residue, const cs_real_t *rhs, cs_real_t *vx, size_t aux_size, void *aux_vectors)
Call PETSc linear equation solver.
Definition: cs_sles_petsc.c:1136
struct _cs_sles_petsc_t cs_sles_petsc_t
Definition: cs_sles_petsc.h:74
void cs_sles_petsc_init(void)
Definition: cs_sles_petsc.c:513
void() cs_sles_petsc_setup_hook_t(void *context, void *ksp)
Function pointer for user settings of a PETSc KSP solver setup.
Definition: cs_sles_petsc.h:69
#define BEGIN_C_DECLS
Definition: cs_defs.h:510
cs_sles_petsc_t * cs_sles_petsc_define(int f_id, const char *name, const char *matrix_type, cs_sles_petsc_setup_hook_t *setup_hook, void *context)
Define and associate a PETSc linear system solver for a given field or equation name.
Definition: cs_sles_petsc.c:551
bool cs_sles_petsc_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 PETSc solver.
Definition: cs_sles_petsc.c:1394
void cs_sles_petsc_library_info(cs_log_t log_type)
Print information on PETSc library.
Definition: cs_sles_petsc.c:1585
double cs_real_t
Floating-point value.
Definition: cs_defs.h:322
struct _cs_matrix_t cs_matrix_t
Definition: cs_matrix.h:110
void cs_sles_petsc_setup(void *context, const char *name, const cs_matrix_t *a, int verbosity)
Setup PETSc linear equation solver.
Definition: cs_sles_petsc.c:743
void * cs_sles_petsc_copy(const void *context)
Create PETSc linear system solver info and context based on existing info and context.
Definition: cs_sles_petsc.c:673
cs_sles_convergence_state_t
Convergence status indicator.
Definition: cs_sles.h:56
struct _cs_sles_t cs_sles_t
Definition: cs_sles.h:68
double precision, save a
Definition: cs_fuel_incl.f90:146
cs_sles_petsc_t * cs_sles_petsc_create(const char *matrix_type, cs_sles_petsc_setup_hook_t *setup_hook, void *context)
Create PETSc linear system solver info and context.
Definition: cs_sles_petsc.c:596
void cs_sles_petsc_free(void *context)
Free PETSc linear equation solver setup context.
Definition: cs_sles_petsc.c:1346
cs_log_t
Definition: cs_log.h:48
const char * cs_sles_petsc_get_mat_type(void *context)
Return matrix type associated with PETSc linear system solver info and context.
Definition: cs_sles_petsc.c:1566
void cs_user_sles_petsc_hook(void *context, void *ksp)
Definition: cs_sles_petsc.c:497
void cs_sles_petsc_log_setup(void *ksp)
Definition: cs_sles_petsc.c:1512
#define END_C_DECLS
Definition: cs_defs.h:511
void cs_sles_petsc_set_cvg_criteria(const void *context, double rtol, int max_iter)
Set some of the convergence criteria for the associated KSP structure.
Definition: cs_sles_petsc.c:1537
void cs_sles_petsc_destroy(void **context)
Destroy PETSc linear system solver info and context.
Definition: cs_sles_petsc.c:697
void cs_sles_petsc_log(const void *context, cs_log_t log_type)
Log sparse linear equation solver info.
Definition: cs_sles_petsc.c:1434