programmer's documentation
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
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-2014 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 
37 /*----------------------------------------------------------------------------*/
38 
40 
41 /*============================================================================
42  * Macro definitions
43  *============================================================================*/
44 
45 /*============================================================================
46  * Type definitions
47  *============================================================================*/
48 
49 /*============================================================================
50  * Global variables
51  *============================================================================*/
52 
53 /*============================================================================
54  * Public function prototypes for Fortran API
55  *============================================================================*/
56 
57 /*----------------------------------------------------------------------------
58  * Build a hierarchy of meshes starting from a fine mesh, for an
59  * ACM (Additive Corrective Multigrid) method.
60  *----------------------------------------------------------------------------*/
61 
62 void CS_PROCF(clmlga, CLMLGA)
63 (
64  const char *cname, /* <-- variable name */
65  const cs_int_t *lname, /* <-- variable name length */
66  const cs_int_t *isym, /* <-- Symmetry indicator:
67  1: symmetric; 2: not symmetric */
68  const cs_int_t *ibsize, /* <-- Matrix block size */
69  const cs_int_t *iesize, /* <-- Matrix extra diag block size */
70  const cs_int_t *nagmax, /* <-- Agglomeration count limit */
71  const cs_int_t *ncpost, /* <-- If > 0, postprocess coarsening, using
72  coarse cell numbers modulo ncpost */
73  const cs_int_t *iwarnp, /* <-- Verbosity level */
74  const cs_int_t *ngrmax, /* <-- Maximum number of grid levels */
75  const cs_int_t *ncegrm, /* <-- Maximum local number of cells on
76  coarsest grid */
77  const cs_real_t *rlxp1, /* <-- P0/P1 relaxation parameter */
78  const cs_real_t *dam, /* <-- Matrix diagonal */
79  const cs_real_t *xam /* <-- Matrix extra-diagonal terms */
80 );
81 
82 /*----------------------------------------------------------------------------
83  * Destroy a hierarchy of meshes starting from a fine mesh, keeping
84  * the corresponding system information for future calls.
85  *----------------------------------------------------------------------------*/
86 
87 void CS_PROCF(dsmlga, DSMLGA)
88 (
89  const char *cname, /* <-- variable name */
90  const cs_int_t *lname /* <-- variable name length */
91 );
92 
93 /*----------------------------------------------------------------------------
94  * General sparse linear system resolution
95  *----------------------------------------------------------------------------*/
96 
97 void CS_PROCF(resmgr, RESMGR)
98 (
99  const char *cname, /* <-- variable name */
100  const cs_int_t *lname, /* <-- variable name length */
101  const cs_int_t *iresds, /* <-- Descent smoother type:
102  0: pcg; 1: Jacobi; 2: cg-stab,
103  200: pcg_single reduction */
104  const cs_int_t *iresas, /* <-- Ascent smoother type:
105  0: pcg; 1: Jacobi; 2: cg-stab,
106  200: pcg_single reduction */
107  const cs_int_t *ireslp, /* <-- Coarse Resolution type:
108  0: pcg; 1: Jacobi; 2: cg-stab,
109  200: pcg_single reduction */
110  const cs_int_t *ipol, /* <-- Preconditioning polynomial degree
111  (0: diagonal, -1: none) */
112  const cs_int_t *ncymxp, /* <-- Max number of cycles */
113  const cs_int_t *nitmds, /* <-- Max number of iterations for descent */
114  const cs_int_t *nitmas, /* <-- Max number of iterations for ascent */
115  const cs_int_t *nitmap, /* <-- Max number of iterations for
116  coarsest solution */
117  const cs_int_t *iinvpe, /* <-- Indicator to cancel increments
118  in rotational periodicity (2) or
119  to exchange them as scalars (1) */
120  const cs_int_t *iwarnp, /* <-- Verbosity level */
121  cs_int_t *ncyclf, /* --> Number of cycles done */
122  cs_int_t *niterf, /* --> Number of iterations done */
123  const cs_real_t *epsilp, /* <-- Precision for iterative resolution */
124  const cs_real_t *rnorm, /* <-- Residue normalization */
125  cs_real_t *residu, /* --> Final non normalized residue */
126  const cs_real_t *rhs, /* <-- System right-hand side */
127  cs_real_t *vx /* <-> System solution */
128 );
129 
130 /*=============================================================================
131  * Public function prototypes
132  *============================================================================*/
133 
134 /*----------------------------------------------------------------------------
135  * Initialize multigrid solver API.
136  *----------------------------------------------------------------------------*/
137 
138 void
140 
141 /*----------------------------------------------------------------------------
142  * Finalize multigrid solver API.
143  *----------------------------------------------------------------------------*/
144 
145 void
147 
148 /*----------------------------------------------------------------------------
149  * Build a hierarchy of meshes starting from a fine mesh, for an
150  * ACM (Additive Corrective Multigrid) method.
151  *
152  * parameters:
153  * var_name <-- variable name
154  * verbosity <-- verbosity level
155  * postprocess_block_size <-- if > 0, postprocess coarsening, using
156  * coarse cell numbers modulo ncpost
157  * aggregation_limit <-- maximum allowed fine cells per coarse cell
158  * n_max_levels <-- maximum number of grid levels
159  * n_g_cells_min <-- global number of cells on coarsest grid
160  * under which no merging occurs
161  * p0p1_relax <-- p0/p1 relaxation_parameter
162  * symmetric <-- indicates if matrix coefficients are symmetric
163  * diag_block_size <-- block sizes for diagonal, or NULL
164  * extra_diag_block_size <-- Block sizes for extra diagonal, or NULL
165  * da <-- diagonal values (NULL if zero)
166  * xa <-- extradiagonal values (NULL if zero)
167  *----------------------------------------------------------------------------*/
168 
169 void
170 cs_multigrid_build(const char *var_name,
171  int verbosity,
172  int postprocess_block_size,
173  int aggregation_limit,
174  int n_max_levels,
175  cs_gnum_t n_g_cells_min,
176  double p0p1_relax,
177  bool symmetric,
178  const int *diag_block_size,
179  const int *extra_diag_block_size,
180  const cs_real_t *da,
181  const cs_real_t *xa);
182 
183 /*----------------------------------------------------------------------------
184  * Destroy a hierarchy of meshes starting from a fine mesh, keeping
185  * the corresponding system and postprocessing information for future calls.
186  *
187  * parameters:
188  * var_name <-- variable name
189  *----------------------------------------------------------------------------*/
190 
191 void
192 cs_multigrid_destroy(const char *var_name);
193 
194 /*----------------------------------------------------------------------------
195  * Sparse linear system resolution using multigrid.
196  *
197  * parameters:
198  * var_name <-- Variable name
199  * descent_smoother_type <-- Type of smoother for descent (PCG, Jacobi, ...)
200  * ascent_smoother_type <-- Type of smoother for ascent (PCG, Jacobi, ...)
201  * coarse_solver_type <-- Type of solver (PCG, Jacobi, ...)
202  * abort_on_divergence <-- Call errorhandler if divergence is detected
203  * poly_degree <-- Preconditioning polynomial degree (0: diagonal)
204  * rotation_mode <-- Halo update option for rotational periodicity
205  * verbosity <-- Verbosity level
206  * n_max_cycles <-- Maximum number of cycles
207  * n_max_iter_descent <-- Maximum nb. of iterations for descent phases
208  * n_max_iter_ascent <-- Maximum nb. of iterations for ascent phases
209  * n_max_iter_coarse <-- Maximum nb. of iterations for coarsest solution
210  * precision <-- Precision limit
211  * r_norm <-- Residue normalization
212  * n_cycles --> Number of cycles
213  * n_iter --> Number of iterations
214  * residue <-> Residue
215  * rhs <-- Right hand side
216  * vx --> System solution
217  * aux_size <-- Number of elements in aux_vectors
218  * aux_vectors --- Optional working area (allocation otherwise)
219  *
220  * returns:
221  * 1 if converged, 0 if not converged, -1 if not converged and maximum
222  * cycle number reached, -2 if divergence is detected.
223  *----------------------------------------------------------------------------*/
224 
225 int
226 cs_multigrid_solve(const char *var_name,
227  cs_sles_type_t descent_smoother_type,
228  cs_sles_type_t ascent_smoother_type,
229  cs_sles_type_t coarse_solver_type,
230  bool abort_on_divergence,
231  int poly_degree,
232  cs_halo_rotation_t rotation_mode,
233  int verbosity,
234  int n_max_cycles,
235  int n_max_iter_descent,
236  int n_max_iter_ascent,
237  int n_max_iter_coarse,
238  double precision,
239  double r_norm,
240  int *n_cycles,
241  int *n_iter,
242  double *residue,
243  const cs_real_t *rhs,
244  cs_real_t *vx,
245  size_t aux_size,
246  void *aux_vectors);
247 
248 /*----------------------------------------------------------------------------*/
249 
251 
252 #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
void resmgr(const char *cname, const cs_int_t *lname, const cs_int_t *iresds, const cs_int_t *iresas, const cs_int_t *ireslp, const cs_int_t *ipol, const cs_int_t *ncymxp, const cs_int_t *nitmds, const cs_int_t *nitmas, const cs_int_t *nitmap, const cs_int_t *iinvpe, const cs_int_t *iwarnp, cs_int_t *ncyclf, cs_int_t *niterf, const cs_real_t *epsilp, const cs_real_t *rnorm, cs_real_t *residu, const cs_real_t *rhs, cs_real_t *vx)
Definition: cs_multigrid.c:1857
void cs_multigrid_destroy(const char *var_name)
Definition: cs_multigrid.c:2282
#define BEGIN_C_DECLS
Definition: cs_defs.h:405
int cs_int_t
Fortran-compatible integer.
Definition: cs_defs.h:295
int cs_multigrid_solve(const char *var_name, cs_sles_type_t descent_smoother_type, cs_sles_type_t ascent_smoother_type, cs_sles_type_t coarse_solver_type, bool abort_on_divergence, int poly_degree, cs_halo_rotation_t rotation_mode, int verbosity, int n_max_cycles, int n_max_iter_descent, int n_max_iter_ascent, int n_max_iter_coarse, double precision, double r_norm, int *n_cycles, int *n_iter, double *residue, const cs_real_t *rhs, cs_real_t *vx, size_t aux_size, void *aux_vectors)
Definition: cs_multigrid.c:2340
void cs_multigrid_finalize(void)
Definition: cs_multigrid.c:1966
void cs_multigrid_initialize(void)
Definition: cs_multigrid.c:1957
#define END_C_DECLS
Definition: cs_defs.h:406
double cs_real_t
Definition: cs_defs.h:296
void cs_multigrid_build(const char *var_name, int verbosity, int postprocess_block_size, int aggregation_limit, int n_max_levels, cs_gnum_t n_g_cells_min, double p0p1_relax, bool symmetric, const int *diag_block_size, const int *extra_diag_block_size, const cs_real_t *da, const cs_real_t *xa)
Definition: cs_multigrid.c:2010
#define CS_PROCF(x, y)
Definition: cs_defs.h:419
void clmlga(const char *cname, const cs_int_t *lname, const cs_int_t *isym, const cs_int_t *ibsize, const cs_int_t *iesize, const cs_int_t *nagmax, const cs_int_t *ncpost, const cs_int_t *iwarnp, const cs_int_t *ngrmax, const cs_int_t *ncegrm, const cs_real_t *rlxp1, const cs_real_t *dam, const cs_real_t *xam)
Definition: cs_multigrid.c:1789
cs_sles_type_t
Definition: cs_sles.h:54
void dsmlga(const char *cname, const cs_int_t *lname)
Definition: cs_multigrid.c:1838