8.3
general documentation
cs_advection_field.h
Go to the documentation of this file.
1#ifndef __CS_ADVECTION_FIELD_H__
2#define __CS_ADVECTION_FIELD_H__
3
4/*============================================================================
5 * Manage the definition/setting of advection fields
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_cdo_connect.h"
35#include "cs_cdo_local.h"
36#include "cs_cdo_quantities.h"
37#include "cs_field.h"
38#include "cs_macfb_builder.h"
39#include "cs_mesh_location.h"
40#include "cs_param_types.h"
41#include "cs_property.h"
42#include "cs_xdef.h"
43#include "cs_xdef_eval.h"
44
45/*----------------------------------------------------------------------------*/
46
48
49/*============================================================================
50 * Macro definitions
51 *============================================================================*/
52
62#define CS_ADVECTION_FIELD_POST_COURANT (1 << 0)
63
66/*============================================================================
67 * Type definitions
68 *============================================================================*/
69
71
122typedef enum {
123
124 /* Category of advection field */
125 /* --------------------------- */
126
128 CS_ADVECTION_FIELD_GWF = 1<<1, /* = 2 */
129 CS_ADVECTION_FIELD_USER = 1<<2, /* = 4 */
130
131 /* Type of advection field */
132 /* ----------------------- */
133
136
137 /* Optional */
138 /* -------- */
139
140 CS_ADVECTION_FIELD_STEADY = 1<<5, /* = 32 */
144
146
151typedef struct {
152
197 int id;
201
206
207 /* We assume that there is only one definition associated to an advection
208 field inside the computational domain */
209
211
212 /* Optional: Definition(s) for the boundary flux */
213
216 short int *bdy_def_ids;
217
219
220/*============================================================================
221 * Global variables
222 *============================================================================*/
223
224/*============================================================================
225 * Inline public function prototypes
226 *============================================================================*/
227
228/*----------------------------------------------------------------------------*/
235/*----------------------------------------------------------------------------*/
236
237static inline void
240{
241 if (adv == NULL)
242 return;
243
244 adv->status = status;
245}
246
247/*----------------------------------------------------------------------------*/
255/*----------------------------------------------------------------------------*/
256
257static inline bool
259{
260 if (adv == NULL)
261 return false;
262
264 return true;
265 else
266 return false;
267}
268
269/*----------------------------------------------------------------------------*/
278/*----------------------------------------------------------------------------*/
279
280static inline bool
282{
283 if (adv == NULL)
284 return false;
285
286 cs_flag_t state = adv->definition->state;
287
288 if (state & CS_FLAG_STATE_UNIFORM)
289 return true;
290 if (state & CS_FLAG_STATE_CELLWISE)
291 return true;
292 else
293 return false;
294}
295
296/*----------------------------------------------------------------------------*/
304/*----------------------------------------------------------------------------*/
305
306static inline const char *
308{
309 if (adv == NULL)
310 return NULL;
311
312 return adv->name;
313}
314
315/*----------------------------------------------------------------------------*/
324/*----------------------------------------------------------------------------*/
325
326static inline cs_xdef_type_t
328{
329 if (adv == NULL)
330 return CS_N_XDEF_TYPES;
331
332 return cs_xdef_get_type(adv->definition);
333}
334
335/*----------------------------------------------------------------------------*/
345/*----------------------------------------------------------------------------*/
346
347static inline cs_field_t *
350{
351 cs_field_t *f = NULL;
352 if (adv == NULL)
353 return f;
354
355 switch (ml_type) {
356
358 if (adv->cell_field_id > -1)
360 else
361 f = NULL;
362 break;
363
365 if (adv->int_field_id > -1)
367 else
368 f = NULL;
369 break;
370
372 if (adv->bdy_field_id > -1)
374 else
375 f = NULL;
376 break;
377
379 if (adv->vtx_field_id > -1)
381 else
382 f = NULL;
383 break;
384
385 default:
386 bft_error(__FILE__, __LINE__, 0,
387 " %s: Invalid mesh location type %d.\n"
388 " Stop retrieving the advection field.\n",
389 __func__, ml_type);
390 }
391
392 return f;
393}
394
395/*============================================================================
396 * Public function prototypes
397 *============================================================================*/
398
399/*----------------------------------------------------------------------------*/
406/*----------------------------------------------------------------------------*/
407
408void
410 const cs_cdo_connect_t *connect);
411
412/*----------------------------------------------------------------------------*/
418/*----------------------------------------------------------------------------*/
419
420int
422
423/*----------------------------------------------------------------------------*/
432/*----------------------------------------------------------------------------*/
433
435cs_advection_field_by_name(const char *name);
436
437/*----------------------------------------------------------------------------*/
446/*----------------------------------------------------------------------------*/
447
450
451/*----------------------------------------------------------------------------*/
459/*----------------------------------------------------------------------------*/
460
462cs_advection_field_add_user(const char *name);
463
464/*----------------------------------------------------------------------------*/
473/*----------------------------------------------------------------------------*/
474
476cs_advection_field_add(const char *name,
478
479/*----------------------------------------------------------------------------*/
485/*----------------------------------------------------------------------------*/
486
487void
489
490/*----------------------------------------------------------------------------*/
499/*----------------------------------------------------------------------------*/
500
501bool
503 const char *ref_name);
504
505/*----------------------------------------------------------------------------*/
509/*----------------------------------------------------------------------------*/
510
511void
513
514/*----------------------------------------------------------------------------*/
521/*----------------------------------------------------------------------------*/
522
523void
525 cs_flag_t post_flag);
526
527/*----------------------------------------------------------------------------*/
534/*----------------------------------------------------------------------------*/
535
536void
538 cs_real_t vector[3]);
539
540/*----------------------------------------------------------------------------*/
548/*----------------------------------------------------------------------------*/
549
550void
552 cs_analytic_func_t *func,
553 void *input);
554
555/*----------------------------------------------------------------------------*/
564/*----------------------------------------------------------------------------*/
565
566void
568 cs_flag_t dof_location,
569 cs_dof_func_t *func,
570 void *input);
571
572/*----------------------------------------------------------------------------*/
586/*----------------------------------------------------------------------------*/
587
588cs_xdef_t *
590 cs_flag_t val_location,
591 cs_real_t *array,
592 bool is_owner);
593
594/*----------------------------------------------------------------------------*/
601/*----------------------------------------------------------------------------*/
602
603void
606
607/*----------------------------------------------------------------------------*/
616/*----------------------------------------------------------------------------*/
617
618void
620 const char *zname,
621 cs_real_t normal_flux);
622
623/*----------------------------------------------------------------------------*/
633/*----------------------------------------------------------------------------*/
634
635void
637 const char *zname,
638 cs_analytic_func_t *func,
639 void *input);
640
641/*----------------------------------------------------------------------------*/
660/*----------------------------------------------------------------------------*/
661
662cs_xdef_t *
664 const char *zname,
665 cs_flag_t val_loc,
666 cs_real_t *array,
667 bool is_owner,
668 bool full_length);
669
670/*----------------------------------------------------------------------------*/
678/*----------------------------------------------------------------------------*/
679
680void
683
684/*----------------------------------------------------------------------------*/
689/*----------------------------------------------------------------------------*/
690
691void
693
694/*----------------------------------------------------------------------------*/
699/*----------------------------------------------------------------------------*/
700
701void
703
704/*----------------------------------------------------------------------------*/
712/*----------------------------------------------------------------------------*/
713
714void
716 const cs_adv_field_t *adv,
717 cs_nvec3_t *vect);
718
719/*----------------------------------------------------------------------------*/
730/*----------------------------------------------------------------------------*/
731
732void
734 const cs_cell_mesh_t *cm,
735 const cs_real_3_t xyz,
736 cs_real_t time_eval,
737 cs_nvec3_t *eval);
738
739/*----------------------------------------------------------------------------*/
748/*----------------------------------------------------------------------------*/
749
750void
752 cs_real_t time_eval,
753 cs_real_t *cell_values);
754
755/*----------------------------------------------------------------------------*/
763/*----------------------------------------------------------------------------*/
764
765void
767 cs_real_t time_eval,
768 cs_real_t *vtx_values);
769
770/*----------------------------------------------------------------------------*/
779/*----------------------------------------------------------------------------*/
780
781void
783 cs_real_t time_eval,
784 cs_real_t *flx_values);
785
786/*----------------------------------------------------------------------------*/
798/*----------------------------------------------------------------------------*/
799
800void
802 const cs_adv_field_t *adv,
803 short int f,
804 cs_real_t time_eval,
805 cs_real_t *fluxes);
806
807/*----------------------------------------------------------------------------*/
819/*----------------------------------------------------------------------------*/
820
823 const short int f,
824 const cs_cell_mesh_t *cm,
825 const cs_adv_field_t *adv);
826
827/*----------------------------------------------------------------------------*/
837/*----------------------------------------------------------------------------*/
838
839void
841 const cs_adv_field_t *adv,
842 cs_real_t time_eval,
843 cs_real_t *fluxes);
844
845/*----------------------------------------------------------------------------*/
855/*----------------------------------------------------------------------------*/
856
857void
859 const cs_adv_field_t *adv,
860 cs_real_t time_eval,
861 cs_real_t *fluxes);
862
863/*----------------------------------------------------------------------------*/
874/*----------------------------------------------------------------------------*/
875
876void
878 const cs_adv_field_t *adv,
879 cs_real_t time_eval,
880 cs_real_t *fluxes);
881
882/*----------------------------------------------------------------------------*/
890/*----------------------------------------------------------------------------*/
891
892void
894 bool cur2prev);
895
896/*----------------------------------------------------------------------------*/
905/*----------------------------------------------------------------------------*/
906
907void
909 const cs_property_t *diff,
910 cs_real_t t_eval,
911 cs_real_t peclet[]);
912
913/*----------------------------------------------------------------------------*/
921/*----------------------------------------------------------------------------*/
922
923void
925 cs_real_t dt_cur,
926 cs_real_t courant[]);
927
928/*----------------------------------------------------------------------------*/
938/*----------------------------------------------------------------------------*/
939
940cs_real_t *
942 cs_real_t t_eval);
943
944/*----------------------------------------------------------------------------*/
945
947
948#endif /* __CS_ADVECTION_FIELD_H__ */
void bft_error(const char *const file_name, const int line_num, const int sys_error_code, const char *const format,...)
Calls the error handler (set by bft_error_handler_set() or default).
Definition: bft_error.cpp:193
void cs_advection_field_macb_dface_flux(const cs_macfb_builder_t *macb, const cs_adv_field_t *adv, cs_real_t time_eval, cs_real_t *fluxes)
Compute the value of the flux of the advection field across the the dual faces of a cell MAC face-bas...
Definition: cs_advection_field.cpp:3066
cs_advection_field_status_bit_t
Bit values for specifying the definition/behavior of an advection field.
Definition: cs_advection_field.h:122
@ CS_ADVECTION_FIELD_STEADY
Definition: cs_advection_field.h:140
@ CS_ADVECTION_FIELD_DEFINE_AT_VERTICES
Definition: cs_advection_field.h:142
@ CS_ADVECTION_FIELD_GWF
Definition: cs_advection_field.h:128
@ CS_ADVECTION_FIELD_DEFINE_AT_BOUNDARY_FACES
Definition: cs_advection_field.h:143
@ CS_ADVECTION_FIELD_USER
Definition: cs_advection_field.h:129
@ CS_ADVECTION_FIELD_NAVSTO
Definition: cs_advection_field.h:127
@ CS_ADVECTION_FIELD_LEGACY_FV
Definition: cs_advection_field.h:141
@ CS_ADVECTION_FIELD_TYPE_VELOCITY_VECTOR
Definition: cs_advection_field.h:134
@ CS_ADVECTION_FIELD_TYPE_SCALAR_FLUX
Definition: cs_advection_field.h:135
void cs_advection_field_get_cell_vector(cs_lnum_t c_id, const cs_adv_field_t *adv, cs_nvec3_t *vect)
Compute the value of the advection field at the cell center.
Definition: cs_advection_field.cpp:1349
static bool cs_advection_field_is_cellwise(const cs_adv_field_t *adv)
returns true if the advection field is uniform in each cell otherwise false
Definition: cs_advection_field.h:281
void cs_advection_field_cw_boundary_f2v_flux(const cs_cell_mesh_t *cm, const cs_adv_field_t *adv, short int f, cs_real_t time_eval, cs_real_t *fluxes)
Compute the value of the normal flux of the advection field across the closure of the dual cell relat...
Definition: cs_advection_field.cpp:2384
cs_adv_field_t * cs_advection_field_add_user(const char *name)
Add and initialize a new user-defined advection field structure.
Definition: cs_advection_field.cpp:440
static cs_field_t * cs_advection_field_get_field(const cs_adv_field_t *adv, cs_mesh_location_type_t ml_type)
Get a cs_field_t structure related to an advection field and a mesh location.
Definition: cs_advection_field.h:348
void cs_advection_field_destroy_all(void)
Free all alllocated cs_adv_field_t structures and its related array.
Definition: cs_advection_field.cpp:546
bool cs_advection_field_check_name(const cs_adv_field_t *adv, const char *ref_name)
Check if the given advection field has the name ref_name.
Definition: cs_advection_field.cpp:588
void cs_advection_field_init_sharing(const cs_cdo_quantities_t *quant, const cs_cdo_connect_t *connect)
Set shared pointers to main domain members.
Definition: cs_advection_field.cpp:355
void cs_advection_field_update(cs_real_t t_eval, bool cur2prev)
For each cs_adv_field_t structures, update the values of the related field(s)
Definition: cs_advection_field.cpp:3159
cs_xdef_t * cs_advection_field_def_by_array(cs_adv_field_t *adv, cs_flag_t val_location, cs_real_t *array, bool is_owner)
Define a cs_adv_field_t structure thanks to an array of values. If an advanced usage of the definitio...
Definition: cs_advection_field.cpp:846
void cs_advection_field_in_cells(const cs_adv_field_t *adv, cs_real_t time_eval, cs_real_t *cell_values)
Compute the mean-value of the vector-valued field related to the advection field inside each cell.
Definition: cs_advection_field.cpp:1527
cs_adv_field_t * cs_advection_field_add(const char *name, cs_advection_field_status_t status)
Add and initialize a new advection field structure.
Definition: cs_advection_field.cpp:457
cs_real_t * cs_advection_field_divergence_at_vertices(const cs_adv_field_t *adv, cs_real_t t_eval)
Compute the divergence of the advection field at vertices Useful for CDO Vertex-based schemes.
Definition: cs_advection_field.cpp:3329
void cs_advection_field_def_by_dof_func(cs_adv_field_t *adv, cs_flag_t dof_location, cs_dof_func_t *func, void *input)
Define a cs_adv_field_t structure using a cs_dof_func_t.
Definition: cs_advection_field.cpp:804
void cs_advection_field_def_by_analytic(cs_adv_field_t *adv, cs_analytic_func_t *func, void *input)
Define a cs_adv_field_t structure thanks to an analytic function.
Definition: cs_advection_field.cpp:768
void cs_advection_field_set_postprocess(cs_adv_field_t *adv, cs_flag_t post_flag)
Set optional post-processings.
Definition: cs_advection_field.cpp:720
void cs_advection_get_courant(const cs_adv_field_t *adv, cs_real_t dt_cur, cs_real_t courant[])
Compute the Courant number in each cell.
Definition: cs_advection_field.cpp:3286
void cs_advection_field_cw_face_flux(const cs_cell_mesh_t *cm, const cs_adv_field_t *adv, cs_real_t time_eval, cs_real_t *fluxes)
Compute the value of the flux of the advection field across the the (primal) faces of a cell.
Definition: cs_advection_field.cpp:2606
void cs_advection_field_def_boundary_flux_by_field(cs_adv_field_t *adv, cs_field_t *field)
Define the value of the boundary normal flux for the given cs_adv_field_t structure using a field str...
Definition: cs_advection_field.cpp:1103
void cs_advection_get_peclet(const cs_adv_field_t *adv, const cs_property_t *diff, cs_real_t t_eval, cs_real_t peclet[])
Compute the Peclet number in each cell.
Definition: cs_advection_field.cpp:3238
static bool cs_advection_field_is_uniform(const cs_adv_field_t *adv)
returns true if the advection field is uniform, otherwise false
Definition: cs_advection_field.h:258
void cs_advection_field_def_boundary_flux_by_analytic(cs_adv_field_t *adv, const char *zname, cs_analytic_func_t *func, void *input)
Define the value of the boundary normal flux for the given cs_adv_field_t structure using an analytic...
Definition: cs_advection_field.cpp:982
cs_real_t cs_advection_field_cw_boundary_face_flux(const cs_real_t time_eval, const short int f, const cs_cell_mesh_t *cm, const cs_adv_field_t *adv)
Compute the value of the normal flux of the advection field across a boundary face f (cellwise versio...
Definition: cs_advection_field.cpp:2151
void cs_advection_field_cw_dface_flux(const cs_cell_mesh_t *cm, const cs_adv_field_t *adv, cs_real_t time_eval, cs_real_t *fluxes)
Compute the value of the flux of the advection field across the the dual faces of a cell.
Definition: cs_advection_field.cpp:2781
cs_flag_t cs_advection_field_status_t
Definition: cs_advection_field.h:70
static void cs_advection_field_set_status(cs_adv_field_t *adv, cs_advection_field_status_t status)
Set a new status for the given advection field structure.
Definition: cs_advection_field.h:238
void cs_advection_field_create_fields(void)
Create all needed cs_field_t structures related to an advection field.
Definition: cs_advection_field.cpp:1144
int cs_advection_field_get_n_fields(void)
Get the number of allocated cs_adv_field_t structures.
Definition: cs_advection_field.cpp:373
cs_xdef_t * cs_advection_field_def_boundary_flux_by_array(cs_adv_field_t *adv, const char *zname, cs_flag_t val_loc, cs_real_t *array, bool is_owner, bool full_length)
Define the value of the boundary normal flux for the given cs_adv_field_t structure using an array of...
Definition: cs_advection_field.cpp:1030
cs_adv_field_t * cs_advection_field_by_name(const char *name)
Search in the array of advection field structures which one has the name given in argument.
Definition: cs_advection_field.cpp:390
void cs_advection_field_def_by_value(cs_adv_field_t *adv, cs_real_t vector[3])
Define the value of a cs_adv_field_t structure.
static cs_xdef_type_t cs_advection_field_get_deftype(const cs_adv_field_t *adv)
Retrieve the type of definition used to set the current advection field structure.
Definition: cs_advection_field.h:327
static const char * cs_advection_field_get_name(const cs_adv_field_t *adv)
Retrieve the name of an advection field.
Definition: cs_advection_field.h:307
void cs_advection_field_finalize_setup(void)
Last stage of the definition of an advection field based on several definitions (i....
Definition: cs_advection_field.cpp:1259
void cs_advection_field_cw_eval_at_xyz(const cs_adv_field_t *adv, const cs_cell_mesh_t *cm, const cs_real_3_t xyz, cs_real_t time_eval, cs_nvec3_t *eval)
Compute the vector-valued interpolation of the advection field at a given location inside a cell.
Definition: cs_advection_field.cpp:1381
void cs_advection_field_at_vertices(const cs_adv_field_t *adv, cs_real_t time_eval, cs_real_t *vtx_values)
Compute the value of the advection field at vertices.
Definition: cs_advection_field.cpp:1719
void cs_advection_field_log_setup(void)
Print all setup information related to cs_adv_field_t structures.
Definition: cs_advection_field.cpp:614
void cs_advection_field_def_by_field(cs_adv_field_t *adv, cs_field_t *field)
Define a cs_adv_field_t structure thanks to a field structure.
Definition: cs_advection_field.cpp:907
void cs_advection_field_def_boundary_flux_by_value(cs_adv_field_t *adv, const char *zname, cs_real_t normal_flux)
Define the value of the boundary normal flux for the given cs_adv_field_t structure.
Definition: cs_advection_field.cpp:946
cs_adv_field_t * cs_advection_field_by_id(int id)
Search in the array of advection field structures which one has the id given in argument.
Definition: cs_advection_field.cpp:417
void cs_advection_field_across_boundary(const cs_adv_field_t *adv, cs_real_t time_eval, cs_real_t *flx_values)
Compute the value of the normal flux of the advection field across the boundary faces.
Definition: cs_advection_field.cpp:1834
#define restrict
Definition: cs_defs.h:145
#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_3_t[3]
vector of 3 floating-point values
Definition: cs_defs.h:359
#define END_C_DECLS
Definition: cs_defs.h:543
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:335
unsigned short int cs_flag_t
Definition: cs_defs.h:344
cs_field_t * cs_field_by_id(int id)
Return a pointer to a field based on its id.
Definition: cs_field.cpp:2465
#define CS_FLAG_STATE_CELLWISE
Definition: cs_flag.h:104
#define CS_FLAG_STATE_UNIFORM
Definition: cs_flag.h:103
cs_mesh_location_type_t
Definition: cs_mesh_location.h:60
@ CS_MESH_LOCATION_CELLS
Definition: cs_mesh_location.h:63
@ CS_MESH_LOCATION_INTERIOR_FACES
Definition: cs_mesh_location.h:64
@ CS_MESH_LOCATION_VERTICES
Definition: cs_mesh_location.h:66
@ CS_MESH_LOCATION_BOUNDARY_FACES
Definition: cs_mesh_location.h:65
void() cs_dof_func_t(cs_lnum_t n_elts, const cs_lnum_t *elt_ids, bool dense_output, void *input, cs_real_t *retval)
Generic function pointer for computing a quantity at predefined locations such as degrees of freedom ...
Definition: cs_param_types.h:154
void() cs_analytic_func_t(cs_real_t time, cs_lnum_t n_elts, const cs_lnum_t *elt_ids, const cs_real_t *coords, bool dense_output, void *input, cs_real_t *retval)
Generic function pointer for an evaluation relying on an analytic function.
Definition: cs_param_types.h:127
cs_xdef_type_t cs_xdef_get_type(const cs_xdef_t *d)
Retrieve the flag dedicated to the state.
Definition: cs_xdef.cpp:863
cs_xdef_type_t
Definition: cs_xdef.h:114
@ CS_N_XDEF_TYPES
Definition: cs_xdef.h:126
Definition: field.f90:27
Definition: cs_advection_field.h:151
cs_advection_field_status_t status
Definition: cs_advection_field.h:199
cs_xdef_t ** bdy_flux_defs
Definition: cs_advection_field.h:215
int cell_field_id
Definition: cs_advection_field.h:203
int bdy_field_id
Definition: cs_advection_field.h:204
int int_field_id
Definition: cs_advection_field.h:205
cs_xdef_t * definition
Definition: cs_advection_field.h:210
cs_flag_t post_flag
Definition: cs_advection_field.h:200
int n_bdy_flux_defs
Definition: cs_advection_field.h:214
char *restrict name
Definition: cs_advection_field.h:198
int vtx_field_id
Definition: cs_advection_field.h:202
short int * bdy_def_ids
Definition: cs_advection_field.h:216
int id
Definition: cs_advection_field.h:197
Definition: cs_cdo_connect.h:61
Definition: cs_cdo_quantities.h:139
Set of local quantities and connectivities related to a mesh cell.
Definition: cs_cdo_local.h:203
Field descriptor.
Definition: cs_field.h:131
Definition: cs_macfb_builder.h:54
Definition: cs_defs.h:400
Structure associated to the definition of a property relying on the cs_xdef_t structure.
Structure storing medata for defining a quantity in a very flexible way.
Definition: cs_xdef.h:160
cs_flag_t state
Definition: cs_xdef.h:199