8.1
general documentation
cs_thermal_model.h
Go to the documentation of this file.
1 #ifndef __CS_THERMAL_MODEL_H__
2 #define __CS_THERMAL_MODEL_H__
3 
4 /*============================================================================
5  * Base thermal model data.
6  *============================================================================*/
7 
8 /*
9  This file is part of code_saturne, a general-purpose CFD tool.
10 
11  Copyright (C) 1998-2023 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_field.h"
36 
37 /*----------------------------------------------------------------------------*/
38 
40 
41 /*=============================================================================
42  * Macro definitions
43  *============================================================================*/
44 
45 /*============================================================================
46  * Type definitions
47  *============================================================================*/
48 
49 /*----------------------------------------------------------------------------
50  * Thermal model type
51  *----------------------------------------------------------------------------*/
52 
53 typedef enum {
54 
61 
63 
64 typedef enum {
65 
69 
71 
72 /* thermal model descriptor */
73 /*--------------------------*/
74 
75 typedef struct {
76 
77  union {
79  int itherm;
80  };
81 
82  union {
83  cs_temperature_scale_t temperature_scale; /* Temperature scale */
84  int itpscl;
85  };
86 
87  /* Has kinetic source terme correction */
89  int cflt; /* compute the thermal cfl condition */
90  int cflp; /* compute the pressure cfl condition */
91  bool has_pdivu;
93 
95 
96 /*============================================================================
97  * Static global variables
98  *============================================================================*/
99 
100 /* Pointer to thermal model structure */
101 
103 
104 /*=============================================================================
105  * Public function prototypes
106  *============================================================================*/
107 
108 /*----------------------------------------------------------------------------
109  * Provide access to cs_glob_thermal_model
110  *
111  * needed to initialize structure with GUI
112  *----------------------------------------------------------------------------*/
113 
116 
117 /*----------------------------------------------------------------------------
118  * Return thermal field (temperature, enthalpy, total energy according to
119  * thermal model).
120  *
121  * returns:
122  * pointer to thermal field
123  *----------------------------------------------------------------------------*/
124 
125 cs_field_t *
127 
128 /*----------------------------------------------------------------------------
129  * Print the thermal model structure to setup.log.
130  *----------------------------------------------------------------------------*/
131 
132 void
134 
135 /*----------------------------------------------------------------------------*/
139 /*----------------------------------------------------------------------------*/
140 
141 void
143 
144 /*----------------------------------------------------------------------------*/
157 /*----------------------------------------------------------------------------*/
158 
159 void
161  const cs_real_t temp[],
162  const cs_real_t pres[],
163  const cs_real_t fracv[],
164  const cs_real_t fracm[],
165  const cs_real_t frace[],
166  cs_real_t dc2[]);
167 
168 /*----------------------------------------------------------------------------*/
177 /*----------------------------------------------------------------------------*/
178 
179 cs_real_t
181  cs_real_t temp,
182  cs_real_t yw);
183 
184 /*----------------------------------------------------------------------------*/
197 /*----------------------------------------------------------------------------*/
198 
199 cs_real_t
201  cs_real_t temp,
202  cs_real_t yw,
203  cs_real_t cpa,
204  cs_real_t cpv,
205  cs_real_t cpl,
206  cs_real_t l00);
207 
208 /*----------------------------------------------------------------------------*/
218 /*----------------------------------------------------------------------------*/
219 
220 void
222  const cs_real_t bmasfl[],
223  const cs_real_t vela[][3],
224  const cs_real_t vel[][3]);
225 
226 /*----------------------------------------------------------------------------*/
233 /*----------------------------------------------------------------------------*/
234 
235 void
237  const cs_real_t cromk[]);
238 
239 /*----------------------------------------------------------------------------*/
245 /*----------------------------------------------------------------------------*/
246 
247 void
249 
250 /*----------------------------------------------------------------------------*/
260 /*----------------------------------------------------------------------------*/
261 
262 void
263 cs_thermal_model_cflp(const cs_real_t croma[],
264  const cs_real_t trav2[][3],
265  const cs_real_t cvara_pr[],
266  const cs_real_t imasfl[],
267  cs_real_t cflp[]);
268 
269 /*----------------------------------------------------------------------------*/
281 /*----------------------------------------------------------------------------*/
282 
283 void
284 cs_thermal_model_cflt(const cs_real_t croma[],
285  const cs_real_t tempk[],
286  const cs_real_t tempka[],
287  const cs_real_t xcvv[],
288  const cs_real_t vel[][3],
289  const cs_real_t imasfl[],
290  cs_real_t cflt[restrict]);
291 
292 /*----------------------------------------------------------------------------*/
298 /*----------------------------------------------------------------------------*/
299 
300 void
302 
303 /*----------------------------------------------------------------------------*/
312 /*----------------------------------------------------------------------------*/
313 
314 void
316  const cs_real_t gradv[][3][3],
317  cs_real_t smbrs[]);
318 
319 /*----------------------------------------------------------------------------*/
333 /*----------------------------------------------------------------------------*/
334 
335 void
336 cs_thermal_model_newton_t(int method,
337  const cs_real_t *pk1,
338  const cs_real_t th_scal[],
339  const cs_real_t cvar_pr[],
340  const cs_real_t cvara_pr[],
341  const cs_real_t yw[],
342  cs_real_t yv[],
343  cs_real_t temp[]);
344 
345 /*----------------------------------------------------------------------------*/
351 /*----------------------------------------------------------------------------*/
352 
353 void
355 
356 /*----------------------------------------------------------------------------*/
357 
359 
360 #endif /* __CS_THERMAL_MODEL_H__ */
#define restrict
Definition: cs_defs.h:139
#define BEGIN_C_DECLS
Definition: cs_defs.h:514
double cs_real_t
Floating-point value.
Definition: cs_defs.h:319
#define END_C_DECLS
Definition: cs_defs.h:515
@ vel
Definition: cs_field_pointer.h:68
@ cp
Definition: cs_field_pointer.h:100
cs_real_t cs_thermal_model_demdt_ecsnt(cs_real_t pres, cs_real_t temp, cs_real_t yw, cs_real_t cpa, cs_real_t cpv, cs_real_t cpl, cs_real_t l00)
Compute the derivative of the internal energy related to the temperature at constant internal energy.
Definition: cs_thermal_model.c:475
cs_thermal_model_t * cs_get_glob_thermal_model(void)
Definition: cs_thermal_model.c:244
void cs_thermal_model_log_setup(void)
Definition: cs_thermal_model.c:256
void cs_thermal_model_cv(cs_real_t *xcvv)
Compute the isochoric heat capacity.
Definition: cs_thermal_model.c:823
void cs_thermal_model_add_kst(cs_real_t smbrs[])
Add the kinetic source term if needed.
Definition: cs_thermal_model.c:656
void cs_thermal_model_cflt(const cs_real_t croma[], const cs_real_t tempk[], const cs_real_t tempka[], const cs_real_t xcvv[], const cs_real_t vel[][3], const cs_real_t imasfl[], cs_real_t cflt[restrict])
Compute the CFL number related to the thermal equation.
Definition: cs_thermal_model.c:1307
void cs_thermal_model_c_square(const cs_real_t cp[], const cs_real_t temp[], const cs_real_t pres[], const cs_real_t fracv[], const cs_real_t fracm[], const cs_real_t frace[], cs_real_t dc2[])
Compute the inverse of the square of sound velocity multiplied by gamma.
Definition: cs_thermal_model.c:328
void cs_thermal_model_kinetic_st_finalize(const cs_real_t cromk1[], const cs_real_t cromk[])
Finalize the computation of the kinetic energy based source term.
Definition: cs_thermal_model.c:624
const cs_thermal_model_t * cs_glob_thermal_model
void cs_thermal_model_pdivu(cs_real_t smbrs[restrict])
Add the term pdivu to the thermal equation rhs.
Definition: cs_thermal_model.c:1120
void cs_thermal_model_cflp(const cs_real_t croma[], const cs_real_t trav2[][3], const cs_real_t cvara_pr[], const cs_real_t imasfl[], cs_real_t cflp[])
Compute the CFL number related to the pressure equation.
Definition: cs_thermal_model.c:687
void cs_thermal_model_kinetic_st_prepare(const cs_real_t imasfl[], const cs_real_t bmasfl[], const cs_real_t vela[][3], const cs_real_t vel[][3])
First pass to compute the contribution of the kinetic energy based source term from the prediction st...
Definition: cs_thermal_model.c:522
void cs_thermal_model_newton_t(int method, const cs_real_t *pk1, const cs_real_t th_scal[], const cs_real_t cvar_pr[], const cs_real_t cvara_pr[], const cs_real_t yw[], cs_real_t yv[], cs_real_t temp[])
Perform the Newton method to compute the temperature from the internal energy.
Definition: cs_thermal_model.c:922
cs_temperature_scale_t
Definition: cs_thermal_model.h:64
@ CS_TEMPERATURE_SCALE_CELSIUS
Definition: cs_thermal_model.h:68
@ CS_TEMPERATURE_SCALE_KELVIN
Definition: cs_thermal_model.h:67
@ CS_TEMPERATURE_SCALE_NONE
Definition: cs_thermal_model.h:66
cs_thermal_model_variable_t
Definition: cs_thermal_model.h:53
@ CS_THERMAL_MODEL_ENTHALPY
Definition: cs_thermal_model.h:57
@ CS_THERMAL_MODEL_INTERNAL_ENERGY
Definition: cs_thermal_model.h:59
@ CS_THERMAL_MODEL_TEMPERATURE
Definition: cs_thermal_model.h:56
@ CS_THERMAL_MODEL_N_TYPES
Definition: cs_thermal_model.h:60
@ CS_THERMAL_MODEL_NONE
Definition: cs_thermal_model.h:55
@ CS_THERMAL_MODEL_TOTAL_ENERGY
Definition: cs_thermal_model.h:58
void cs_thermal_model_init(void)
Initialize thermal variables if needed.
Definition: cs_thermal_model.c:305
cs_real_t cs_thermal_model_demdt(cs_real_t pres, cs_real_t temp, cs_real_t yw)
Compute the derivative of the internal energy related to the temperature at constant pressure.
Definition: cs_thermal_model.c:428
void cs_thermal_model_dissipation(const cs_real_t vistot[], const cs_real_t gradv[][3][3], cs_real_t smbrs[])
Compute and add the dissipation term of the thermal equation to its right hand side.
Definition: cs_thermal_model.c:876
cs_field_t * cs_thermal_model_field(void)
Definition: cs_thermal_model.c:216
real(c_double), pointer, save l00
Latent heat.
Definition: cstphy.f90:251
Field descriptor.
Definition: cs_field.h:131
Thermal model descriptor.
Definition: cs_thermal_model.h:75
int itpscl
Definition: cs_thermal_model.h:84
cs_thermal_model_variable_t thermal_variable
Definition: cs_thermal_model.h:78
cs_temperature_scale_t temperature_scale
Definition: cs_thermal_model.h:83
int itherm
Definition: cs_thermal_model.h:79
int cflp
Definition: cs_thermal_model.h:90
bool has_dissipation
Definition: cs_thermal_model.h:92
int cflt
Definition: cs_thermal_model.h:89
int has_kinetic_st
Definition: cs_thermal_model.h:88
bool has_pdivu
Definition: cs_thermal_model.h:91