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-2021 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 * PETSc headers
32 *----------------------------------------------------------------------------*/
33
34 #include <petscversion.h>
35 #include <petscconf.h>
36 #include <petscviewer.h>
37 #include <petscksp.h>
38
39 /*----------------------------------------------------------------------------
40 * Local headers
41 *----------------------------------------------------------------------------*/
42
43 #include "cs_base.h"
44 #include "cs_matrix.h"
45 #include "cs_sles.h"
46
47 /*----------------------------------------------------------------------------*/
48
49 BEGIN_C_DECLS
50
51 /*============================================================================
52 * Macro definitions
53 *============================================================================*/
54
55 /*============================================================================
56 * Type definitions
57 *============================================================================*/
58
59 /*----------------------------------------------------------------------------
60 * Function pointer for user settings of a PETSc KSP solver setup.
61 *
62 * This function is called at the end of the setup stage for a KSP solver.
63 *
64 * Note that using the advanced KSPSetPostSolve and KSPSetPreSolve functions,
65 * this also allows setting further function pointers for pre and post-solve
66 * operations (see the PETSc documentation).
67 *
68 * Note: if the context pointer is non-NULL, it must point to valid data
69 * when the selection function is called so that value or structure should
70 * not be temporary (i.e. local);
71 *
72 * parameters:
73 * context <-> pointer to optional (untyped) value or structure
74 * ksp <-> pointer to PETSc KSP context
75 *----------------------------------------------------------------------------*/
76
77 typedef void
78 (cs_sles_petsc_setup_hook_t) (void *context,
79 KSP ksp);
80
81 /* Iterative linear solver context (opaque) */
82
83 typedef struct _cs_sles_petsc_t cs_sles_petsc_t;
84
85 /*============================================================================
86 * Global variables
87 *============================================================================*/
88
89 /*=============================================================================
90 * User function prototypes
91 *============================================================================*/
92
93 /*----------------------------------------------------------------------------
94 * Function pointer for user settings of a PETSc KSP solver setup.
95 *
96 * This function is called at the end of the setup stage for a KSP solver.
97 *
98 * Note that using the advanced KSPSetPostSolve and KSPSetPreSolve functions,
99 * this also allows setting further function pointers for pre and post-solve
100 * operations (see the PETSc documentation).
101 *
102 * Note: if the context pointer is non-NULL, it must point to valid data
103 * when the selection function is called so that value or structure should
104 * not be temporary (i.e. local);
105 *
106 * parameters:
107 * context <-> pointer to optional (untyped) value or structure
108 * ksp <-> pointer to PETSc KSP context
109 *----------------------------------------------------------------------------*/
110
111 void
112 cs_user_sles_petsc_hook(void *context,
113 KSP ksp);
114
115 /*=============================================================================
116 * Static inline public function prototypes
117 *============================================================================*/
118
119 /*----------------------------------------------------------------------------
120 * Initialize PETSc if needed
121 *----------------------------------------------------------------------------*/
122
123 void
124 cs_sles_petsc_init(void);
125
126 /*----------------------------------------------------------------------------
127 * \brief Output the settings of a KSP structure
128 *
129 * \param[in] ksp Krylov SubSpace structure
130 *----------------------------------------------------------------------------*/
131
132 static inline void
cs_sles_petsc_log_setup(KSP ksp)133 cs_sles_petsc_log_setup(KSP ksp)
134 {
135 PetscViewer v;
136
137 PetscViewerCreate(PETSC_COMM_WORLD, &v);
138 PetscViewerSetType(v, PETSCVIEWERASCII);
139 PetscViewerFileSetMode(v, FILE_MODE_APPEND);
140 PetscViewerFileSetName(v, "petsc_setup.log");
141
142 KSPView(ksp, v);
143 PetscViewerDestroy(&v);
144 }
145
146 /*=============================================================================
147 * Public function prototypes
148 *============================================================================*/
149
150 /*----------------------------------------------------------------------------
151 * Define and associate a PETSc linear system solver
152 * for a given field or equation name.
153 *
154 * If this system did not previously exist, it is added to the list of
155 * "known" systems. Otherwise, its definition is replaced by the one
156 * defined here.
157 *
158 * This is a utility function: if finer control is needed, see
159 * cs_sles_define() and cs_sles_petsc_create().
160 *
161 * In case of rotational periodicity for a block (non-scalar) matrix,
162 * the matrix type will be forced to MATSHELL ("shell") regardless
163 * of the option used.
164 *
165 * Note that this function returns a pointer directly to the iterative solver
166 * management structure. This may be used to set further options.
167 * If needed, cs_sles_find() may be used to obtain a pointer to the matching
168 * cs_sles_t container.
169 *
170 * parameters:
171 * f_id <-- associated field id, or < 0
172 * name <-- associated name if f_id < 0, or NULL
173 * matrix_type <-- PETSc matrix type
174 * setup_hook <-- pointer to optional setup epilogue function
175 * context <-> pointer to optional (untyped) value or structure
176 * for setup_hook, or NULL
177 *
178 * returns:
179 * pointer to newly created iterative solver info object.
180 *----------------------------------------------------------------------------*/
181
182 cs_sles_petsc_t *
183 cs_sles_petsc_define(int f_id,
184 const char *name,
185 MatType matrix_type,
186 cs_sles_petsc_setup_hook_t *setup_hook,
187 void *context);
188
189 /*----------------------------------------------------------------------------
190 * Create PETSc linear system solver info and context.
191 *
192 * In case of rotational periodicity for a block (non-scalar) matrix,
193 * the matrix type will be forced to MATSHELL ("shell") regardless
194 * of the option used.
195 *
196 * parameters:
197 * matrix_type <-- PETSc matrix type
198 * setup_hook <-- pointer to optional setup epilogue function
199 * context <-> pointer to optional (untyped) value or structure
200 * for setup_hook, or NULL
201 *
202 * returns:
203 * pointer to newly created solver info object.
204 *----------------------------------------------------------------------------*/
205
206 cs_sles_petsc_t *
207 cs_sles_petsc_create(MatType matrix_type,
208 cs_sles_petsc_setup_hook_t *setup_hook,
209 void *context);
210
211 /*----------------------------------------------------------------------------
212 * Create PETSc linear system solver info and context
213 * based on existing info and context.
214 *
215 * parameters:
216 * context <-- pointer to reference info and context
217 * (actual type: cs_sles_petsc_t *)
218 *
219 * returns:
220 * pointer to newly created solver info object
221 * (actual type: cs_sles_petsc_t *)
222 *----------------------------------------------------------------------------*/
223
224 void *
225 cs_sles_petsc_copy(const void *context);
226
227 /*----------------------------------------------------------------------------
228 * Destroy PETSc linear system solver info and context.
229 *
230 * parameters:
231 * context <-> pointer to PETSc linear solver info
232 * (actual type: cs_sles_petsc_t **)
233 *----------------------------------------------------------------------------*/
234
235 void
236 cs_sles_petsc_destroy(void **context);
237
238 /*----------------------------------------------------------------------------
239 * Setup PETSc linear equation solver.
240 *
241 * parameters:
242 * context <-> pointer to PETSc linear solver info
243 * (actual type: cs_sles_petsc_t *)
244 * name <-- pointer to system name
245 * a <-- associated matrix
246 * verbosity <-- verbosity level
247 *----------------------------------------------------------------------------*/
248
249 void
250 cs_sles_petsc_setup(void *context,
251 const char *name,
252 const cs_matrix_t *a,
253 int verbosity);
254
255 /*----------------------------------------------------------------------------
256 * Call PETSc linear equation solver.
257 *
258 * parameters:
259 * context <-> pointer to PETSc linear solver info
260 * (actual type: cs_sles_petsc_t *)
261 * name <-- pointer to system name
262 * a <-- matrix
263 * verbosity <-- verbosity level
264 * precision <-- solver precision
265 * r_norm <-- residue normalization
266 * n_iter --> number of iterations
267 * residue --> residue
268 * rhs <-- right hand side
269 * vx <-> system solution
270 * aux_size <-- number of elements in aux_vectors (in bytes)
271 * aux_vectors --- optional working area (internal allocation if NULL)
272 *
273 * returns:
274 * convergence state
275 *----------------------------------------------------------------------------*/
276
277 cs_sles_convergence_state_t
278 cs_sles_petsc_solve(void *context,
279 const char *name,
280 const cs_matrix_t *a,
281 int verbosity,
282 double precision,
283 double r_norm,
284 int *n_iter,
285 double *residue,
286 const cs_real_t *rhs,
287 cs_real_t *vx,
288 size_t aux_size,
289 void *aux_vectors);
290
291 /*----------------------------------------------------------------------------
292 * Free PETSc linear equation solver setup context.
293 *
294 * This function frees resolution-related data, such as
295 * buffers and preconditioning but does not free the whole context,
296 * as info used for logging (especially performance data) is maintained.
297
298 * parameters:
299 * context <-> pointer to PETSc linear solver info
300 * (actual type: cs_sles_petsc_t *)
301 *----------------------------------------------------------------------------*/
302
303 void
304 cs_sles_petsc_free(void *context);
305
306 /*----------------------------------------------------------------------------*/
307 /*!
308 * \brief Error handler for PETSc solver.
309 *
310 * In case of divergence or breakdown, this error handler outputs an error
311 * message
312 * It does nothing in case the maximum iteration count is reached.
313
314 * \param[in, out] sles pointer to solver object
315 * \param[in] state convergence state
316 * \param[in] a matrix
317 * \param[in] rhs right hand side
318 * \param[in, out] vx system solution
319 *
320 * \return false (do not attempt new solve)
321 */
322 /*----------------------------------------------------------------------------*/
323
324 bool
325 cs_sles_petsc_error_post_and_abort(cs_sles_t *sles,
326 cs_sles_convergence_state_t state,
327 const cs_matrix_t *a,
328 const cs_real_t *rhs,
329 cs_real_t *vx);
330
331 /*----------------------------------------------------------------------------
332 * Log sparse linear equation solver info.
333 *
334 * parameters:
335 * context <-> pointer to PETSc linear solver info
336 * (actual type: cs_sles_petsc_t *)
337 * log_type <-- log type
338 *----------------------------------------------------------------------------*/
339
340 void
341 cs_sles_petsc_log(const void *context,
342 cs_log_t log_type);
343
344 /*----------------------------------------------------------------------------*/
345 /*!
346 * \brief Print information on PETSc library.
347 *
348 * \param[in] log_type log type
349 */
350 /*----------------------------------------------------------------------------*/
351
352 void
353 cs_sles_petsc_library_info(cs_log_t log_type);
354
355 /*----------------------------------------------------------------------------*/
356
357 END_C_DECLS
358
359 #endif /* __CS_SLES_PETSC_H__ */
360