1 /*============================================================================
2 * Base gas mix data.
3 *============================================================================*/
4
5 /*
6 This file is part of Code_Saturne, a general-purpose CFD tool.
7
8 Copyright (C) 1998-2021 EDF S.A.
9
10 This program is free software; you can redistribute it and/or modify it under
11 the terms of the GNU General Public License as published by the Free Software
12 Foundation; either version 2 of the License, or (at your option) any later
13 version.
14
15 This program is distributed in the hope that it will be useful, but WITHOUT
16 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
18 details.
19
20 You should have received a copy of the GNU General Public License along with
21 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
22 Street, Fifth Floor, Boston, MA 02110-1301, USA.
23 */
24
25 /*----------------------------------------------------------------------------*/
26
27 #include "cs_defs.h"
28
29 /*----------------------------------------------------------------------------*/
30
31 /*----------------------------------------------------------------------------
32 * Standard C library headers
33 *----------------------------------------------------------------------------*/
34
35 #include <assert.h>
36 #include <stdio.h>
37 #include <stdlib.h>
38 #include <string.h>
39
40 /*----------------------------------------------------------------------------
41 * Local headers
42 *----------------------------------------------------------------------------*/
43
44 #include "cs_log.h"
45
46 #include "bft_mem.h"
47 #include "bft_error.h"
48 #include "bft_printf.h"
49
50 #include "cs_field.h"
51 #include "cs_parall.h"
52 #include "cs_physical_model.h"
53 #include "cs_post.h"
54
55 /*----------------------------------------------------------------------------
56 * Header for the current file
57 *----------------------------------------------------------------------------*/
58
59 #include "cs_gas_mix.h"
60
61 /*----------------------------------------------------------------------------*/
62
63 BEGIN_C_DECLS
64
65 /*=============================================================================
66 * Additional doxygen documentation
67 *============================================================================*/
68
69 /*!
70 \file cs_gas_mix.c
71 Base gas mix data.
72 */
73 /*----------------------------------------------------------------------------*/
74
75 /*! \struct cs_gas_mix_t
76
77 \brief Gas mix descriptor.
78
79 Members of this structure are publicly accessible, to allow for
80 concise syntax, as they are expected to be used in many places.
81
82 \var cs_gas_mix_t::n_species
83 number of species in the gas mix
84 */
85
86 /*----------------------------------------------------------------------------*/
87
88 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
89
90 /*=============================================================================
91 * Macro definitions
92 *============================================================================*/
93
94 /*============================================================================
95 * Type definitions
96 *============================================================================*/
97
98 /*============================================================================
99 * Static global variables
100 *============================================================================*/
101
102 /* main gas mix options and associated pointer */
103
104 static cs_gas_mix_t _gas_mix = {
105 .n_species = 0,
106 .n_species_solved = 0,
107 .species_to_field_id = NULL
108 };
109
110 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
111
112 /*============================================================================
113 * Static global variables
114 *============================================================================*/
115
116 /* Default variable compute options */
117
118 static cs_gas_mix_species_prop_t _gas_mix_species_prop =
119 {
120 -1., /* molar mass */
121 -1., /* specific heat */
122 -1., /* volume diffusion */
123 -1., /* dynamic viscosity a */
124 -1., /* dynamic viscosity b */
125 -1., /* thermal conductivity a */
126 -1., /* thermal conductivity b */
127 -1., /* reference viscosity (Sutherland) */
128 -1., /* reference conductivity (Sutherland) */
129 -1., /* reference temperature for viscosity */
130 -1., /* reference temperature for conductivity */
131 -1., /* Sutherland temperature for viscosity */
132 -1., /* Sutherland temperature for conductivity */
133 };
134
135 /*============================================================================
136 * Global variables
137 *============================================================================*/
138
139 const cs_gas_mix_t *cs_glob_gas_mix = &_gas_mix;
140
141 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
142
143 /*============================================================================
144 * Prototypes for functions intended for use only by Fortran wrappers.
145 * (descriptions follow, with function bodies).
146 *============================================================================*/
147
148 void
149 cs_f_gas_mix_get_pointers(int **nscasp);
150
151 int
152 cs_f_gas_mix_species_to_field_id(int sp_id);
153
154 /*============================================================================
155 * Private function definitions
156 *============================================================================*/
157
158 /* Log values of the structure */
159
160 static void
_log_func_gas_mix_species_prop(const void * t)161 _log_func_gas_mix_species_prop(const void *t)
162 {
163 const char fmt[] = N_(" %-19s %-12.3g\n");
164 const cs_gas_mix_species_prop_t *_t = (const void *)t;
165 cs_log_printf(CS_LOG_SETUP, fmt, "mol_mas ", _t->mol_mas);
166 cs_log_printf(CS_LOG_SETUP, fmt, "cp ", _t->cp);
167 cs_log_printf(CS_LOG_SETUP, fmt, "vol_diff", _t->vol_dif);
168 cs_log_printf(CS_LOG_SETUP, fmt, "mu_a ", _t->mu_a);
169 cs_log_printf(CS_LOG_SETUP, fmt, "mu_b ", _t->mu_b);
170 cs_log_printf(CS_LOG_SETUP, fmt, "lambda_a", _t->lambda_a);
171 cs_log_printf(CS_LOG_SETUP, fmt, "lambda_b", _t->lambda_b);
172 cs_log_printf(CS_LOG_SETUP, fmt, "muref ", _t->muref);
173 cs_log_printf(CS_LOG_SETUP, fmt, "lamref ", _t->lamref);
174 cs_log_printf(CS_LOG_SETUP, fmt, "trefmu ", _t->trefmu);
175 cs_log_printf(CS_LOG_SETUP, fmt, "treflam ", _t->treflam);
176 cs_log_printf(CS_LOG_SETUP, fmt, "smu ", _t->smu);
177 cs_log_printf(CS_LOG_SETUP, fmt, "slam ", _t->slam);
178 }
179
180 /* Log default values of the structure */
181
182 static void
_log_func_default_gas_mix_species_prop(const void * t)183 _log_func_default_gas_mix_species_prop(const void *t)
184 {
185 const char fmt[] = " %-19s %-12.3g %s\n";
186 const cs_gas_mix_species_prop_t *_t = (const void *)t;
187 cs_log_printf(CS_LOG_SETUP, fmt, "mol_mas ", _t->mol_mas,
188 _("Molar mass"));
189 cs_log_printf(CS_LOG_SETUP, fmt, "cp ", _t->cp,
190 _("Specific heat"));
191 cs_log_printf(CS_LOG_SETUP, fmt, "vol_diff", _t->vol_dif,
192 _("Volume diffusion"));
193 cs_log_printf(CS_LOG_SETUP, fmt, "mu_a ", _t->mu_a,
194 _("Dynamic viscosity a"));
195 cs_log_printf(CS_LOG_SETUP, fmt, "mu_b ", _t->mu_b,
196 _("Dynamic viscosity b"));
197 cs_log_printf(CS_LOG_SETUP, fmt, "lambda_a", _t->lambda_a,
198 _("Thermal conductivity a"));
199 cs_log_printf(CS_LOG_SETUP, fmt, "lambda_b", _t->lambda_b,
200 _("Thermal conductivity b"));
201 cs_log_printf(CS_LOG_SETUP, fmt, "muref ", _t->muref,
202 _("Reference thermal viscosity (Sutherland)"));
203 cs_log_printf(CS_LOG_SETUP, fmt, "lamref ", _t->lamref,
204 _("Reference thermal conductivity (Sutherland)"));
205 cs_log_printf(CS_LOG_SETUP, fmt, "trefmu ", _t->trefmu,
206 _("Reference temperature (Sutherland for viscosity)"));
207 cs_log_printf(CS_LOG_SETUP, fmt, "treflam ", _t->treflam,
208 _("Reference temperature (Sutherland conductivity)"));
209 cs_log_printf(CS_LOG_SETUP, fmt, "smu ", _t->smu,
210 _("Sutherland temperature for viscosity"));
211 cs_log_printf(CS_LOG_SETUP, fmt, "slam ", _t->slam,
212 _("Sutherland temperature for conductivity"));
213 }
214
215 /*----------------------------------------------------------------------------*/
216 /*!
217 * \brief Add a species field to the gas mix (set of fields).
218 *
219 * The species to field mapping is usually done in the order of addition
220 * of species, but non-resolved (deduced) fields are placed last in the
221 * list; the order of solved species relative to each other is maintained,
222 * as is that of deduced species. In most cases, we have at most one deduced
223 * species, which is mapped last.
224 *
225 * \param[in] f pointer to the field
226 */
227 /*----------------------------------------------------------------------------*/
228
229 static void
_map_field(const cs_field_t * f)230 _map_field(const cs_field_t *f)
231 {
232 /* Check if already mapped */
233
234 for (int i = 0; i < _gas_mix.n_species; i++) {
235 if (_gas_mix.species_to_field_id[i] == f->id)
236 return;
237 }
238
239 int is_solved = f->type & CS_FIELD_VARIABLE;
240 int insert_id = (is_solved) ?
241 _gas_mix.n_species_solved : _gas_mix.n_species;
242
243 int n_species_ini = _gas_mix.n_species;
244 _gas_mix.n_species++;
245
246 if (is_solved)
247 _gas_mix.n_species_solved += 1;
248
249 BFT_REALLOC(_gas_mix.species_to_field_id, _gas_mix.n_species, int);
250
251 /* If we need to insert a solved variable with non-solved fields
252 already mapped, shift the non-solved field map positions */
253
254 for (int i = n_species_ini; i > insert_id; i--)
255 _gas_mix.species_to_field_id[i] = _gas_mix.species_to_field_id[i-1];
256
257 _gas_mix.species_to_field_id[insert_id] = f->id;
258 }
259
260 /*----------------------------------------------------------------------------*/
261 /*!
262 * \brief Set predefined gas properties for a field based on its name.
263 *
264 * This only applies to fields for which properties have not already been
265 * set (i.e. for which \ref cs_gas_mix_add_species_with_properties
266 * has not been called), so as to complete the property definitions
267 * based on predefined values.
268 *
269 * \param[in, out] f_id field id of an already created scalar model field
270 */
271 /*----------------------------------------------------------------------------*/
272
273 static void
_set_predefined_property(cs_field_t * f)274 _set_predefined_property(cs_field_t *f)
275 {
276 int k_id = cs_gas_mix_get_field_key();
277
278 cs_gas_mix_species_prop_t gmp;
279 cs_field_get_key_struct(f, k_id, &gmp);
280
281 if (gmp.mol_mas > 0)
282 return;
283
284 if (strcmp(f->name, "y_o2") == 0) {
285 gmp.mol_mas= 0.032;
286 gmp.cp = 930.;
287 gmp.vol_dif = 19.70;
288 gmp.mu_a = 5.086522e-8;
289 gmp.mu_b = 5.512391e-6;
290 gmp.lambda_a = 6.2e-5;
291 gmp.lambda_b = 8.1e-3;
292 gmp.muref = 1.919e-5;
293 gmp.lamref = 0.0244;
294 gmp.trefmu = 273.;
295 gmp.treflam = 273.;
296 gmp.smu = 139.;
297 gmp.slam = 240.;
298 }
299
300 else if (strcmp(f->name, "y_n2") == 0) {
301 gmp.mol_mas = 0.028;
302 gmp.cp = 1042.;
303 gmp.vol_dif = 19.70;
304 gmp.mu_a = 4.210130e-8;
305 gmp.mu_b = 5.494348e-6;
306 gmp.lambda_a = 6.784141e-5;
307 gmp.lambda_b = 5.564317e-3;
308 gmp.muref = 1.663e-5;
309 gmp.lamref = 0.0242;
310 gmp.trefmu = 273.;
311 gmp.treflam = 273.;
312 gmp.smu = 107.;
313 gmp.slam = 150.;
314 }
315
316 else if (strcmp(f->name, "y_he") == 0) {
317 gmp.mol_mas = 0.004;
318 gmp.cp = 5194.;
319 gmp.vol_dif = 2.67;
320 gmp.mu_a = 18.5752e-6;
321 gmp.mu_b = 0.0;
322 gmp.lambda_a = 0.144;
323 gmp.lambda_b = 0.0;
324 gmp.muref = 1.874e-5;
325 gmp.lamref = 0.647;
326 gmp.trefmu = 273.;
327 gmp.treflam = 273.;
328 gmp.smu = 78.;
329 gmp.slam = 78.;
330 }
331
332 else if (strcmp(f->name, "y_h2") == 0) {
333 gmp.mol_mas = 0.002;
334 gmp.cp = 14560.;
335 gmp.vol_dif = 6.12;
336 gmp.mu_a = 1.93e-9;
337 gmp.mu_b = 8.40e-6;
338 gmp.lambda_a = 4.431e-4;
339 gmp.lambda_b = 5.334e-2;
340 gmp.muref = 8.411e-6;
341 gmp.lamref = 0.0168;
342 gmp.trefmu = 273.;
343 gmp.treflam = 273.;
344 gmp.smu = 97.;
345 gmp.slam = 120.;
346 }
347
348 else if (strcmp(f->name, "y_h2o_g") == 0) {
349 gmp.mol_mas = 0.018;
350 gmp.cp = 2060.;
351 gmp.vol_dif = 13.1;
352 gmp.mu_a = 3.8496e-8;
353 gmp.mu_b = 8.2997e-6;
354 gmp.lambda_a = 7.6209e-5;
355 gmp.lambda_b = 0.016949;
356 gmp.muref = 1.12e-5;
357 gmp.lamref = 0.0181;
358 gmp.trefmu = 350.;
359 gmp.treflam = 300.;
360 gmp.smu = 1064.;
361 gmp.slam = 2200.;
362 }
363
364 else {
365 bft_error(__FILE__, __LINE__, 0,
366 _("%s: no predefined properties for field %s."),
367 __func__, f->name);
368
369 }
370
371 cs_field_set_key_struct(f, k_id, &gmp);
372 }
373
374 /*============================================================================
375 * Fortran wrapper function definitions
376 *============================================================================*/
377
378 /*----------------------------------------------------------------------------
379 * Get pointers to members of the global gas mix structure.
380 *
381 * This function is intended for use by Fortran wrappers, and
382 * enables mapping to Fortran global pointers.
383 *
384 * parameters:
385 * nscasp --> pointer to cs_glob_gas_mix->n_species
386 *----------------------------------------------------------------------------*/
387
388 void
cs_f_gas_mix_get_pointers(int ** nscasp)389 cs_f_gas_mix_get_pointers(int **nscasp)
390 {
391 *nscasp = &(_gas_mix.n_species_solved);
392 }
393
394 /*----------------------------------------------------------------------------
395 * Return field id matching a given species id.
396 *
397 * This function is intended for use by Fortran wrappers.
398 *
399 * parameters:
400 * sp_id --> matching species id (1-based,
401 * deduced species if sp_id == nscasp + 1)
402 *----------------------------------------------------------------------------*/
403
404 int
cs_f_gas_mix_species_to_field_id(int sp_id)405 cs_f_gas_mix_species_to_field_id(int sp_id)
406 {
407 int retval = -1;
408
409 int _sp_id = sp_id-1;
410
411 if (_sp_id >= 0 && _sp_id < _gas_mix.n_species)
412 retval = _gas_mix.species_to_field_id[_sp_id];
413
414 return retval;
415 }
416
417 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
418
419 /*=============================================================================
420 * Public function definitions
421 *============================================================================*/
422
423 /*----------------------------------------------------------------------------*/
424 /*!
425 * \brief Get the field key for gas mix properties.
426 *
427 * \return field key id for gas mix properties
428 */
429 /*----------------------------------------------------------------------------*/
430
431 int
cs_gas_mix_get_field_key(void)432 cs_gas_mix_get_field_key(void)
433 {
434 static int k_id = -1;
435
436 /* Structure containing physical properties relative to
437 species scalars used by the gas mixture modelling */
438
439 if (k_id < 0) {
440
441 const char key[] = "gas_mix_species_prop";
442
443 cs_field_define_key_struct(key,
444 &_gas_mix_species_prop,
445 _log_func_gas_mix_species_prop,
446 _log_func_default_gas_mix_species_prop,
447 NULL,
448 sizeof(cs_gas_mix_species_prop_t),
449 0);
450
451 k_id = cs_field_key_id_try(key);
452
453 }
454
455 return k_id;
456 }
457
458 /*----------------------------------------------------------------------------*/
459 /*!
460 * \brief Add a species field to the gas mix (set of fields).
461 *
462 * \param[in] f_id field id of an already created scalar model field
463 */
464 /*----------------------------------------------------------------------------*/
465
466 void
cs_gas_mix_add_species(int f_id)467 cs_gas_mix_add_species(int f_id)
468 {
469 if (cs_glob_physical_model_flag[CS_GAS_MIX] == -1)
470 bft_error(__FILE__, __LINE__, 0,
471 _("No gas species can be added."
472 " The gas mix model is not enabled.\n"));
473
474 cs_field_t *f = cs_field_by_id(f_id);
475
476 if ( strcmp(f->name, "y_o2") != 0
477 && strcmp(f->name, "y_n2") != 0
478 && strcmp(f->name, "y_he") != 0
479 && strcmp(f->name, "y_h2") != 0)
480 bft_error(__FILE__, __LINE__, 0,
481 _("Only the species having the following field names "
482 "can be added to a gas mix:\n"
483 "y_o2, y_n2, y_he, y_h2\n"));
484
485 _map_field(f);
486
487 _set_predefined_property(f);
488 }
489
490 /*----------------------------------------------------------------------------*/
491 /*!
492 * \brief Add a species field to the gas mix (set of fields).
493 *
494 * \param[in] f_id id of field representing species mixture fraction.
495 * \param[in] mol_mass molar mass
496 * \param[in] cp specific heat
497 * \param[in] col_diff volume diffusion
498 * \param[in] mu_a dynamic viscosity a
499 * \param[in] mu_b dynamic viscosity b
500 * \param[in] lambda_a thermal conductivity a
501 * \param[in] lambda_b thermal conductivity b
502 * \param[in] mu_ref reference viscosity (Sutherland)
503 * \param[in] lambda_ref reference conductivity (Sutherland)
504 * \param[in] tref_mu reference temperature for viscosity
505 * \param[in] tref_lambda reference temperature for conductivity
506 * \param[in] s_mu Sutherland temperature for viscosity
507 * \param[in] s_lambda Sutherland temperature for conductivity
508 */
509 /*----------------------------------------------------------------------------*/
510
511 void
cs_gas_mix_add_species_with_properties(int f_id,cs_real_t mol_mass,cs_real_t cp,cs_real_t vol_diff,cs_real_t mu_a,cs_real_t mu_b,cs_real_t lambda_a,cs_real_t lambda_b,cs_real_t mu_ref,cs_real_t lambda_ref,cs_real_t tref_mu,cs_real_t tref_lambda,cs_real_t s_mu,cs_real_t s_lambda)512 cs_gas_mix_add_species_with_properties(int f_id,
513 cs_real_t mol_mass,
514 cs_real_t cp,
515 cs_real_t vol_diff,
516 cs_real_t mu_a,
517 cs_real_t mu_b,
518 cs_real_t lambda_a,
519 cs_real_t lambda_b,
520 cs_real_t mu_ref,
521 cs_real_t lambda_ref,
522 cs_real_t tref_mu,
523 cs_real_t tref_lambda,
524 cs_real_t s_mu,
525 cs_real_t s_lambda)
526 {
527 cs_field_t *f = cs_field_by_id(f_id);
528
529 _map_field(f);
530
531 cs_gas_mix_species_prop_t gmp
532 = {.mol_mas = mol_mass,
533 .cp = cp,
534 .vol_dif = vol_diff,
535 .mu_a = mu_a,
536 .mu_b = mu_b,
537 .lambda_a = lambda_a,
538 .lambda_b = lambda_b,
539 .muref = mu_ref,
540 .lamref = lambda_ref,
541 .trefmu = tref_mu,
542 .treflam = tref_lambda,
543 .smu = s_mu,
544 .slam = s_lambda};
545
546 int k_id = cs_gas_mix_get_field_key();
547
548 cs_field_set_key_struct(f, k_id, &gmp);
549
550 /* Set model flag.
551 Reserve lower values (0-5 currently used) for predefined cases;
552 TODO: only 0/1 should be needed here, using specific functions
553 or an ENUM for pre-definition */
554
555 cs_glob_physical_model_flag[CS_GAS_MIX] = CS_GAS_MIX_USER;
556 }
557
558 /*----------------------------------------------------------------------------*/
559 /*!
560 * \brief Add property fields specific to a gas mix.
561 */
562 /*----------------------------------------------------------------------------*/
563
564 void
cs_gas_mix_add_property_fields(void)565 cs_gas_mix_add_property_fields(void)
566 {
567 cs_field_t *f;
568
569 int field_type = CS_FIELD_INTENSIVE | CS_FIELD_PROPERTY;
570 const int klbl = cs_field_key_id("label");
571 const int keyvis = cs_field_key_id("post_vis");
572 const int keylog = cs_field_key_id("log");
573
574 bool add_deduced = true;
575
576 char name[32], label[32];
577
578 switch(cs_glob_physical_model_flag[CS_GAS_MIX]) {
579
580 case CS_GAS_MIX_AIR_HELIUM:
581 strncpy(name, "y_he", 32);
582 strncpy(label, "Y_he", 32);
583 break;
584
585 case CS_GAS_MIX_AIR_HYDROGEN:
586 strncpy(name, "y_h2", 32);
587 strncpy(label, "Y_H2", 32);
588 break;
589
590 case CS_GAS_MIX_AIR_STEAM:
591 case CS_GAS_MIX_AIR_HELIUM_STEAM:
592 case CS_GAS_MIX_AIR_HYDROGEN_STEAM:
593 strncpy(name, "y_h2o_g", 32);
594 strncpy(label, "Y_H20_g", 32);
595 break;
596
597 case CS_GAS_MIX_HELIUM_AIR:
598 strncpy(name, "y_o2", 32);
599 strncpy(label, "Y_O2", 32);
600 break;
601
602 default:
603 add_deduced = false; /* should be added by user */
604
605 }
606
607 if (add_deduced) {
608
609 const int post_flag = CS_POST_ON_LOCATION | CS_POST_MONITOR;
610
611 f = cs_field_create(name,
612 field_type,
613 CS_MESH_LOCATION_CELLS,
614 1, /* dim */
615 true);
616
617 cs_field_set_key_int(f, keyvis, post_flag);
618 cs_field_set_key_int(f, keylog, 1);
619 cs_field_set_key_str(f, klbl, label);
620
621 _map_field(f);
622 _set_predefined_property(f);
623
624 }
625
626 /* Add molar mass property */
627
628 f = cs_field_create("mix_mol_mas",
629 field_type,
630 CS_MESH_LOCATION_CELLS,
631 1, /* dim */
632 false);
633
634 cs_field_set_key_int(f, keyvis, 0);
635 cs_field_set_key_int(f, keylog, 1);
636 }
637
638 /*----------------------------------------------------------------------------*/
639 /*!
640 * \brief Free array mapping gas mix species ids to field ids.
641 */
642 /*----------------------------------------------------------------------------*/
643
644 void
cs_gas_mix_finalize(void)645 cs_gas_mix_finalize(void)
646 {
647 BFT_FREE(_gas_mix.species_to_field_id);
648 _gas_mix.n_species = 0;
649 }
650
651 /*----------------------------------------------------------------------------*/
652
653 END_C_DECLS
654