6.2
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-2020 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 
57 typedef enum {
58 
65 
66 /* Multigrid linear solver context (opaque) */
67 
68 typedef struct _cs_multigrid_t cs_multigrid_t;
69 
70 /*============================================================================
71  * Global variables
72  *============================================================================*/
73 
74 /* Names for multigrid types */
75 
76 extern const char *cs_multigrid_type_name[];
77 
78 /*=============================================================================
79  * Public function prototypes
80  *============================================================================*/
81 
82 /*----------------------------------------------------------------------------
83  * Initialize multigrid solver API.
84  *----------------------------------------------------------------------------*/
85 
86 void
88 
89 /*----------------------------------------------------------------------------
90  * Finalize multigrid solver API.
91  *----------------------------------------------------------------------------*/
92 
93 void
95 
96 /*----------------------------------------------------------------------------
97  * Define and associate a multigrid sparse linear system solver
98  * for a given field or equation name.
99  *
100  * If this system did not previously exist, it is added to the list of
101  * "known" systems. Otherwise, its definition is replaced by the one
102  * defined here.
103  *
104  * This is a utility function: if finer control is needed, see
105  * cs_sles_define() and cs_multigrid_create().
106  *
107  * Note that this function returns a pointer directly to the multigrid solver
108  * management structure. This may be used to set further options, for
109  * example calling cs_multigrid_set_coarsening_options() and
110  * cs_multigrid_set_solver_options().
111  * If needed, cs_sles_find() may be used to obtain a pointer to the
112  * matching cs_sles_t container.
113  *
114  * \param[in] f_id associated field id, or < 0
115  * \param[in] name associated name if f_id < 0, or NULL
116  * \param[in] mg_type type of multigrid algorithm to use
117  *
118  * \return pointer to new multigrid info and context
119  */
120 /*----------------------------------------------------------------------------*/
121 
123 cs_multigrid_define(int f_id,
124  const char *name,
125  cs_multigrid_type_t mg_type);
126 
127 /*----------------------------------------------------------------------------*/
137 /*----------------------------------------------------------------------------*/
138 
141 
142 /*----------------------------------------------------------------------------
143  * Destroy multigrid linear system solver info and context.
144  *
145  * parameters:
146  * context <-> pointer to multigrid linear solver info
147  * (actual type: cs_multigrid_t **)
148  *----------------------------------------------------------------------------*/
149 
150 void
151 cs_multigrid_destroy(void **context);
152 
153 /*----------------------------------------------------------------------------
154  * Create multigrid sparse linear system solver info and context
155  * based on existing info and context.
156  *
157  * parameters:
158  * context <-- pointer to reference info and context
159  * (actual type: cs_multigrid_t *)
160  *
161  * returns:
162  * pointer to newly created solver info object
163  * (actual type: cs_multigrid_t *)
164  *----------------------------------------------------------------------------*/
165 
166 void *
167 cs_multigrid_copy(const void *context);
168 
169 /*----------------------------------------------------------------------------
170  * Set multigrid coarsening parameters.
171  *
172  * parameters:
173  * mg <-> pointer to multigrid info and context
174  * aggregation_limit <-- maximum allowed fine rows per coarse cell
175  * coarsening_type <-- coarsening type:
176  * 0: algebraic, natural face traversal;
177  * 1: algebraic, face traveral by criteria;
178  * 2: algebraic, Hilbert face traversal;
179  * n_max_levels <-- maximum number of grid levels
180  * min_g_rows <-- global number of rows on coarse grids
181  * under which no coarsening occurs
182  * p0p1_relax <-- p0/p1 relaxation_parameter
183  * postprocess_block_size <-- if > 0, postprocess coarsening
184  * (using coarse cell numbers modulo this value)
185  *----------------------------------------------------------------------------*/
186 
187 void
189  int aggregation_limit,
190  int coarsening_type,
191  int n_max_levels,
192  cs_gnum_t min_g_rows,
193  double p0p1_relax,
194  int postprocess_block_size);
195 
196 /*----------------------------------------------------------------------------
197  * Set multigrid parameters for associated iterative solvers.
198  *
199  * parameters:
200  * mg <-> pointer to multigrid info and context
201  * descent_smoother_type <-- type of smoother for descent
202  * ascent_smoother_type <-- type of smoother for ascent
203  * coarse_solver_type <-- type of solver
204  * n_max_cycles <-- maximum number of cycles
205  * n_max_iter_descent <-- maximum iterations per descent phase
206  * n_max_iter_ascent <-- maximum iterations per descent phase
207  * n_max_iter_coarse <-- maximum iterations per coarsest solution
208  * poly_degree_descent <-- preconditioning polynomial degree
209  * for descent phases (0: diagonal)
210  * poly_degree_ascent <-- preconditioning polynomial degree
211  * for ascent phases (0: diagonal)
212  * poly_degree_coarse <-- preconditioning polynomial degree
213  * for coarse solver (0: diagonal)
214  * precision_mult_descent <-- precision multiplier for descent phases
215  * (levels >= 1)
216  * precision_mult_ascent <-- precision multiplier for ascent phases
217  * precision_mult_coarse <-- precision multiplier for coarsest grid
218  *----------------------------------------------------------------------------*/
219 
220 void
222  cs_sles_it_type_t descent_smoother_type,
223  cs_sles_it_type_t ascent_smoother_type,
224  cs_sles_it_type_t coarse_solver_type,
225  int n_max_cycles,
226  int n_max_iter_descent,
227  int n_max_iter_ascent,
228  int n_max_iter_coarse,
229  int poly_degree_descent,
230  int poly_degree_ascent,
231  int poly_degree_coarse,
232  double precision_mult_descent,
233  double precision_mult_ascent,
234  double precision_mult_coarse);
235 
236 /*----------------------------------------------------------------------------
237  * Return solver type used on fine mesh.
238  *
239  * parameters:
240  * mg <-- pointer to multigrid info and context
241  *
242  * returns:
243  * type of smoother for descent (used for fine mesh)
244  *----------------------------------------------------------------------------*/
245 
248 
249 /*----------------------------------------------------------------------------
250  * Setup multigrid sparse linear equation solver.
251  *
252  * parameters:
253  * context <-> pointer to multigrid info and context
254  * (actual type: cs_multigrid_t *)
255  * name <-- pointer to name of linear system
256  * a <-- associated matrix
257  * verbosity <-- associated verbosity
258  *----------------------------------------------------------------------------*/
259 
260 void
261 cs_multigrid_setup(void *context,
262  const char *name,
263  const cs_matrix_t *a,
264  int verbosity);
265 
266 /*----------------------------------------------------------------------------
267  * Setup multigrid sparse linear equation solver with separate
268  * convection-diffusion matrixes
269  *
270  * parameters:
271  * context <-> pointer to multigrid info and context
272  * (actual type: cs_multigrid_t *)
273  * name <-- pointer to name of linear system
274  * a <-- associated matrix
275  * a_conv <-- associated matrix (convection)
276  * a_diff <-- associated matrix (diffusion)
277  * verbosity <-- associated verbosity
278  *----------------------------------------------------------------------------*/
279 
280 void
281 cs_multigrid_setup_conv_diff(void *context,
282  const char *name,
283  const cs_matrix_t *a,
284  const cs_matrix_t *a_conv,
285  const cs_matrix_t *a_diff,
286  int verbosity);
287 
288 /*----------------------------------------------------------------------------
289  * Call multigrid sparse linear equation solver.
290  *
291  * parameters:
292  * context <-> pointer to iterative sparse linear solver info
293  * (actual type: cs_multigrid_t *)
294  * name <-- pointer to name of linear system
295  * a <-- matrix
296  * verbosity <-- associated verbosity
297  * rotation_mode <-- halo update option for rotational periodicity
298  * precision <-- solver precision
299  * r_norm <-- residue normalization
300  * n_iter --> number of iterations
301  * residue --> residue
302  * rhs <-- right hand side
303  * vx <-> system solution
304  * aux_size <-- number of elements in aux_vectors
305  * aux_vectors --- optional working area (internal allocation if NULL)
306  *
307  * returns:
308  * convergence state
309  *----------------------------------------------------------------------------*/
310 
312 cs_multigrid_solve(void *context,
313  const char *name,
314  const cs_matrix_t *a,
315  int verbosity,
316  cs_halo_rotation_t rotation_mode,
317  double precision,
318  double r_norm,
319  int *n_iter,
320  double *residue,
321  const cs_real_t *rhs,
322  cs_real_t *vx,
323  size_t aux_size,
324  void *aux_vectors);
325 
326 /*----------------------------------------------------------------------------
327  * Free iterative sparse linear equation solver setup context.
328  *
329  * Note that this function should free resolution-related data, such as
330  * buffers and preconditioning but does not free the whole context,
331  * as info used for logging (especially performance data) is maintained.
332 
333  * parameters:
334  * context <-> pointer to iterative sparse linear solver info
335  * (actual type: cs_multigrid_t *)
336  *----------------------------------------------------------------------------*/
337 
338 void
339 cs_multigrid_free(void *context);
340 
341 /*----------------------------------------------------------------------------
342  * Log sparse linear equation solver info.
343  *
344  * parameters:
345  * context <-> pointer to iterative sparse linear solver info
346  * (actual type: cs_multigrid_t *)
347  * log_type <-- log type
348  *----------------------------------------------------------------------------*/
349 
350 void
351 cs_multigrid_log(const void *context,
352  cs_log_t log_type);
353 
354 /*----------------------------------------------------------------------------*/
362 /*----------------------------------------------------------------------------*/
363 
364 cs_sles_pc_t *
366 
367 /*----------------------------------------------------------------------------
368  * Error handler for multigrid sparse linear equation solver.
369  *
370  * In case of divergence or breakdown, this error handler outputs
371  * postprocessing data to assist debugging, then aborts the run.
372  * It does nothing in case the maximum iteration count is reached.
373  *
374  * parameters:
375  * sles <-> pointer to solver object
376  * state <-- convergence status
377  * a <-- matrix
378  * rotation_mode <-- halo update option for rotational periodicity
379  * rhs <-- right hand side
380  * vx <-> system solution
381  *
382  * returns:
383  * false (do not attempt new solve)
384  *----------------------------------------------------------------------------*/
385 
386 bool
389  const cs_matrix_t *a,
390  cs_halo_rotation_t rotation_mode,
391  const cs_real_t *rhs,
392  cs_real_t *vx);
393 
394 /*----------------------------------------------------------------------------
395  * Set plotting options for multigrid.
396  *
397  * parameters:
398  * mg <-> pointer to multigrid info and context
399  * base_name <-- base plot name to activate, NULL otherwise
400  * use_iteration <-- if true, use iteration as time stamp
401  * otherwise, use wall clock time
402  *----------------------------------------------------------------------------*/
403 
404 void
406  const char *base_name,
407  bool use_iteration);
408 
409 /*----------------------------------------------------------------------------*/
421 /*----------------------------------------------------------------------------*/
422 
423 void
425  int *rank_stride,
426  int *rows_mean_threshold,
427  cs_gnum_t *rows_glob_threshold);
428 
429 /*----------------------------------------------------------------------------*/
441 /*----------------------------------------------------------------------------*/
442 
443 void
445  int rank_stride,
446  int rows_mean_threshold,
447  cs_gnum_t rows_glob_threshold);
448 
449 /*----------------------------------------------------------------------------*/
450 
452 
453 #endif /* __CS_MULTIGRID_H__ */
Definition: cs_multigrid.h:61
unsigned long cs_gnum_t
global mesh entity number
Definition: cs_defs.h:286
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.c:4792
cs_halo_rotation_t
Definition: cs_halo.h:60
cs_sles_it_type_t
Definition: cs_sles_it.h:57
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.c:3949
Definition: cs_multigrid.h:62
void cs_multigrid_setup_conv_diff(void *context, const char *name, const cs_matrix_t *a, const cs_matrix_t *a_conv, const cs_matrix_t *a_diff, int verbosity)
Setup multigrid sparse linear equation solver.
Definition: cs_multigrid.c:4175
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.c:4746
cs_multigrid_t * cs_multigrid_create(cs_multigrid_type_t mg_type)
Create multigrid linear system solver info and context.
Definition: cs_multigrid.c:3781
Definition: cs_multigrid.h:59
#define BEGIN_C_DECLS
Definition: cs_defs.h:495
void cs_multigrid_destroy(void **context)
Destroy multigrid linear system solver info and context.
Definition: cs_multigrid.c:3896
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.c:4147
const char * cs_multigrid_type_name[]
struct _cs_sles_pc_t cs_sles_pc_t
Definition: cs_sles_pc.h:66
Definition: cs_multigrid.h:60
double cs_real_t
Floating-point value.
Definition: cs_defs.h:307
void cs_multigrid_log(const void *context, cs_log_t log_type)
Log multigrid solver info.
Definition: cs_multigrid.c:3977
struct _cs_matrix_t cs_matrix_t
Definition: cs_matrix.h:94
cs_sles_convergence_state_t
Convergence status indicator.
Definition: cs_sles.h:56
void cs_multigrid_set_coarsening_options(cs_multigrid_t *mg, int aggregation_limit, int 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.c:4009
struct _cs_sles_t cs_sles_t
Definition: cs_sles.h:68
void cs_multigrid_finalize(void)
Finalize multigrid solver API.
Definition: cs_multigrid.c:3711
cs_multigrid_type_t
Definition: cs_multigrid.h:57
cs_sles_pc_t * cs_multigrid_pc_create(cs_multigrid_type_t mg_type)
Create a multigrid preconditioner.
Definition: cs_multigrid.c:4500
double precision, save a
Definition: cs_fuel_incl.f90:146
cs_sles_convergence_state_t cs_multigrid_solve(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)
Call multigrid sparse linear equation solver.
Definition: cs_multigrid.c:4265
struct _cs_multigrid_t cs_multigrid_t
Definition: cs_multigrid.h:68
bool cs_multigrid_error_post_and_abort(cs_sles_t *sles, cs_sles_convergence_state_t state, const cs_matrix_t *a, cs_halo_rotation_t rotation_mode, const cs_real_t *rhs, cs_real_t *vx)
void cs_multigrid_free(void *context)
Free multigrid sparse linear equation solver setup context.
Definition: cs_multigrid.c:4432
cs_log_t
Definition: cs_log.h:48
cs_sles_it_type_t cs_multigrid_get_fine_solver_type(const cs_multigrid_t *mg)
Return solver type used on fine mesh.
Definition: cs_multigrid.c:4124
void cs_multigrid_initialize(void)
Initialize multigrid solver API.
Definition: cs_multigrid.c:3700
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.c:4829
#define END_C_DECLS
Definition: cs_defs.h:496
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.c:4062
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.c:3744