8.3
general documentation
cs_equation_iterative_solve.h
Go to the documentation of this file.
1#ifndef __CS_EQUATION_ITERATIVE_SOLVE_H__
2#define __CS_EQUATION_ITERATIVE_SOLVE_H__
3
4/*============================================================================
5 * Scalar convection diffusion equation solver.
6 *============================================================================*/
7
8/*
9 This file is part of code_saturne, a general-purpose CFD tool.
10
11 Copyright (C) 1998-2024 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_defs.h"
35#include "cs_parameters.h"
36
38
39/*=============================================================================
40 * Public function prototypes
41 *============================================================================*/
42
43/*----------------------------------------------------------------------------*/
49/*----------------------------------------------------------------------------*/
50
51/*----------------------------------------------------------------------------*/
52/*
53 * \brief Solve an advection diffusion equation with source
54 * terms for one time step for the variable \f$ a \f$.
55 *
56 * The equation reads:
57 *
58 * \f[
59 * f_s^{imp}(a^{n+1}-a^n)
60 * + \divs \left( a^{n+1} \rho \vect{u} - \mu \grad a^{n+1} \right)
61 * = Rhs
62 * \f]
63 *
64 * This equation is rewritten as:
65 *
66 * \f[
67 * f_s^{imp} \delta a
68 * + \divs \left( \delta a \rho \vect{u} - \mu \grad \delta a \right)
69 * = Rhs^1
70 * \f]
71 *
72 * where \f$ \delta a = a^{n+1} - a^n\f$ and
73 * \f$ Rhs^1 = Rhs - \divs( a^n \rho \vect{u} - \mu \grad a^n)\f$
74 *
75 *
76 * It is in fact solved with the following iterative process:
77 *
78 * \f[
79 * f_s^{imp} \delta a^k
80 * + \divs \left(\delta a^k \rho \vect{u}-\mu\grad\delta a^k \right)
81 * = Rhs^k
82 * \f]
83 *
84 * where \f$Rhs^k=Rhs-f_s^{imp}(a^k-a^n)
85 * - \divs \left( a^k\rho\vect{u}-\mu\grad a^k \right)\f$
86 *
87 * Be careful, it is forbidden to modify \f$ f_s^{imp} \f$ here!
88 *
89 * Please refer to the
90 * <a href="../../theory.pdf#codits"><b>codits</b></a> section of the
91 * theory guide for more informations.
92 *
93 * \param[in] idtvar indicator of the temporal scheme
94 * \param[in] iterns external sub-iteration number
95 * \param[in] f_id field id (or -1)
96 * \param[in] name associated name if f_id < 0, or NULL
97 * \param[in] iescap compute the predictor indicator if 1
98 * \param[in] imucpp indicator
99 * - 0 do not multiply the convectiv term by Cp
100 * - 1 do multiply the convectiv term by Cp
101 * \param[in] normp Reference norm to solve the system (optional)
102 * if negative: recomputed here
103 * \param[in] eqp pointer to a cs_equation_param_t structure which
104 * contains variable calculation options
105 * \param[in] pvara variable at the previous time step
106 * \f$ a^n \f$
107 * \param[in] pvark variable at the previous sub-iteration
108 * \f$ a^k \f$.
109 * If you sub-iter on Navier-Stokes, then
110 * it allows to initialize by something else than
111 * pvara (usually pvar=pvara)
112 * \param[in] bc_coeffs boundary condition structure for the variable
113 * \param[in] i_massflux mass flux at interior faces
114 * \param[in] b_massflux mass flux at boundary faces
115 * \param[in] i_viscm \f$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} \f$
116 * at interior faces for the matrix
117 * \param[in] b_viscm \f$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} \f$
118 * at boundary faces for the matrix
119 * \param[in] i_visc \f$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} \f$
120 * at interior faces for the r.h.s.
121 * \param[in] b_visc \f$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} \f$
122 * at boundary faces for the r.h.s.
123 * \param[in] viscel symmetric cell tensor \f$ \tens{\mu}_\celli \f$
124 * \param[in] weighf internal face weight between cells i j in case
125 * of tensor diffusion
126 * \param[in] weighb boundary face weight for cells i in case
127 * of tensor diffusion
128 * \param[in] icvflb global indicator of boundary convection flux
129 * - 0 upwind scheme at all boundary faces
130 * - 1 imposed flux at some boundary faces
131 * \param[in] icvfli boundary face indicator array of convection flux
132 * - 0 upwind scheme
133 * - 1 imposed flux
134 * \param[in] rovsdt \f$ f_s^{imp} \f$
135 * \param[in,out] smbrp Right hand side \f$ Rhs^k \f$
136 * \param[in,out] pvar current variable
137 * \param[out] dpvar last variable increment
138 * \param[in] xcpp array of specific heat (Cp)
139 * \param[out] eswork prediction-stage error estimator
140 * (if iescap > 0)
141 */
142/*----------------------------------------------------------------------------*/
143
144void
146 int iterns,
147 int f_id,
148 const char *name,
149 int iescap,
150 int imucpp,
151 cs_real_t normp,
153 const cs_real_t pvara[],
154 const cs_real_t pvark[],
155 const cs_field_bc_coeffs_t *bc_coeffs,
156 const cs_real_t i_massflux[],
157 const cs_real_t b_massflux[],
158 const cs_real_t i_viscm[],
159 const cs_real_t b_viscm[],
160 const cs_real_t i_visc[],
161 const cs_real_t b_visc[],
162 cs_real_6_t viscel[],
163 const cs_real_2_t weighf[],
164 const cs_real_t weighb[],
165 int icvflb,
166 const int icvfli[],
167 const cs_real_t rovsdt[],
168 cs_real_t smbrp[],
169 cs_real_t pvar[],
170 cs_real_t dpvar[],
171 const cs_real_t xcpp[],
172 cs_real_t eswork[]);
173
174/*----------------------------------------------------------------------------*/
175/*
176 * \brief This function solves an advection diffusion equation with source terms
177 * for one time step for the vector variable \f$ \vect{a} \f$.
178 *
179 * The equation reads:
180 *
181 * \f[
182 * \tens{f_s}^{imp}(\vect{a}^{n+1}-\vect{a}^n)
183 * + \divv \left( \vect{a}^{n+1} \otimes \rho \vect {u}
184 * - \mu \gradt \vect{a}^{n+1}\right)
185 * = \vect{Rhs}
186 * \f]
187 *
188 * This equation is rewritten as:
189 *
190 * \f[
191 * \tens{f_s}^{imp} \delta \vect{a}
192 * + \divv \left( \delta \vect{a} \otimes \rho \vect{u}
193 * - \mu \gradt \delta \vect{a} \right)
194 * = \vect{Rhs}^1
195 * \f]
196 *
197 * where \f$ \delta \vect{a} = \vect{a}^{n+1} - \vect{a}^n\f$ and
198 * \f$ \vect{Rhs}^1 = \vect{Rhs}
199 * - \divv \left( \vect{a}^n \otimes \rho \vect{u}
200 * - \mu \gradt \vect{a}^n \right)\f$
201 *
202 *
203 * It is in fact solved with the following iterative process:
204 *
205 * \f[
206 * \tens{f_s}^{imp} \delta \vect{a}^k
207 * + \divv \left( \delta \vect{a}^k \otimes \rho \vect{u}
208 * - \mu \gradt \delta \vect{a}^k \right)
209 * = \vect{Rhs}^k
210 * \f]
211 *
212 * where \f$ \vect{Rhs}^k = \vect{Rhs}
213 * - \tens{f_s}^{imp} \left(\vect{a}^k-\vect{a}^n \right)
214 * - \divv \left( \vect{a}^k \otimes \rho \vect{u}
215 * - \mu \gradt \vect{a}^k \right)\f$
216 *
217 * Be careful, it is forbidden to modify \f$ \tens{f_s}^{imp} \f$ here!
218 *
219 * \param[in] idtvar indicator of the temporal scheme
220 * \param[in] iterns external sub-iteration number
221 * \param[in] f_id field id (or -1)
222 * \param[in] name associated name if f_id < 0, or NULL
223 * \param[in] ivisep indicator to take \f$ \divv
224 * \left(\mu \gradt \transpose{\vect{a}} \right)
225 * -2/3 \grad\left( \mu \dive \vect{a} \right)\f$
226 * - 1 take into account,
227 * - 0 otherwise
228 * \param[in] iescap compute the predictor indicator if >= 1
229 * \param[in] eqp pointer to a cs_equation_param_t structure which
230 * contains variable calculation options
231 * \param[in] pvara variable at the previous time step
232 * \f$ \vect{a}^n \f$
233 * \param[in] pvark variable at the previous sub-iteration
234 * \f$ \vect{a}^k \f$.
235 * If you sub-iter on Navier-Stokes, then
236 * it allows to initialize by something else than
237 * \c pvara (usually \c pvar= \c pvara)
238 * \param[in] bc_coeffs_v boundary condition structure for the variable
239 * \param[in] i_massflux mass flux at interior faces
240 * \param[in] b_massflux mass flux at boundary faces
241 * \param[in] i_viscm \f$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} \f$
242 * at interior faces for the matrix
243 * \param[in] b_viscm \f$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} \f$
244 * at boundary faces for the matrix
245 * \param[in] i_visc \f$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} \f$
246 * at interior faces for the r.h.s.
247 * \param[in] b_visc \f$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} \f$
248 * at boundary faces for the r.h.s.
249 * \param[in] i_secvis secondary viscosity at interior faces
250 * \param[in] b_secvis secondary viscosity at boundary faces
251 * \param[in] viscel symmetric cell tensor \f$ \tens{\mu}_\celli \f$
252 * \param[in] weighf internal face weight between cells i j in case
253 * of tensor diffusion
254 * \param[in] weighb boundary face weight for cells i in case
255 * of tensor diffusion
256 * \param[in] icvflb global indicator of boundary convection flux
257 * - 0 upwind scheme at all boundary faces
258 * - 1 imposed flux at some boundary faces
259 * \param[in] icvfli boundary face indicator array of convection flux
260 * - 0 upwind scheme
261 * - 1 imposed flux
262 * \param[in, out] fimp \f$ \tens{f_s}^{imp} \f$
263 * \param[in,out] smbrp Right hand side \f$ \vect{Rhs}^k \f$
264 * \param[in, out] pvar current variable
265 * \param[out] eswork prediction-stage error estimator
266 * (if iescap >= 0)
267 */
268/*----------------------------------------------------------------------------*/
269
270void
272 int iterns,
273 int f_id,
274 const char *name,
275 int ivisep,
276 int iescap,
278 const cs_real_t pvara[][3],
279 const cs_real_t pvark[][3],
280 const cs_field_bc_coeffs_t *bc_coeffs_v,
281 const cs_real_t i_massflux[],
282 const cs_real_t b_massflux[],
283 const cs_real_t i_viscm[],
284 const cs_real_t b_viscm[],
285 const cs_real_t i_visc[],
286 const cs_real_t b_visc[],
287 const cs_real_t i_secvis[],
288 const cs_real_t b_secvis[],
289 cs_real_t viscel[][6],
290 const cs_real_2_t weighf[],
291 const cs_real_t weighb[],
292 int icvflb,
293 const int icvfli[],
294 cs_real_t fimp[][3][3],
295 cs_real_t smbrp[][3],
296 cs_real_t pvar[][3],
297 cs_real_t eswork[][3]);
298
299/*----------------------------------------------------------------------------*/
300/*
301 * \brief This function solves an advection diffusion equation with source
302 * terms for one time step for the symmetric tensor variable
303 * \f$ \tens{\variat} \f$.
304 *
305 * The equation reads:
306 *
307 * \f[
308 * \tens{f_s}^{imp}(\tens{\variat}^{n+1}-\tens{\variat}^n)
309 * + \divt \left( \tens{\variat}^{n+1} \otimes \rho \vect {u}
310 * - \mu \gradtt \tens{\variat}^{n+1}\right)
311 * = \tens{Rhs}
312 * \f]
313 *
314 * This equation is rewritten as:
315 *
316 * \f[
317 * \tens{f_s}^{imp} \delta \tens{\variat}
318 * + \divt \left( \delta \tens{\variat} \otimes \rho \vect{u}
319 * - \mu \gradtt \delta \tens{\variat} \right)
320 * = \tens{Rhs}^1
321 * \f]
322 *
323 * where \f$ \delta \tens{\variat} = \tens{\variat}^{n+1} - \tens{\variat}^n\f$
324 * and \f$ \tens{Rhs}^1 = \tens{Rhs}
325 * - \divt \left( \tens{\variat}^n \otimes \rho \vect{u}
326 * - \mu \gradtt \tens{\variat}^n \right)\f$
327 *
328 *
329 * It is in fact solved with the following iterative process:
330 *
331 * \f[
332 * \tens{f_s}^{imp} \delta \tens{\variat}^k
333 * + \divt \left( \delta \tens{\variat}^k \otimes \rho \vect{u}
334 * - \mu \gradtt \delta \tens{\variat}^k \right)
335 * = \tens{Rhs}^k
336 * \f]
337 *
338 * where \f$ \tens{Rhs}^k = \tens{Rhs}
339 * - \tens{f_s}^{imp} \left(\tens{\variat}^k-\tens{\variat}^n \right)
340 * - \divt \left( \tens{\variat}^k \otimes \rho \vect{u}
341 * - \mu \gradtt \tens{\variat}^k \right)\f$
342 *
343 * Be careful, it is forbidden to modify \f$ \tens{f_s}^{imp} \f$ here!
344 *
345 * \param[in] idtvar indicator of the temporal scheme
346 * \param[in] f_id field id (or -1)
347 * \param[in] name associated name if f_id < 0, or NULL
348 * \param[in] eqp pointer to a cs_equation_param_t structure which
349 * contains variable calculation options
350 * \param[in] pvara variable at the previous time step
351 * \f$ \vect{a}^n \f$
352 * \param[in] pvark variable at the previous sub-iteration
353 * \f$ \vect{a}^k \f$.
354 * If you sub-iter on Navier-Stokes, then
355 * it allows to initialize by something else than
356 * pvara (usually pvar=pvara)
357 * \param[in] bc_coeffs_ts boundary condition structure for the variable
358 * \param[in] i_massflux mass flux at interior faces
359 * \param[in] b_massflux mass flux at boundary faces
360 * \param[in] i_viscm \f$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} \f$
361 * at interior faces for the matrix
362 * \param[in] b_viscm \f$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} \f$
363 * at boundary faces for the matrix
364 * \param[in] i_visc \f$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} \f$
365 * at interior faces for the r.h.s.
366 * \param[in] b_visc \f$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} \f$
367 * at boundary faces for the r.h.s.
368 * \param[in] viscel symmetric cell tensor \f$ \tens{\mu}_\celli \f$
369 * \param[in] weighf internal face weight between cells i j in case
370 * of tensor diffusion
371 * \param[in] weighb boundary face weight for cells i in case
372 * of tensor diffusion
373 * \param[in] icvflb global indicator of boundary convection flux
374 * - 0 upwind scheme at all boundary faces
375 * - 1 imposed flux at some boundary faces
376 * \param[in] icvfli boundary face indicator array of convection flux
377 * - 0 upwind scheme
378 * - 1 imposed flux
379 * \param[in,out] fimp \f$ \tens{f_s}^{imp} \f$
380 * \param[in,out] smbrp Right hand side \f$ \vect{Rhs}^k \f$
381 * \param[in,out] pvar current variable
382 */
383/*----------------------------------------------------------------------------*/
384
385void
387 int f_id,
388 const char *name,
390 const cs_real_t pvara[][6],
391 const cs_real_t pvark[][6],
392 const cs_field_bc_coeffs_t *bc_coeffs_ts,
393 const cs_real_t i_massflux[],
394 const cs_real_t b_massflux[],
395 const cs_real_t i_viscm[],
396 const cs_real_t b_viscm[],
397 const cs_real_t i_visc[],
398 const cs_real_t b_visc[],
399 cs_real_t viscel[][6],
400 const cs_real_2_t weighf[],
401 const cs_real_t weighb[],
402 int icvflb,
403 const int icvfli[],
404 cs_real_t fimp[][6][6],
405 cs_real_t smbrp[][6],
406 cs_real_t pvar[][6]);
407
409
410#endif /* __CS_EQUATION_ITERATIVE_SOLVE_H__ */
#define BEGIN_C_DECLS
Definition: cs_defs.h:542
double cs_real_t
Floating-point value.
Definition: cs_defs.h:342
cs_real_t cs_real_2_t[2]
vector of 2 floating-point values
Definition: cs_defs.h:358
cs_real_t cs_real_6_t[6]
vector of 6 floating-point values
Definition: cs_defs.h:361
#define END_C_DECLS
Definition: cs_defs.h:543
void cs_equation_iterative_solve_vector(int idtvar, int iterns, int f_id, const char *name, int ivisep, int iescap, cs_equation_param_t *eqp, const cs_real_t pvara[][3], const cs_real_t pvark[][3], const cs_field_bc_coeffs_t *bc_coeffs_v, const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_viscm[], const cs_real_t b_viscm[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t i_secvis[], const cs_real_t b_secvis[], cs_real_t viscel[][6], const cs_real_2_t weighf[], const cs_real_t weighb[], int icvflb, const int icvfli[], cs_real_t fimp[][3][3], cs_real_t smbrp[][3], cs_real_t pvar[][3], cs_real_t eswork[][3])
This function solves an advection diffusion equation with source terms for one time step for the vect...
Definition: cs_equation_iterative_solve.cpp:2082
void cs_equation_iterative_solve_tensor(int idtvar, int f_id, const char *name, cs_equation_param_t *eqp, const cs_real_t pvara[][6], const cs_real_t pvark[][6], const cs_field_bc_coeffs_t *bc_coeffs_ts, const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_viscm[], const cs_real_t b_viscm[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_t viscel[][6], const cs_real_2_t weighf[], const cs_real_t weighb[], int icvflb, const int icvfli[], cs_real_t fimp[][6][6], cs_real_t smbrp[][6], cs_real_t pvar[][6])
This function solves an advection diffusion equation with source terms for one time step for the symm...
Definition: cs_equation_iterative_solve.cpp:2226
void cs_equation_iterative_solve_scalar(int idtvar, int iterns, int f_id, const char *name, int iescap, int imucpp, cs_real_t normp, cs_equation_param_t *eqp, const cs_real_t pvara[], const cs_real_t pvark[], const cs_field_bc_coeffs_t *bc_coeffs, const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_viscm[], const cs_real_t b_viscm[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_6_t viscel[], const cs_real_2_t weighf[], const cs_real_t weighb[], int icvflb, const int icvfli[], const cs_real_t rovsdt[], cs_real_t smbrp[], cs_real_t pvar[], cs_real_t dpvar[], const cs_real_t xcpp[], cs_real_t eswork[])
Solve an advection diffusion equation with source terms for one time step for the variable .
Definition: cs_equation_iterative_solve.cpp:1207
integer(c_int), pointer, save idtvar
option for a variable time step
Definition: optcal.f90:193
Set of parameters to handle an unsteady convection-diffusion-reaction equation with term sources.
Definition: cs_equation_param.h:192
Field boundary condition descriptor (for variables)
Definition: cs_field.h:104