programmer's documentation
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
cs_sles.h
Go to the documentation of this file.
1 #ifndef __CS_SLES_H__
2 #define __CS_SLES_H__
3 
4 /*============================================================================
5  * Sparse Linear Equation Solvers
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_halo_perio.h"
36 #include "cs_matrix.h"
37 
38 /*----------------------------------------------------------------------------*/
39 
41 
42 /*============================================================================
43  * Macro definitions
44  *============================================================================*/
45 
46 /*============================================================================
47  * Type definitions
48  *============================================================================*/
49 
50 /*----------------------------------------------------------------------------
51  * Solver types
52  *----------------------------------------------------------------------------*/
53 
54 typedef enum {
55 
56  CS_SLES_PCG, /* Preconditionned conjugate gradient */
57  CS_SLES_PCG_SR, /* Preconditionned conjugate gradient, single reduction*/
58  CS_SLES_JACOBI, /* Jacobi */
59  CS_SLES_BICGSTAB, /* Bi-conjugate gradient stabilized */
60  CS_SLES_GMRES, /* Generalized minimal residual */
61  CS_SLES_N_TYPES /* Number of resolution algorithms */
62 
64 
65 /*============================================================================
66  * Global variables
67  *============================================================================*/
68 
69 /* Short names for matrix types */
70 
71 extern const char *cs_sles_type_name[];
72 
73 /*=============================================================================
74  * Public function prototypes for Fortran API
75  *============================================================================*/
76 
77 /*----------------------------------------------------------------------------
78  * General sparse linear system resolution
79  *----------------------------------------------------------------------------*/
80 
81 void CS_PROCF(reslin, RESLIN)
82 (
83  const char *cname, /* <-- variable name */
84  const cs_int_t *lname, /* <-- variable name length */
85  const cs_int_t *ncelet, /* <-- Number of cells, halo included */
86  const cs_int_t *ncel, /* <-- Number of local cells */
87  const cs_int_t *nfac, /* <-- Number of faces */
88  const cs_int_t *isym, /* <-- Symmetry indicator:
89  1: symmetric; 2: not symmetric */
90  const cs_int_t *ilved, /* <-- Interleaved indicator */
91  /* 1: interleaved; 2: not interleaved */
92  const cs_int_t *ibsize, /* <-- Block size of element ii,ii */
93  const cs_int_t *iesize, /* <-- Block size of element ij */
94  const cs_int_t *ireslp, /* <-- Resolution type:
95  0: pcg; 1: Jacobi; 2: cg-stab, 3: gmres,
96  200: pcg_sr */
97  const cs_int_t *ipol, /* <-- Preconditioning polynomial degree
98  (0: diagonal; -1: non-preconditionned) */
99  const cs_int_t *nitmap, /* <-- Number of max iterations */
100  const cs_int_t *iinvpe, /* <-- Indicator to cancel increments
101  in rotational periodicty (2) or
102  to exchange them as scalars (1) */
103  const cs_int_t *iwarnp, /* <-- Verbosity level */
104  cs_int_t *niterf, /* --> Number of iterations done */
105  const cs_real_t *epsilp, /* <-- Precision for iterative resolution */
106  const cs_real_t *rnorm, /* <-- Residue normalization */
107  cs_real_t *residu, /* --> Final non normalized residue */
108  const cs_int_t *ifacel, /* <-- Face -> cell connectivity */
109  const cs_real_t *dam, /* <-- Matrix diagonal */
110  const cs_real_t *xam, /* <-- Matrix extra-diagonal terms */
111  const cs_real_t *smbrp, /* <-- System right-hand side */
112  cs_real_t *vx /* <-> System solution */
113 );
114 
115 /*=============================================================================
116  * Public function prototypes
117  *============================================================================*/
118 
119 /*----------------------------------------------------------------------------
120  * Initialize sparse linear equation solver API.
121  *----------------------------------------------------------------------------*/
122 
123 void
124 cs_sles_initialize(void);
125 
126 /*----------------------------------------------------------------------------
127  * Finalize sparse linear equation solver API.
128  *----------------------------------------------------------------------------*/
129 
130 void
131 cs_sles_finalize(void);
132 
133 /*----------------------------------------------------------------------------
134  * Update sparse linear equation solver API in case of mesh modification.
135  *----------------------------------------------------------------------------*/
136 
137 void
138 cs_sles_update_mesh(void);
139 
140 #if defined(HAVE_MPI)
141 
142 /*----------------------------------------------------------------------------
143  * Set MPI communicator for dot products.
144  *----------------------------------------------------------------------------*/
145 
146 void
147 cs_sles_set_mpi_reduce_comm(MPI_Comm comm);
148 
149 #endif /* defined(HAVE_MPI) */
150 
151 /*----------------------------------------------------------------------------
152  * Test if a general sparse linear system needs solving or if the right-hand
153  * side is already zero within convergence criteria.
154  *
155  * The computed residue is also updated;
156  *
157  * parameters:
158  * var_name <-- Variable name
159  * solver_name <-- Name of solver
160  * n_rows <-- Number of (non ghost) rows in rhs
161  * verbosity <-- Verbosity level
162  * r_norm <-- Residue normalization
163  * residue <-> Residue
164  * rhs <-- Right hand side
165  *
166  * returns:
167  * 1 if solving is required, 0 if the rhs is already zero within tolerance
168  * criteria (precision of residue normalization)
169  *----------------------------------------------------------------------------*/
170 
171 int
172 cs_sles_needs_solving(const char *var_name,
173  const char *solver_name,
174  cs_int_t n_rows,
175  int verbosity,
176  double r_norm,
177  double *residue,
178  const cs_real_t *rhs);
179 
180 /*----------------------------------------------------------------------------
181  * General sparse linear system resolution.
182  *
183  * parameters:
184  * var_name <-- Variable name
185  * solver_type <-- Type of solver (PCG, Jacobi, ...)
186  * update_stats <-- Automatic solver statistics indicator
187  * symmetric <-- Symmetric coefficients indicator
188  * a <-- Matrix
189  * poly_degree <-- Preconditioning polynomial degree
190  * (0: diagonal; -1: non-preconditioned)
191  * rotation_mode <-- Halo update option for rotational periodicity
192  * verbosity <-- Verbosity level
193  * n_max_iter <-- Maximum number of iterations
194  * precision <-- Precision limit
195  * r_norm <-- Residue normalization
196  * n_iter --> Number of iterations
197  * residue <-> Residue
198  * rhs <-- Right hand side
199  * vx --> System solution
200  * aux_size <-- Number of elements in aux_vectors
201  * aux_vectors --- Optional working area (allocation otherwise)
202  *
203  * returns:
204  * 1 if converged, 0 if not converged, -1 if not converged and maximum
205  * iteration number reached, -2 if divergence is detected.
206  *----------------------------------------------------------------------------*/
207 
208 int
209 cs_sles_solve(const char *var_name,
210  cs_sles_type_t solver_type,
211  bool update_stats,
212  const cs_matrix_t *a,
213  int poly_degree,
214  cs_halo_rotation_t rotation_mode,
215  int verbosity,
216  int n_max_iter,
217  double precision,
218  double r_norm,
219  int *n_iter,
220  double *residue,
221  const cs_real_t *rhs,
222  cs_real_t *vx,
223  size_t aux_size,
224  void *aux_vectors);
225 
226 /*----------------------------------------------------------------------------
227  * Output default post-processing data for failed system convergence.
228  *
229  * parameters:
230  * var_name <-- Variable name
231  * mesh_id <-- id of error output mesh, or 0 if none
232  * rotation_mode <-- Halo update option for rotational periodicity
233  * a <-- Linear equation matrix
234  * rhs <-- Right hand side
235  * vx <-> Current system solution
236  *----------------------------------------------------------------------------*/
237 
238 void
239 cs_sles_post_error_output_def(const char *var_name,
240  int mesh_id,
241  cs_halo_rotation_t rotation_mode,
242  const cs_matrix_t *a,
243  const cs_real_t *rhs,
244  cs_real_t *vx);
245 
246 /*----------------------------------------------------------------------------
247  * Output post-processing variable for failed system convergence.
248  *
249  * parameters:
250  * var_name <-- Variable name
251  * diag_block_size <-- Block size for diagonal
252  * mesh_id <-- id of error output mesh, or 0 if none
253  * var <-- Variable values
254  *----------------------------------------------------------------------------*/
255 
256 void
257 cs_sles_post_error_output_var(const char *var_name,
258  int mesh_id,
259  int diag_block_size,
260  cs_real_t *var);
261 
262 /*----------------------------------------------------------------------------*/
263 
265 
266 #endif /* __CS_SLES_H__ */