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