programmer's 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-2015 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_sles.h"
36 #include "cs_sles_it.h"
37 #include "cs_time_plot.h"
38 
39 /*----------------------------------------------------------------------------*/
40 
42 
43 /*============================================================================
44  * Macro definitions
45  *============================================================================*/
46 
47 /*============================================================================
48  * Type definitions
49  *============================================================================*/
50 
51 /* Multigrid linear solver context (opaque) */
52 
53 typedef struct _cs_multigrid_t cs_multigrid_t;
54 
55 /*============================================================================
56  * Global variables
57  *============================================================================*/
58 
59 /*=============================================================================
60  * Public function prototypes
61  *============================================================================*/
62 
63 /*----------------------------------------------------------------------------
64  * Initialize multigrid solver API.
65  *----------------------------------------------------------------------------*/
66 
67 void
69 
70 /*----------------------------------------------------------------------------
71  * Finalize multigrid solver API.
72  *----------------------------------------------------------------------------*/
73 
74 void
76 
77 /*----------------------------------------------------------------------------
78  * Indicate if multigrid solver API is used for at least one system.
79  *
80  * returns:
81  * true if at least one system uses a multigrid solver, false otherwise
82  *----------------------------------------------------------------------------*/
83 
84 bool
86 
87 /*----------------------------------------------------------------------------
88  * Define and associate a multigrid sparse linear system solver
89  * for a given field or equation name.
90  *
91  * If this system did not previously exist, it is added to the list of
92  * "known" systems. Otherwise, its definition is replaced by the one
93  * defined here.
94  *
95  * This is a utility function: if finer control is needed, see
96  * cs_sles_define() and cs_multigrid_create().
97  *
98  * Note that this function returns a pointer directly to the multigrid solver
99  * management structure. This may be used to set further options, for
100  * example calling cs_multigrid_set_coarsening_options() and
101  * cs_multigrid_set_solver_options().
102  * If needed, cs_sles_find() may be used to obtain a pointer to the
103  * matching cs_sles_t container.
104  *
105  * parameters:
106  * f_id <-- associated field id, or < 0
107  * name <-- associated name if f_id < 0, or NULL
108  *
109  * \return pointer to new multigrid info and context
110  */
111 /*----------------------------------------------------------------------------*/
112 
114 cs_multigrid_define(int f_id,
115  const char *name);
116 
117 /*----------------------------------------------------------------------------
118  * Create multigrid linear system solver info and context.
119  *
120  * The multigrid variant is an ACM (Additive Corrective Multigrid) method.
121  *
122  * returns:
123  * pointer to new multigrid info and context
124  *----------------------------------------------------------------------------*/
125 
127 cs_multigrid_create(void);
128 
129 /*----------------------------------------------------------------------------
130  * Destroy multigrid linear system solver info and context.
131  *
132  * parameters:
133  * context <-> pointer to multigrid linear solver info
134  * (actual type: cs_multigrid_t **)
135  *----------------------------------------------------------------------------*/
136 
137 void
138 cs_multigrid_destroy(void **context);
139 
140 /*----------------------------------------------------------------------------
141  * Create multigrid sparse linear system solver info and context
142  * based on existing info and context.
143  *
144  * parameters:
145  * context <-- pointer to reference info and context
146  * (actual type: cs_multigrid_t *)
147  *
148  * returns:
149  * pointer to newly created solver info object
150  * (actual type: cs_multigrid_t *)
151  *----------------------------------------------------------------------------*/
152 
153 void *
154 cs_multigrid_copy(const void *context);
155 
156 /*----------------------------------------------------------------------------
157  * Set multigrid coarsening parameters.
158  *
159  * parameters:
160  * mg <-> pointer to multigrid info and context
161  * aggregation_limit <-- maximum allowed fine cells per coarse cell
162  * coarsening_type <-- coarsening type:
163  * 0: algebraic, natural face traversal;
164  * 1: algebraic, face traveral by criteria;
165  * 2: algebraic, Hilbert face traversal;
166  * n_max_levels <-- maximum number of grid levels
167  * min_g_cells <-- global number of cells on coarse grids
168  * under which no coarsening occurs
169  * p0p1_relax <-- p0/p1 relaxation_parameter
170  * postprocess_block_size <-- if > 0, postprocess coarsening
171  * (using coarse cell numbers modulo this value)
172  *----------------------------------------------------------------------------*/
173 
174 void
176  int aggregation_limit,
177  int coarsening_type,
178  int n_max_levels,
179  cs_gnum_t min_g_cells,
180  double p0p1_relax,
181  int postprocess_block_size);
182 
183 /*----------------------------------------------------------------------------
184  * Set multigrid parameters for associated iterative solvers.
185  *
186  * parameters:
187  * mg <-> pointer to multigrid info and context
188  * descent_smoother_type <-- type of smoother for descent
189  * ascent_smoother_type <-- type of smoother for ascent
190  * coarse_solver_type <-- type of solver
191  * n_max_cycles <-- maximum number of cycles
192  * n_max_iter_descent <-- maximum iterations per descent phase
193  * n_max_iter_ascent <-- maximum iterations per descent phase
194  * n_max_iter_coarse <-- maximum iterations per coarsest solution
195  * poly_degree_descent <-- preconditioning polynomial degree
196  * for descent phases (0: diagonal)
197  * poly_degree_ascent <-- preconditioning polynomial degree
198  * for ascent phases (0: diagonal)
199  * poly_degree_coarse <-- preconditioning polynomial degree
200  * for coarse solver (0: diagonal)
201  * precision_mult_descent <-- precision multiplier for descent phases
202  * (levels >= 1)
203  * precision_mult_ascent <-- precision multiplier for ascent phases
204  * precision_mult_coarse <-- precision multiplier for coarsest grid
205  *----------------------------------------------------------------------------*/
206 
207 void
209  cs_sles_it_type_t descent_smoother_type,
210  cs_sles_it_type_t ascent_smoother_type,
211  cs_sles_it_type_t coarse_solver_type,
212  int n_max_cycles,
213  int n_max_iter_descent,
214  int n_max_iter_ascent,
215  int n_max_iter_coarse,
216  int poly_degree_descent,
217  int poly_degree_ascent,
218  int poly_degree_coarse,
219  double precision_mult_descent,
220  double precision_mult_ascent,
221  double precision_mult_coarse);
222 
223 /*----------------------------------------------------------------------------
224  * Setup multigrid sparse linear equation solver.
225  *
226  * parameters:
227  * context <-> pointer to multigrid info and context
228  * (actual type: cs_multigrid_t *)
229  * name <-- pointer to name of linear system
230  * a <-- associated matrix
231  * verbosity <-- associated verbosity
232  *----------------------------------------------------------------------------*/
233 
234 void
235 cs_multigrid_setup(void *context,
236  const char *name,
237  const cs_matrix_t *a,
238  int verbosity);
239 
240 /*----------------------------------------------------------------------------
241  * Call multigrid sparse linear equation solver.
242  *
243  * parameters:
244  * context <-> pointer to iterative sparse linear solver info
245  * (actual type: cs_multigrid_t *)
246  * name <-- pointer to name of linear system
247  * a <-- matrix
248  * verbosity <-- associated verbosity
249  * rotation_mode <-- halo update option for rotational periodicity
250  * precision <-- solver precision
251  * r_norm <-- residue normalization
252  * n_iter --> number of iterations
253  * residue --> residue
254  * rhs <-- right hand side
255  * vx <-> system solution
256  * aux_size <-- number of elements in aux_vectors
257  * aux_vectors --- optional working area (internal allocation if NULL)
258  *
259  * returns:
260  * convergence state
261  *----------------------------------------------------------------------------*/
262 
264 cs_multigrid_solve(void *context,
265  const char *name,
266  const cs_matrix_t *a,
267  int verbosity,
268  cs_halo_rotation_t rotation_mode,
269  double precision,
270  double r_norm,
271  int *n_iter,
272  double *residue,
273  const cs_real_t *rhs,
274  cs_real_t *vx,
275  size_t aux_size,
276  void *aux_vectors);
277 
278 /*----------------------------------------------------------------------------
279  * Free iterative sparse linear equation solver setup context.
280  *
281  * Note that this function should free resolution-related data, such as
282  * buffers and preconditioning but doesd not free the whole context,
283  * as info used for logging (especially performance data) is maintained.
284 
285  * parameters:
286  * context <-> pointer to iterative sparse linear solver info
287  * (actual type: cs_multigrid_t *)
288  *----------------------------------------------------------------------------*/
289 
290 void
291 cs_multigrid_free(void *context);
292 
293 /*----------------------------------------------------------------------------
294  * Log sparse linear equation solver info.
295  *
296  * parameters:
297  * context <-> pointer to iterative sparse linear solver info
298  * (actual type: cs_multigrid_t *)
299  * log_type <-- log type
300  *----------------------------------------------------------------------------*/
301 
302 void
303 cs_multigrid_log(const void *context,
304  cs_log_t log_type);
305 
306 /*----------------------------------------------------------------------------
307  * Error handler for multigrid sparse linear equation solver.
308  *
309  * In case of divergence or breakdown, this error handler outputs
310  * postprocessing data to assist debugging, then aborts the run.
311  * It does nothing in case the maximum iteration count is reached.
312  *
313  * parameters:
314  * context <-> pointer to multigrid sparse linear system solver info
315  * (actual type: cs_multigrid_t *)
316  * name <-- pointer to name of linear system
317  * state <-- convergence status
318  * a <-- matrix
319  * rotation_mode <-- halo update option for rotational periodicity
320  * rhs <-- right hand side
321  * vx <-> system solution
322  *----------------------------------------------------------------------------*/
323 
324 void
327  const char *name,
328  const cs_matrix_t *a,
329  cs_halo_rotation_t rotation_mode,
330  const cs_real_t *rhs,
331  cs_real_t *vx);
332 
333 /*----------------------------------------------------------------------------
334  * Set plotting options for multigrid.
335  *
336  * parameters:
337  * mg <-> pointer to multigrid info and context
338  * base_name <-- base plot name to activate, NULL otherwise
339  * use_iteration <-- if true, use iteration as time stamp
340  * otherwise, use wall clock time
341  *----------------------------------------------------------------------------*/
342 
343 void
345  const char *base_name,
346  bool use_iteration);
347 
348 /*----------------------------------------------------------------------------*/
349 
351 
352 #endif /* __CS_MULTIGRID_H__ */
unsigned long cs_gnum_t
global mesh entity number
Definition: cs_defs.h:280
cs_halo_rotation_t
Definition: cs_halo.h:59
cs_sles_it_type_t
Iterative solver types.
Definition: cs_sles_it.h:56
cs_multigrid_t * cs_multigrid_define(int f_id, const char *name)
Define and associate a multigrid sparse linear system solver for a given field or equation name...
Definition: cs_multigrid.c:1926
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:2069
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:2893
#define BEGIN_C_DECLS
Definition: cs_defs.h:419
void cs_multigrid_destroy(void **context)
Destroy multigrid linear system solver info and context.
Definition: cs_multigrid.c:2013
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:2236
bool cs_multigrid_needed(void)
Indicate if multigrid solver API is used for at least one system.
Definition: cs_multigrid.c:1894
void cs_multigrid_log(const void *context, cs_log_t log_type)
Log multigrid solver info.
Definition: cs_multigrid.c:2097
struct _cs_matrix_t cs_matrix_t
Definition: cs_matrix.h:86
void cs_multigrid_error_post_and_abort(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)
cs_multigrid_t * cs_multigrid_create(void)
Create multigrid linear system solver info and context.
Definition: cs_multigrid.c:1960
cs_sles_convergence_state_t
Convergence status indicator.
Definition: cs_sles.h:55
void cs_multigrid_finalize(void)
Finalize multigrid solver API.
Definition: cs_multigrid.c:1878
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:2516
struct _cs_multigrid_t cs_multigrid_t
Definition: cs_multigrid.h:53
void cs_multigrid_free(void *context)
Free multigrid sparse linear equation solver setup context.
Definition: cs_multigrid.c:2655
cs_log_t
Definition: cs_log.h:47
void cs_multigrid_initialize(void)
Initialize multigrid solver API.
Definition: cs_multigrid.c:1867
#define END_C_DECLS
Definition: cs_defs.h:420
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_cells, double p0p1_relax, int postprocess_block_size)
Set multigrid coarsening parameters.
Definition: cs_multigrid.c:2131
double cs_real_t
Definition: cs_defs.h:296
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:2184