8.0
general documentation
Loading...
Searching...
No Matches
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-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_cdo_connect.h"
35#include "cs_cdo_local.h"
36#include "cs_cdo_quantities.h"
37#include "cs_field.h"
38#include "cs_mesh_location.h"
39#include "cs_param_types.h"
40#include "cs_property.h"
41#include "cs_xdef.h"
42#include "cs_xdef_eval.h"
43
44/*----------------------------------------------------------------------------*/
45
47
48/*============================================================================
49 * Macro definitions
50 *============================================================================*/
51
58
60
61#define CS_ADVECTION_FIELD_POST_COURANT (1 << 0)
62
64
65/*============================================================================
66 * Type definitions
67 *============================================================================*/
68
70
120
121typedef enum {
122
123 /* Category of advection field */
124 /* --------------------------- */
125
127 CS_ADVECTION_FIELD_GWF = 1<<1, /* = 2 */
128 CS_ADVECTION_FIELD_USER = 1<<2, /* = 4 */
129
130 /* Type of advection field */
131 /* ----------------------- */
132
135
136 /* Optional */
137 /* -------- */
138
139 CS_ADVECTION_FIELD_STEADY = 1<<5, /* = 32 */
143
145
149
150typedef struct {
151
195
196 int id;
200
205
206 /* We assume that there is only one definition associated to an advection
207 field inside the computational domain */
208
210
211 /* Optional: Definition(s) for the boundary flux */
212
215 short int *bdy_def_ids;
216
218
219/*============================================================================
220 * Global variables
221 *============================================================================*/
222
223/*============================================================================
224 * Inline public function prototypes
225 *============================================================================*/
226
227/*----------------------------------------------------------------------------*/
234/*----------------------------------------------------------------------------*/
235
236static inline void
239{
240 if (adv == NULL)
241 return;
242
243 adv->status = status;
244}
245
246/*----------------------------------------------------------------------------*/
254/*----------------------------------------------------------------------------*/
255
256static inline bool
258{
259 if (adv == NULL)
260 return false;
261
263 return true;
264 else
265 return false;
266}
267
268/*----------------------------------------------------------------------------*/
277/*----------------------------------------------------------------------------*/
278
279static inline bool
281{
282 if (adv == NULL)
283 return false;
284
285 cs_flag_t state = adv->definition->state;
286
287 if (state & CS_FLAG_STATE_UNIFORM)
288 return true;
289 if (state & CS_FLAG_STATE_CELLWISE)
290 return true;
291 else
292 return false;
293}
294
295/*----------------------------------------------------------------------------*/
303/*----------------------------------------------------------------------------*/
304
305static inline const char *
307{
308 if (adv == NULL)
309 return NULL;
310
311 return adv->name;
312}
313
314/*----------------------------------------------------------------------------*/
323/*----------------------------------------------------------------------------*/
324
325static inline cs_xdef_type_t
327{
328 if (adv == NULL)
329 return CS_N_XDEF_TYPES;
330
331 return cs_xdef_get_type(adv->definition);
332}
333
334/*----------------------------------------------------------------------------*/
344/*----------------------------------------------------------------------------*/
345
346static inline cs_field_t *
349{
350 cs_field_t *f = NULL;
351 if (adv == NULL)
352 return f;
353
354 switch (ml_type) {
355
357 if (adv->cell_field_id > -1)
359 else
360 f = NULL;
361 break;
362
364 if (adv->int_field_id > -1)
366 else
367 f = NULL;
368 break;
369
371 if (adv->bdy_field_id > -1)
373 else
374 f = NULL;
375 break;
376
378 if (adv->vtx_field_id > -1)
380 else
381 f = NULL;
382 break;
383
384 default:
385 bft_error(__FILE__, __LINE__, 0,
386 " %s: Invalid mesh location type %d.\n"
387 " Stop retrieving the advection field.\n",
388 __func__, ml_type);
389 }
390
391 return f;
392}
393
394/*============================================================================
395 * Public function prototypes
396 *============================================================================*/
397
398/*----------------------------------------------------------------------------*/
405/*----------------------------------------------------------------------------*/
406
407void
409 const cs_cdo_connect_t *connect);
410
411/*----------------------------------------------------------------------------*/
417/*----------------------------------------------------------------------------*/
418
419int
421
422/*----------------------------------------------------------------------------*/
431/*----------------------------------------------------------------------------*/
432
434cs_advection_field_by_name(const char *name);
435
436/*----------------------------------------------------------------------------*/
445/*----------------------------------------------------------------------------*/
446
449
450/*----------------------------------------------------------------------------*/
458/*----------------------------------------------------------------------------*/
459
461cs_advection_field_add_user(const char *name);
462
463/*----------------------------------------------------------------------------*/
472/*----------------------------------------------------------------------------*/
473
475cs_advection_field_add(const char *name,
477
478/*----------------------------------------------------------------------------*/
484/*----------------------------------------------------------------------------*/
485
486void
488
489/*----------------------------------------------------------------------------*/
498/*----------------------------------------------------------------------------*/
499
500bool
502 const char *ref_name);
503
504/*----------------------------------------------------------------------------*/
508/*----------------------------------------------------------------------------*/
509
510void
512
513/*----------------------------------------------------------------------------*/
520/*----------------------------------------------------------------------------*/
521
522void
524 cs_flag_t post_flag);
525
526/*----------------------------------------------------------------------------*/
533/*----------------------------------------------------------------------------*/
534
535void
537 cs_real_t vector[3]);
538
539/*----------------------------------------------------------------------------*/
547/*----------------------------------------------------------------------------*/
548
549void
551 cs_analytic_func_t *func,
552 void *input);
553
554/*----------------------------------------------------------------------------*/
563/*----------------------------------------------------------------------------*/
564
565void
567 cs_flag_t loc_flag,
568 cs_dof_func_t *func,
569 void *input);
570
571/*----------------------------------------------------------------------------*/
585/*----------------------------------------------------------------------------*/
586
587cs_xdef_t *
589 cs_flag_t val_location,
590 cs_real_t *array,
591 bool is_owner);
592
593/*----------------------------------------------------------------------------*/
600/*----------------------------------------------------------------------------*/
601
602void
605
606/*----------------------------------------------------------------------------*/
615/*----------------------------------------------------------------------------*/
616
617void
619 const char *zname,
620 cs_real_t normal_flux);
621
622/*----------------------------------------------------------------------------*/
632/*----------------------------------------------------------------------------*/
633
634void
636 const char *zname,
637 cs_analytic_func_t *func,
638 void *input);
639
640/*----------------------------------------------------------------------------*/
659/*----------------------------------------------------------------------------*/
660
661cs_xdef_t *
663 const char *zname,
664 cs_flag_t val_loc,
665 cs_real_t *array,
666 bool is_owner,
667 bool full_length);
668
669/*----------------------------------------------------------------------------*/
677/*----------------------------------------------------------------------------*/
678
679void
682
683/*----------------------------------------------------------------------------*/
688/*----------------------------------------------------------------------------*/
689
690void
692
693/*----------------------------------------------------------------------------*/
698/*----------------------------------------------------------------------------*/
699
700void
702
703/*----------------------------------------------------------------------------*/
711/*----------------------------------------------------------------------------*/
712
713void
715 const cs_adv_field_t *adv,
716 cs_nvec3_t *vect);
717
718/*----------------------------------------------------------------------------*/
729/*----------------------------------------------------------------------------*/
730
731void
733 const cs_cell_mesh_t *cm,
734 const cs_real_3_t xyz,
735 cs_real_t time_eval,
736 cs_nvec3_t *eval);
737
738/*----------------------------------------------------------------------------*/
747/*----------------------------------------------------------------------------*/
748
749void
751 cs_real_t time_eval,
752 cs_real_t *cell_values);
753
754/*----------------------------------------------------------------------------*/
762/*----------------------------------------------------------------------------*/
763
764void
766 cs_real_t time_eval,
767 cs_real_t *vtx_values);
768
769/*----------------------------------------------------------------------------*/
778/*----------------------------------------------------------------------------*/
779
780void
782 cs_real_t time_eval,
783 cs_real_t *flx_values);
784
785/*----------------------------------------------------------------------------*/
797/*----------------------------------------------------------------------------*/
798
799void
801 const cs_adv_field_t *adv,
802 short int f,
803 cs_real_t time_eval,
804 cs_real_t *fluxes);
805
806/*----------------------------------------------------------------------------*/
818/*----------------------------------------------------------------------------*/
819
822 const short int f,
823 const cs_cell_mesh_t *cm,
824 const cs_adv_field_t *adv);
825
826/*----------------------------------------------------------------------------*/
836/*----------------------------------------------------------------------------*/
837
838void
840 const cs_adv_field_t *adv,
841 cs_real_t time_eval,
842 cs_real_t *fluxes);
843
844/*----------------------------------------------------------------------------*/
854/*----------------------------------------------------------------------------*/
855
856void
858 const cs_adv_field_t *adv,
859 cs_real_t time_eval,
860 cs_real_t *fluxes);
861
862/*----------------------------------------------------------------------------*/
870/*----------------------------------------------------------------------------*/
871
872void
874 bool cur2prev);
875
876/*----------------------------------------------------------------------------*/
885/*----------------------------------------------------------------------------*/
886
887void
889 const cs_property_t *diff,
890 cs_real_t t_eval,
891 cs_real_t peclet[]);
892
893/*----------------------------------------------------------------------------*/
901/*----------------------------------------------------------------------------*/
902
903void
905 cs_real_t dt_cur,
906 cs_real_t courant[]);
907
908/*----------------------------------------------------------------------------*/
918/*----------------------------------------------------------------------------*/
919
920cs_real_t *
922 cs_real_t t_eval);
923
924/*----------------------------------------------------------------------------*/
925
927
928#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.c:193
cs_advection_field_status_bit_t
Bit values for specifying the definition/behavior of an advection field.
Definition cs_advection_field.h:121
@ CS_ADVECTION_FIELD_STEADY
Definition cs_advection_field.h:139
@ CS_ADVECTION_FIELD_LEGACY_FV
Definition cs_advection_field.h:140
@ CS_ADVECTION_FIELD_TYPE_SCALAR_FLUX
Definition cs_advection_field.h:134
@ CS_ADVECTION_FIELD_DEFINE_AT_BOUNDARY_FACES
Definition cs_advection_field.h:142
@ CS_ADVECTION_FIELD_NAVSTO
Definition cs_advection_field.h:126
@ CS_ADVECTION_FIELD_TYPE_VELOCITY_VECTOR
Definition cs_advection_field.h:133
@ CS_ADVECTION_FIELD_GWF
Definition cs_advection_field.h:127
@ CS_ADVECTION_FIELD_USER
Definition cs_advection_field.h:128
@ CS_ADVECTION_FIELD_DEFINE_AT_VERTICES
Definition cs_advection_field.h:141
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.c:1358
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:280
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.c:2419
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.c:449
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:347
void cs_advection_field_destroy_all(void)
Free all alllocated cs_adv_field_t structures and its related array.
Definition cs_advection_field.c:555
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.c:595
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.c:363
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.c:3143
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.c:849
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.c:1545
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.c:466
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.c:3317
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.c:772
void cs_advection_field_set_postprocess(cs_adv_field_t *adv, cs_flag_t post_flag)
Set optional post-processings.
Definition cs_advection_field.c:723
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.c:3273
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.c:2648
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.c:1109
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.c:3225
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:257
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.c:985
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.c:2169
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.c:2839
cs_flag_t cs_advection_field_status_t
Definition cs_advection_field.h:69
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:237
void cs_advection_field_def_by_dof_func(cs_adv_field_t *adv, cs_flag_t loc_flag, cs_dof_func_t *func, void *input)
Define a cs_adv_field_t structure using a cs_dof_func_t.
Definition cs_advection_field.c:807
void cs_advection_field_create_fields(void)
Create all needed cs_field_t structures related to an advection field.
Definition cs_advection_field.c:1150
int cs_advection_field_get_n_fields(void)
Get the number of allocated cs_adv_field_t structures.
Definition cs_advection_field.c:381
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.c:1035
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.c:398
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:326
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:306
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.c:1266
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.c:1390
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.c:1741
void cs_advection_field_log_setup(void)
Print all setup information related to cs_adv_field_t structures.
Definition cs_advection_field.c:621
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.c:909
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.c:949
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.c:426
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.c:1859
#define restrict
Definition cs_defs.h:139
#define BEGIN_C_DECLS
Definition cs_defs.h:509
double cs_real_t
Floating-point value.
Definition cs_defs.h:319
#define END_C_DECLS
Definition cs_defs.h:510
int cs_lnum_t
local mesh entity id
Definition cs_defs.h:313
unsigned short int cs_flag_t
Definition cs_defs.h:321
cs_real_t cs_real_3_t[3]
vector of 3 floating-point values
Definition cs_defs.h:332
cs_field_t * cs_field_by_id(int id)
Return a pointer to a field based on its id.
Definition cs_field.c:2316
#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_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
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
cs_xdef_type_t cs_xdef_get_type(const cs_xdef_t *d)
Retrieve the flag dedicated to the state.
Definition cs_xdef.c:845
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:150
cs_xdef_t ** bdy_flux_defs
Definition cs_advection_field.h:214
int n_bdy_flux_defs
Definition cs_advection_field.h:213
int bdy_field_id
Definition cs_advection_field.h:203
int vtx_field_id
Definition cs_advection_field.h:201
int id
Definition cs_advection_field.h:196
cs_flag_t post_flag
Definition cs_advection_field.h:199
short int * bdy_def_ids
Definition cs_advection_field.h:215
cs_advection_field_status_t status
Definition cs_advection_field.h:198
int cell_field_id
Definition cs_advection_field.h:202
cs_xdef_t * definition
Definition cs_advection_field.h:209
int int_field_id
Definition cs_advection_field.h:204
char *restrict name
Definition cs_advection_field.h:197
Definition cs_cdo_connect.h:61
Definition cs_cdo_quantities.h:137
Set of local quantities and connectivities related to a mesh cell.
Definition cs_cdo_local.h:203
Field descriptor.
Definition cs_field.h:130
Definition cs_defs.h:367
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