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-2013 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  const cs_int_t *iresas, /* <-- Ascent smoother type:
104  0: pcg; 1: Jacobi; 2: cg-stab */
105  const cs_int_t *ireslp, /* <-- Coarse Resolution type:
106  0: pcg; 1: Jacobi; 2: cg-stab */
107  const cs_int_t *ipol, /* <-- Preconditioning polynomial degree
108  (0: diagonal) */
109  const cs_int_t *ncymxp, /* <-- Max number of cycles */
110  const cs_int_t *nitmds, /* <-- Max number of iterations for descent */
111  const cs_int_t *nitmas, /* <-- Max number of iterations for ascent */
112  const cs_int_t *nitmap, /* <-- Max number of iterations for
113  coarsest solution */
114  const cs_int_t *iinvpe, /* <-- Indicator to cancel increments
115  in rotational periodicity (2) or
116  to exchange them as scalars (1) */
117  const cs_int_t *iwarnp, /* <-- Verbosity level */
118  cs_int_t *ncyclf, /* --> Number of cycles done */
119  cs_int_t *niterf, /* --> Number of iterations done */
120  const cs_real_t *epsilp, /* <-- Precision for iterative resolution */
121  const cs_real_t *rnorm, /* <-- Residue normalization */
122  cs_real_t *residu, /* --> Final non normalized residue */
123  const cs_real_t *rhs, /* <-- System right-hand side */
124  cs_real_t *vx /* <-> System solution */
125 );
126 
127 /*=============================================================================
128  * Public function prototypes
129  *============================================================================*/
130 
131 /*----------------------------------------------------------------------------
132  * Initialize multigrid solver API.
133  *----------------------------------------------------------------------------*/
134 
135 void
137 
138 /*----------------------------------------------------------------------------
139  * Finalize multigrid solver API.
140  *----------------------------------------------------------------------------*/
141 
142 void
144 
145 /*----------------------------------------------------------------------------
146  * Build a hierarchy of meshes starting from a fine mesh, for an
147  * ACM (Additive Corrective Multigrid) method.
148  *
149  * parameters:
150  * var_name <-- variable name
151  * verbosity <-- verbosity level
152  * postprocess_block_size <-- if > 0, postprocess coarsening, using
153  * coarse cell numbers modulo ncpost
154  * aggregation_limit <-- maximum allowed fine cells per coarse cell
155  * n_max_levels <-- maximum number of grid levels
156  * n_g_cells_min <-- global number of cells on coarsest grid
157  * under which no merging occurs
158  * p0p1_relax <-- p0/p1 relaxation_parameter
159  * symmetric <-- indicates if matrix coefficients are symmetric
160  * diag_block_size <-- block sizes for diagonal, or NULL
161  * extra_diag_block_size <-- Block sizes for extra diagonal, or NULL
162  * da <-- diagonal values (NULL if zero)
163  * xa <-- extradiagonal values (NULL if zero)
164  *----------------------------------------------------------------------------*/
165 
166 void
167 cs_multigrid_build(const char *var_name,
168  int verbosity,
169  int postprocess_block_size,
170  int aggregation_limit,
171  int n_max_levels,
172  cs_gnum_t n_g_cells_min,
173  double p0p1_relax,
174  bool symmetric,
175  const int *diag_block_size,
176  const int *extra_diag_block_size,
177  const cs_real_t *da,
178  const cs_real_t *xa);
179 
180 /*----------------------------------------------------------------------------
181  * Destroy a hierarchy of meshes starting from a fine mesh, keeping
182  * the corresponding system and postprocessing information for future calls.
183  *
184  * parameters:
185  * var_name <-- variable name
186  *----------------------------------------------------------------------------*/
187 
188 void
189 cs_multigrid_destroy(const char *var_name);
190 
191 /*----------------------------------------------------------------------------
192  * Sparse linear system resolution using multigrid.
193  *
194  * parameters:
195  * var_name <-- Variable name
196  * descent_smoother_type <-- Type of smoother for descent (PCG, Jacobi, ...)
197  * ascent_smoother_type <-- Type of smoother for ascent (PCG, Jacobi, ...)
198  * coarse_solver_type <-- Type of solver (PCG, Jacobi, ...)
199  * abort_on_divergence <-- Call errorhandler if divergence is detected
200  * poly_degree <-- Preconditioning polynomial degree (0: diagonal)
201  * rotation_mode <-- Halo update option for rotational periodicity
202  * verbosity <-- Verbosity level
203  * n_max_cycles <-- Maximum number of cycles
204  * n_max_iter_descent <-- Maximum nb. of iterations for descent phases
205  * n_max_iter_ascent <-- Maximum nb. of iterations for ascent phases
206  * n_max_iter_coarse <-- Maximum nb. of iterations for coarsest solution
207  * precision <-- Precision limit
208  * r_norm <-- Residue normalization
209  * n_cycles --> Number of cycles
210  * n_iter --> Number of iterations
211  * residue <-> Residue
212  * rhs <-- Right hand side
213  * vx --> System solution
214  * aux_size <-- Number of elements in aux_vectors
215  * aux_vectors --- Optional working area (allocation otherwise)
216  *
217  * returns:
218  * 1 if converged, 0 if not converged, -1 if not converged and maximum
219  * cycle number reached, -2 if divergence is detected.
220  *----------------------------------------------------------------------------*/
221 
222 int
223 cs_multigrid_solve(const char *var_name,
224  cs_sles_type_t descent_smoother_type,
225  cs_sles_type_t ascent_smoother_type,
226  cs_sles_type_t coarse_solver_type,
227  bool abort_on_divergence,
228  int poly_degree,
229  cs_halo_rotation_t rotation_mode,
230  int verbosity,
231  int n_max_cycles,
232  int n_max_iter_descent,
233  int n_max_iter_ascent,
234  int n_max_iter_coarse,
235  double precision,
236  double r_norm,
237  int *n_cycles,
238  int *n_iter,
239  double *residue,
240  const cs_real_t *rhs,
241  cs_real_t *vx,
242  size_t aux_size,
243  void *aux_vectors);
244 
245 /*----------------------------------------------------------------------------*/
246 
248 
249 #endif /* __CS_MULTIGRID_H__ */