8.3
general documentation
cs_multigrid.h
Go to the documentation of this file.
1#ifndef __CS_MULTIGRID_H__
2#define __CS_MULTIGRID_H__
3
4/*============================================================================
5 * Multigrid solver.
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_grid.h"
36#include "cs_sles.h"
37#include "cs_sles_it.h"
38#include "cs_sles_pc.h"
39#include "cs_time_plot.h"
40
41/*----------------------------------------------------------------------------*/
42
44
45/*============================================================================
46 * Macro definitions
47 *============================================================================*/
48
49/*============================================================================
50 * Type definitions
51 *============================================================================*/
52
53/*----------------------------------------------------------------------------
54 * Multigrid types
55 *----------------------------------------------------------------------------*/
56
57typedef enum {
58
65
66/* Multigrid linear solver context (opaque) */
67
68typedef struct _cs_multigrid_t cs_multigrid_t;
69
70/*============================================================================
71 * Global variables
72 *============================================================================*/
73
74/* Names for multigrid types */
75
76extern const char *cs_multigrid_type_name[];
77
78/*=============================================================================
79 * Public function prototypes
80 *============================================================================*/
81
82/*----------------------------------------------------------------------------
83 * Initialize multigrid solver API.
84 *----------------------------------------------------------------------------*/
85
86void
88
89/*----------------------------------------------------------------------------
90 * Finalize multigrid solver API.
91 *----------------------------------------------------------------------------*/
92
93void
95
96/*----------------------------------------------------------------------------*
97 * Set maximum grid level which should run on device (i.e. GPU).
98 *----------------------------------------------------------------------------*/
99
100void
102
103/*----------------------------------------------------------------------------
104 * Define and associate a multigrid sparse linear system solver
105 * for a given field or equation name.
106 *
107 * If this system did not previously exist, it is added to the list of
108 * "known" systems. Otherwise, its definition is replaced by the one
109 * defined here.
110 *
111 * This is a utility function: if finer control is needed, see
112 * cs_sles_define() and cs_multigrid_create().
113 *
114 * Note that this function returns a pointer directly to the multigrid solver
115 * management structure. This may be used to set further options, for
116 * example calling cs_multigrid_set_coarsening_options() and
117 * cs_multigrid_set_solver_options().
118 * If needed, cs_sles_find() may be used to obtain a pointer to the
119 * matching cs_sles_t container.
120 *
121 * \param[in] f_id associated field id, or < 0
122 * \param[in] name associated name if f_id < 0, or NULL
123 * \param[in] mg_type type of multigrid algorithm to use
124 *
125 * \return pointer to new multigrid info and context
126 */
127/*----------------------------------------------------------------------------*/
128
130cs_multigrid_define(int f_id,
131 const char *name,
132 cs_multigrid_type_t mg_type);
133
134/*----------------------------------------------------------------------------*/
144/*----------------------------------------------------------------------------*/
145
148
149/*----------------------------------------------------------------------------
150 * Destroy multigrid linear system solver info and context.
151 *
152 * parameters:
153 * context <-> pointer to multigrid linear solver info
154 * (actual type: cs_multigrid_t **)
155 *----------------------------------------------------------------------------*/
156
157void
158cs_multigrid_destroy(void **context);
159
160/*----------------------------------------------------------------------------
161 * Create multigrid sparse linear system solver info and context
162 * based on existing info and context.
163 *
164 * parameters:
165 * context <-- pointer to reference info and context
166 * (actual type: cs_multigrid_t *)
167 *
168 * returns:
169 * pointer to newly created solver info object
170 * (actual type: cs_multigrid_t *)
171 *----------------------------------------------------------------------------*/
172
173void *
174cs_multigrid_copy(const void *context);
175
176/*----------------------------------------------------------------------------
177 * Set multigrid coarsening parameters.
178 *
179 * parameters:
180 * mg <-> pointer to multigrid info and context
181 * aggregation_limit <-- maximum allowed fine rows per coarse cell
182 * coarsening_type <-- coarsening type; see cs_grid_coarsening_t
183 * n_max_levels <-- maximum number of grid levels
184 * min_g_rows <-- global number of rows on coarse grids
185 * under which no coarsening occurs
186 * p0p1_relax <-- p0/p1 relaxation_parameter
187 * postprocess_block_size <-- if > 0, postprocess coarsening
188 * (using coarse cell numbers modulo this value)
189 *----------------------------------------------------------------------------*/
190
191void
193 (cs_multigrid_t *mg,
194 int aggregation_limit,
195 cs_grid_coarsening_t coarsening_type,
196 int n_max_levels,
197 cs_gnum_t min_g_rows,
198 double p0p1_relax,
199 int postprocess_block_size);
200
201/*----------------------------------------------------------------------------*/
202/*
203 * \brief Set specific multigrid coarsening parameters for fine grids.
204 *
205 * \param[in, out] mg pointer to multigrid info and context
206 * \param[in] f_settings_threshold grids of this level or higher use
207 * standard (coarse grid) settings.
208 * \param[in] aggregation_limit maximum allowed fine rows
209 * per coarse row
210 * \param[in] coarsening_type coarsening type;
211 * see \ref cs_grid_coarsening_t
212 */
213/*----------------------------------------------------------------------------*/
214
215void
217 (cs_multigrid_t *mg,
218 int f_settings_threshold,
219 int aggregation_limit,
220 cs_grid_coarsening_t coarsening_type);
221
222/*----------------------------------------------------------------------------
223 * Set multigrid parameters for associated iterative solvers.
224 *
225 * On a GPU, some parameters may be replaced by the closest GPU equivalents.
226 * For finer control, use \ref cs_multigrid_set_coarsening_options_d after
227 * calling this function to modify parameters for solvers running on device.
228 *
229 * parameters:
230 * mg <-> pointer to multigrid info and context
231 * descent_smoother_type <-- type of smoother for descent
232 * ascent_smoother_type <-- type of smoother for ascent
233 * coarse_solver_type <-- type of solver
234 * n_max_cycles <-- maximum number of cycles
235 * n_max_iter_descent <-- maximum iterations per descent phase
236 * n_max_iter_ascent <-- maximum iterations per descent phase
237 * n_max_iter_coarse <-- maximum iterations per coarsest solution
238 * poly_degree_descent <-- preconditioning polynomial degree
239 * for descent phases (0: diagonal)
240 * poly_degree_ascent <-- preconditioning polynomial degree
241 * for ascent phases (0: diagonal)
242 * poly_degree_coarse <-- preconditioning polynomial degree
243 * for coarse solver (0: diagonal)
244 * precision_mult_descent <-- precision multiplier for descent phases
245 * (levels >= 1)
246 * precision_mult_ascent <-- precision multiplier for ascent phases
247 * precision_mult_coarse <-- precision multiplier for coarsest grid
248 *----------------------------------------------------------------------------*/
249
250void
252 cs_sles_it_type_t descent_smoother_type,
253 cs_sles_it_type_t ascent_smoother_type,
254 cs_sles_it_type_t coarse_solver_type,
255 int n_max_cycles,
256 int n_max_iter_descent,
257 int n_max_iter_ascent,
258 int n_max_iter_coarse,
259 int poly_degree_descent,
260 int poly_degree_ascent,
261 int poly_degree_coarse,
262 double precision_mult_descent,
263 double precision_mult_ascent,
264 double precision_mult_coarse);
265
266/*----------------------------------------------------------------------------*/
267/*
268 * \brief Set multigrid device solver parameters for associated
269 * iterative solvers.
270 *
271 * \param[in, out] mg pointer to multigrid info
272 * and context
273 * \param[in] descent_smoother_type type of smoother for descent
274 * \param[in] ascent_smoother_type type of smoother for ascent
275 * \param[in] coarse_solver_type type of solver for coarsest grid
276 * \param[in] n_max_iter_descent maximum iterations
277 * per descent smoothing
278 * \param[in] n_max_iter_ascent maximum iterations
279 * per ascent smoothing
280 * \param[in] n_max_iter_coarse maximum iterations
281 * per coarsest solution
282 * \param[in] poly_degree_descent preconditioning polynomial degree
283 * for descent phases (0: diagonal)
284 * \param[in] poly_degree_ascent preconditioning polynomial degree
285 * for ascent phases (0: diagonal)
286 * \param[in] poly_degree_coarse preconditioning polynomial degree
287 * for coarse solver (0: diagonal)
288 */
289/*----------------------------------------------------------------------------*/
290
291void
293 cs_sles_it_type_t descent_smoother_type,
294 cs_sles_it_type_t ascent_smoother_type,
295 cs_sles_it_type_t coarse_solver_type,
296 int n_max_iter_descent,
297 int n_max_iter_ascent,
298 int n_max_iter_coarse,
299 int poly_degree_descent,
300 int poly_degree_ascent,
301 int poly_degree_coarse);
302
303/*----------------------------------------------------------------------------*/
304/*
305 * \brief Set the max. number of cycles for a multigrid
306 *
307 * \param[in, out] mg pointer to multigrid info and context
308 * \param[in] n_max_cycles maximum number of cycles
309 */
310/*----------------------------------------------------------------------------*/
311
312void
314 int n_max_cycles);
315
316/*----------------------------------------------------------------------------*/
317/*
318 * \brief Indicate if a multigrid solver requires an MSR matrix input.
319 *
320 * \param[in] mg pointer to multigrid info and context
321 *
322 * \return true if MSR is needed, false otherwise.
323 */
324/*----------------------------------------------------------------------------*/
325
326bool
328
329/*----------------------------------------------------------------------------
330 * Setup multigrid sparse linear equation solver.
331 *
332 * parameters:
333 * context <-> pointer to multigrid info and context
334 * (actual type: cs_multigrid_t *)
335 * name <-- pointer to name of linear system
336 * a <-- associated matrix
337 * verbosity <-- associated verbosity
338 *----------------------------------------------------------------------------*/
339
340void
341cs_multigrid_setup(void *context,
342 const char *name,
343 const cs_matrix_t *a,
344 int verbosity);
345
346/*----------------------------------------------------------------------------
347 * Setup multigrid sparse linear equation solver with separate
348 * convection-diffusion matrixes
349 *
350 * parameters:
351 * context <-> pointer to multigrid info and context
352 * (actual type: cs_multigrid_t *)
353 * name <-- pointer to name of linear system
354 * a <-- associated matrix
355 * conv_diff <-- convection-diffusion mode
356 * verbosity <-- associated verbosity
357 *----------------------------------------------------------------------------*/
358
359void
360cs_multigrid_setup_conv_diff(void *context,
361 const char *name,
362 const cs_matrix_t *a,
363 bool conv_diff,
364 int verbosity);
365
366/*----------------------------------------------------------------------------
367 * Call multigrid sparse linear equation solver.
368 *
369 * parameters:
370 * context <-> pointer to iterative sparse linear solver info
371 * (actual type: cs_multigrid_t *)
372 * name <-- pointer to name of linear system
373 * a <-- matrix
374 * verbosity <-- associated verbosity
375 * precision <-- solver precision
376 * r_norm <-- residual normalization
377 * n_iter --> number of iterations
378 * residual --> residual
379 * rhs <-- right hand side
380 * vx_ini <-- initial solution (vx if nonzero, nullptr if zero)
381 * vx <-> system solution
382 * aux_size <-- number of elements in aux_vectors
383 * aux_vectors --- optional working area (internal allocation if NULL)
384 *
385 * returns:
386 * convergence state
387 *----------------------------------------------------------------------------*/
388
390cs_multigrid_solve(void *context,
391 const char *name,
392 const cs_matrix_t *a,
393 int verbosity,
394 double precision,
395 double r_norm,
396 int *n_iter,
397 double *residual,
398 const cs_real_t *rhs,
399 cs_real_t *vx_ini,
400 cs_real_t *vx,
401 size_t aux_size,
402 void *aux_vectors);
403
404/*----------------------------------------------------------------------------
405 * Free iterative sparse linear equation solver setup context.
406 *
407 * Note that this function should free resolution-related data, such as
408 * buffers and preconditioning but does not free the whole context,
409 * as info used for logging (especially performance data) is maintained.
410
411 * parameters:
412 * context <-> pointer to iterative sparse linear solver info
413 * (actual type: cs_multigrid_t *)
414 *----------------------------------------------------------------------------*/
415
416void
417cs_multigrid_free(void *context);
418
419/*----------------------------------------------------------------------------
420 * Log sparse linear equation solver info.
421 *
422 * parameters:
423 * context <-> pointer to iterative sparse linear solver info
424 * (actual type: cs_multigrid_t *)
425 * log_type <-- log type
426 *----------------------------------------------------------------------------*/
427
428void
429cs_multigrid_log(const void *context,
430 cs_log_t log_type);
431
432/*----------------------------------------------------------------------------*/
440/*----------------------------------------------------------------------------*/
441
444
445/*----------------------------------------------------------------------------
446 * Error handler for multigrid sparse linear equation solver.
447 *
448 * In case of divergence or breakdown, this error handler outputs
449 * postprocessing data to assist debugging, then aborts the run.
450 * It does nothing in case the maximum iteration count is reached.
451 *
452 * parameters:
453 * sles <-> pointer to solver object
454 * state <-- convergence status
455 * a <-- matrix
456 * rhs <-- right hand side
457 * vx <-> system solution
458 *
459 * returns:
460 * false (do not attempt new solve)
461 *----------------------------------------------------------------------------*/
462
463bool
466 const cs_matrix_t *a,
467 const cs_real_t *rhs,
468 cs_real_t *vx);
469
470/*----------------------------------------------------------------------------
471 * Set plotting options for multigrid.
472 *
473 * parameters:
474 * mg <-> pointer to multigrid info and context
475 * base_name <-- base plot name to activate, NULL otherwise
476 * use_iteration <-- if true, use iteration as time stamp
477 * otherwise, use wall clock time
478 *----------------------------------------------------------------------------*/
479
480void
482 const char *base_name,
483 bool use_iteration);
484
485/*----------------------------------------------------------------------------*/
497/*----------------------------------------------------------------------------*/
498
499void
501 int *rank_stride,
502 int *rows_mean_threshold,
503 cs_gnum_t *rows_glob_threshold);
504
505/*----------------------------------------------------------------------------*/
517/*----------------------------------------------------------------------------*/
518
519void
521 int rank_stride,
522 int rows_mean_threshold,
523 cs_gnum_t rows_glob_threshold);
524
525/*----------------------------------------------------------------------------*/
538/*----------------------------------------------------------------------------*/
539
540const cs_grid_t *
542 int level);
543
544/*----------------------------------------------------------------------------*/
545
547
548#endif /* __CS_MULTIGRID_H__ */
#define BEGIN_C_DECLS
Definition: cs_defs.h:542
double cs_real_t
Floating-point value.
Definition: cs_defs.h:342
uint64_t cs_gnum_t
global mesh entity number
Definition: cs_defs.h:325
#define END_C_DECLS
Definition: cs_defs.h:543
cs_grid_coarsening_t
Definition: cs_grid.h:55
struct _cs_grid_t cs_grid_t
Definition: cs_grid.h:68
cs_log_t
Definition: cs_log.h:48
struct _cs_matrix_t cs_matrix_t
Definition: cs_matrix.h:110
void cs_multigrid_setup(void *context, const char *name, const cs_matrix_t *a, int verbosity)
Setup multigrid sparse linear equation solver.
Definition: cs_multigrid.cpp:5056
void cs_multigrid_set_coarsening_options_fine_grid(cs_multigrid_t *mg, int f_settings_threshold, int aggregation_limit, cs_grid_coarsening_t coarsening_type)
Set specific multigrid coarsening parameters for fine grids.
Definition: cs_multigrid.cpp:4778
void cs_multigrid_setup_conv_diff(void *context, const char *name, const cs_matrix_t *a, bool conv_diff, int verbosity)
Setup multigrid sparse linear equation solver.
Definition: cs_multigrid.cpp:5082
void cs_multigrid_set_solver_options_d(cs_multigrid_t *mg, cs_sles_it_type_t descent_smoother_type, cs_sles_it_type_t ascent_smoother_type, cs_sles_it_type_t coarse_solver_type, int n_max_iter_descent, int n_max_iter_ascent, int n_max_iter_coarse, int poly_degree_descent, int poly_degree_ascent, int poly_degree_coarse)
Set multigrid device solver parameters for associated iterative solvers.
Definition: cs_multigrid.cpp:4947
struct _cs_multigrid_t cs_multigrid_t
Definition: cs_multigrid.h:68
cs_multigrid_type_t
Definition: cs_multigrid.h:57
@ CS_MULTIGRID_K_CYCLE
Definition: cs_multigrid.h:60
@ CS_MULTIGRID_N_TYPES
Definition: cs_multigrid.h:62
@ CS_MULTIGRID_K_CYCLE_HPC
Definition: cs_multigrid.h:61
@ CS_MULTIGRID_V_CYCLE
Definition: cs_multigrid.h:59
void cs_multigrid_get_merge_options(const cs_multigrid_t *mg, int *rank_stride, int *rows_mean_threshold, cs_gnum_t *rows_glob_threshold)
Query the global multigrid parameters for parallel grid merging.
Definition: cs_multigrid.cpp:5764
void cs_multigrid_set_max_grid_level_for_device(int level)
Set maximum grid level which should run on device (i.e. GPU).
Definition: cs_multigrid.cpp:4431
cs_multigrid_t * cs_multigrid_define(int f_id, const char *name, cs_multigrid_type_t mg_type)
Define and associate a multigrid sparse linear system solver for a given field or equation name.
Definition: cs_multigrid.cpp:4464
void * cs_multigrid_copy(const void *context)
Create multigrid sparse linear system solver info and context based on existing info and context.
Definition: cs_multigrid.cpp:4676
void cs_multigrid_set_coarsening_options(cs_multigrid_t *mg, int aggregation_limit, cs_grid_coarsening_t coarsening_type, int n_max_levels, cs_gnum_t min_g_rows, double p0p1_relax, int postprocess_block_size)
Set multigrid coarsening parameters.
Definition: cs_multigrid.cpp:4741
cs_sles_pc_t * cs_multigrid_pc_create(cs_multigrid_type_t mg_type)
Create a multigrid preconditioner.
Definition: cs_multigrid.cpp:5476
cs_sles_convergence_state_t cs_multigrid_solve(void *context, const char *name, const cs_matrix_t *a, int verbosity, double precision, double r_norm, int *n_iter, double *residual, const cs_real_t *rhs, cs_real_t *vx_ini, cs_real_t *vx, size_t aux_size, void *aux_vectors)
Call multigrid sparse linear equation solver.
Definition: cs_multigrid.cpp:5167
void cs_multigrid_set_merge_options(cs_multigrid_t *mg, int rank_stride, int rows_mean_threshold, cs_gnum_t rows_glob_threshold)
Set global multigrid parameters for parallel grid merging behavior.
Definition: cs_multigrid.cpp:5801
void cs_multigrid_set_plot_options(cs_multigrid_t *mg, const char *base_name, bool use_iteration)
Set plotting options for multigrid.
Definition: cs_multigrid.cpp:5718
cs_multigrid_t * cs_multigrid_create(cs_multigrid_type_t mg_type)
Create multigrid linear system solver info and context.
Definition: cs_multigrid.cpp:4501
bool cs_multigrid_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)
const cs_grid_t * cs_multigrid_get_grid(const cs_multigrid_t *mg, int level)
Return a pointer to a grid associated with a given multigrid setup and level.
Definition: cs_multigrid.cpp:5831
void cs_multigrid_set_solver_options(cs_multigrid_t *mg, cs_sles_it_type_t descent_smoother_type, cs_sles_it_type_t ascent_smoother_type, cs_sles_it_type_t coarse_solver_type, int n_max_cycles, int n_max_iter_descent, int n_max_iter_ascent, int n_max_iter_coarse, int poly_degree_descent, int poly_degree_ascent, int poly_degree_coarse, double precision_mult_descent, double precision_mult_ascent, double precision_mult_coarse)
Set multigrid parameters for associated iterative solvers.
Definition: cs_multigrid.cpp:4827
void cs_multigrid_free(void *context)
Free multigrid sparse linear equation solver setup context.
Definition: cs_multigrid.cpp:5406
bool cs_multigrid_need_msr(const cs_multigrid_t *mg)
Indicate if a multigrid solver requires an MSR matrix input.
Definition: cs_multigrid.cpp:5020
void cs_multigrid_initialize(void)
Initialize multigrid solver API.
Definition: cs_multigrid.cpp:4408
void cs_multigrid_log(const void *context, cs_log_t log_type)
Log multigrid solver info.
Definition: cs_multigrid.cpp:4708
void cs_multigrid_set_max_cycles(cs_multigrid_t *mg, int n_max_cycles)
Set the max. number of cycles for a multigrid.
Definition: cs_multigrid.cpp:4998
void cs_multigrid_finalize(void)
Finalize multigrid solver API.
Definition: cs_multigrid.cpp:4419
const char * cs_multigrid_type_name[]
void cs_multigrid_destroy(void **context)
Destroy multigrid linear system solver info and context.
Definition: cs_multigrid.cpp:4623
cs_sles_convergence_state_t
Definition: cs_sles.h:56
struct _cs_sles_t cs_sles_t
Definition: cs_sles.h:68
cs_sles_it_type_t
Definition: cs_sles_it.h:55
struct _cs_sles_pc_t cs_sles_pc_t
Definition: cs_sles_pc.h:66