My Project
programmer's 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-2019 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_mesh_location.h"
38#include "cs_param.h"
39#include "cs_property.h"
40#include "cs_xdef.h"
41#include "cs_xdef_eval.h"
42
43/*----------------------------------------------------------------------------*/
44
46
47/*============================================================================
48 * Macro definitions
49 *============================================================================*/
50
59#define CS_ADVECTION_FIELD_POST_COURANT (1 << 0)
61#define CS_ADVECTION_FIELD_STEADY (1 << 1) /* steady-state field */
62
65/*============================================================================
66 * Type definitions
67 *============================================================================*/
68
96
124
144
149typedef struct {
150
197 int id;
202
207
208 /* We assume that there is only one definition associated to an advection
209 field inside the computational domain */
211
212 /* Optional: Definition(s) for the boundary flux */
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->type = type;
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
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/*----------------------------------------------------------------------------*/
565/*----------------------------------------------------------------------------*/
566
567void
569 cs_flag_t loc,
570 cs_real_t *array,
571 bool is_owner,
572 cs_lnum_t *index);
573
574/*----------------------------------------------------------------------------*/
581/*----------------------------------------------------------------------------*/
582
583void
586
587/*----------------------------------------------------------------------------*/
596/*----------------------------------------------------------------------------*/
597
598void
600 const char *zname,
601 cs_real_t normal_flux);
602
603/*----------------------------------------------------------------------------*/
613/*----------------------------------------------------------------------------*/
614
615void
617 const char *zname,
618 cs_analytic_func_t *func,
619 void *input);
620
621/*----------------------------------------------------------------------------*/
634/*----------------------------------------------------------------------------*/
635
636void
638 const char *zname,
639 cs_flag_t loc,
640 cs_real_t *array,
641 bool is_owner,
642 cs_lnum_t *index);
643
644/*----------------------------------------------------------------------------*/
649/*----------------------------------------------------------------------------*/
650
651void
653
654/*----------------------------------------------------------------------------*/
659/*----------------------------------------------------------------------------*/
660
661void
663
664/*----------------------------------------------------------------------------*/
672/*----------------------------------------------------------------------------*/
673
674void
676 const cs_adv_field_t *adv,
677 cs_nvec3_t *vect);
678
679/*----------------------------------------------------------------------------*/
690/*----------------------------------------------------------------------------*/
691
692void
694 const cs_cell_mesh_t *cm,
695 const cs_real_3_t xyz,
696 cs_real_t time_eval,
697 cs_nvec3_t *eval);
698
699/*----------------------------------------------------------------------------*/
707/*----------------------------------------------------------------------------*/
708
709void
711 cs_real_t time_eval,
712 cs_real_t *cell_values);
713
714/*----------------------------------------------------------------------------*/
722/*----------------------------------------------------------------------------*/
723
724void
726 cs_real_t time_eval,
727 cs_real_t *vtx_values);
728
729/*----------------------------------------------------------------------------*/
738/*----------------------------------------------------------------------------*/
739
740void
742 cs_real_t time_eval,
743 cs_real_t *flx_values);
744
745/*----------------------------------------------------------------------------*/
757/*----------------------------------------------------------------------------*/
758
759void
761 const cs_adv_field_t *adv,
762 short int f,
763 cs_real_t time_eval,
764 cs_real_t *fluxes);
765
766/*----------------------------------------------------------------------------*/
778/*----------------------------------------------------------------------------*/
779
782 const short int f,
783 const cs_cell_mesh_t *cm,
784 const cs_adv_field_t *adv);
785
786/*----------------------------------------------------------------------------*/
796/*----------------------------------------------------------------------------*/
797
798void
800 const cs_adv_field_t *adv,
801 cs_real_t time_eval,
802 cs_real_t *fluxes);
803
804/*----------------------------------------------------------------------------*/
814/*----------------------------------------------------------------------------*/
815
816void
818 const cs_adv_field_t *adv,
819 cs_real_t time_eval,
820 cs_real_t *fluxes);
821
822/*----------------------------------------------------------------------------*/
830/*----------------------------------------------------------------------------*/
831
832void
834 bool cur2prev);
835
836/*----------------------------------------------------------------------------*/
845/*----------------------------------------------------------------------------*/
846
847void
849 const cs_property_t *diff,
850 cs_real_t t_eval,
851 cs_real_t peclet[]);
852
853/*----------------------------------------------------------------------------*/
861/*----------------------------------------------------------------------------*/
862
863void
865 cs_real_t dt_cur,
866 cs_real_t courant[]);
867
868/*----------------------------------------------------------------------------*/
878/*----------------------------------------------------------------------------*/
879
880cs_real_t *
882 cs_real_t t_eval);
883
884/*----------------------------------------------------------------------------*/
885
887
888#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
static void cs_advection_field_set_type(cs_adv_field_t *adv, cs_advection_field_type_t type)
Set the type of advection field for the given structure.
Definition cs_advection_field.h:237
void cs_advection_field_set_option(cs_adv_field_t *adv, cs_advection_field_key_t key)
Set optional parameters related to a cs_adv_field_t structure.
Definition cs_advection_field.c:568
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:1037
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:1923
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:310
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_set_shared_pointers(const cs_cdo_quantities_t *quant, const cs_cdo_connect_t *connect)
Set shared pointers to main domain members.
Definition cs_advection_field.c:225
void cs_advection_field_destroy_all(void)
Free all alllocated cs_adv_field_t structures and its related array.
Definition cs_advection_field.c:392
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:432
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:2633
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 advection field inside each cell.
Definition cs_advection_field.c:1133
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:327
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:2801
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:637
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:2756
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:2156
void cs_advection_field_def_boundary_flux_by_array(cs_adv_field_t *adv, const char *zname, cs_flag_t loc, cs_real_t *array, bool is_owner, cs_lnum_t *index)
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:821
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:2708
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:779
cs_advection_field_status_t
Type of advection field.
Definition cs_advection_field.h:115
@ CS_ADVECTION_FIELD_GWF
Definition cs_advection_field.h:119
@ CS_N_ADVECTION_FIELD_STATUS
Definition cs_advection_field.h:121
@ CS_ADVECTION_FIELD_USER
Definition cs_advection_field.h:120
@ CS_ADVECTION_FIELD_NAVSTO
Definition cs_advection_field.h:117
@ CS_ADVECTION_FIELD_LEGACY_NAVSTO
Definition cs_advection_field.h:118
void cs_advection_field_def_by_array(cs_adv_field_t *adv, cs_flag_t loc, cs_real_t *array, bool is_owner, cs_lnum_t *index)
Define a cs_adv_field_t structure thanks to an array of values.
Definition cs_advection_field.c:671
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:1711
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:2343
cs_advection_field_type_t
Status of the advection field. The advection field stands for what.
Definition cs_advection_field.h:137
@ CS_ADVECTION_FIELD_TYPE_VELOCITY
Definition cs_advection_field.h:139
@ CS_N_ADVECTION_FIELD_TYPES
Definition cs_advection_field.h:141
@ CS_ADVECTION_FIELD_TYPE_FLUX
Definition cs_advection_field.h:140
void cs_advection_field_create_fields(void)
Create all needed cs_field_t structures related to an advection field.
Definition cs_advection_field.c:864
int cs_advection_field_get_n_fields(void)
Get the number of allocated cs_adv_field_t structures.
Definition cs_advection_field.c:242
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:259
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.
cs_advection_field_key_t
List of available keys for setting an advection field.
Definition cs_advection_field.h:87
@ CS_ADVKEY_N_KEYS
Definition cs_advection_field.h:93
@ CS_ADVKEY_POST_COURANT
Definition cs_advection_field.h:91
@ CS_ADVKEY_STATE_STEADY
Definition cs_advection_field.h:92
@ CS_ADVKEY_DEFINE_AT_VERTICES
Definition cs_advection_field.h:89
@ CS_ADVKEY_DEFINE_AT_BOUNDARY_FACES
Definition cs_advection_field.h:90
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:972
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 value of the advection field at a specific location inside a cell.
Definition cs_advection_field.c:1068
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:1281
void cs_advection_field_log_setup(void)
Print all setup information related to cs_adv_field_t structures.
Definition cs_advection_field.c:458
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 an array of values.
Definition cs_advection_field.c:707
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:743
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:287
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:1452
#define restrict
Definition cs_defs.h:127
#define BEGIN_C_DECLS
Definition cs_defs.h:467
double cs_real_t
Floating-point value.
Definition cs_defs.h:302
cs_real_t cs_real_3_t[3]
vector of 3 floating-point values
Definition cs_defs.h:315
#define END_C_DECLS
Definition cs_defs.h:468
int cs_lnum_t
local mesh entity id
Definition cs_defs.h:298
unsigned short int cs_flag_t
Definition cs_defs.h:304
cs_field_t * cs_field_by_id(int id)
Return a pointer to a field based on its id.
Definition cs_field.c:2307
#define CS_FLAG_STATE_CELLWISE
Definition cs_flag.h:75
#define CS_FLAG_STATE_UNIFORM
Definition cs_flag.h:74
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 compact, void *input, cs_real_t *retval)
Generic function pointer for an analytic function elt_ids is optional. If not NULL,...
Definition cs_param.h:66
cs_xdef_type_t cs_xdef_get_type(const cs_xdef_t *d)
Retrieve the flag dedicated to the state.
Definition cs_xdef.c:604
cs_xdef_type_t
Definition cs_xdef.h:82
@ CS_N_XDEF_TYPES
Definition cs_xdef.h:92
static int input(void)
Definition field.f90:27
Definition cs_advection_field.h:149
cs_advection_field_status_t status
Definition cs_advection_field.h:199
cs_xdef_t ** bdy_flux_defs
Definition cs_advection_field.h:214
int cell_field_id
Definition cs_advection_field.h:204
int bdy_field_id
Definition cs_advection_field.h:205
int int_field_id
Definition cs_advection_field.h:206
cs_xdef_t * definition
Definition cs_advection_field.h:210
cs_flag_t flag
Definition cs_advection_field.h:201
cs_advection_field_type_t type
Definition cs_advection_field.h:200
int n_bdy_flux_defs
Definition cs_advection_field.h:213
char *restrict name
Definition cs_advection_field.h:198
int vtx_field_id
Definition cs_advection_field.h:203
short int * bdy_def_ids
Definition cs_advection_field.h:215
int id
Definition cs_advection_field.h:197
Definition cs_cdo_connect.h:74
Definition cs_cdo_quantities.h:94
Set of local quantities and connectivities related to a mesh cell This is a key structure for all cel...
Definition cs_cdo_local.h:146
Field descriptor.
Definition cs_field.h:124
Definition cs_defs.h:343
Definition cs_property.h:104
Structure storing medata for defining a quantity in a very flexible way.
Definition cs_xdef.h:126
cs_flag_t state
Definition cs_xdef.h:164