8.3
general documentation
cs_cdofb_navsto.h
Go to the documentation of this file.
1#ifndef __CS_CDOFB_NAVSTO_H__
2#define __CS_CDOFB_NAVSTO_H__
3
4/*============================================================================
5 * Functions shared among all face-based schemes for the discretization of the
6 * Navier--Stokes system
7 *============================================================================*/
8
9/*
10 This file is part of code_saturne, a general-purpose CFD tool.
11
12 Copyright (C) 1998-2024 EDF S.A.
13
14 This program is free software; you can redistribute it and/or modify it under
15 the terms of the GNU General Public License as published by the Free Software
16 Foundation; either version 2 of the License, or (at your option) any later
17 version.
18
19 This program is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
21 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
22 details.
23
24 You should have received a copy of the GNU General Public License along with
25 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
26 Street, Fifth Floor, Boston, MA 02110-1301, USA.
27*/
28
29/*----------------------------------------------------------------------------*/
30
31#include "cs_defs.h"
32
33/*----------------------------------------------------------------------------
34 * Standard C library headers
35 *----------------------------------------------------------------------------*/
36
37/*----------------------------------------------------------------------------
38 * Local headers
39 *----------------------------------------------------------------------------*/
40
41#include "cs_base.h"
42#include "cs_cdo_connect.h"
43#include "cs_cdo_quantities.h"
44#include "cs_field.h"
45#include "cs_iter_algo.h"
46#include "cs_math.h"
47#include "cs_matrix.h"
48#include "cs_mesh.h"
49#include "cs_navsto_param.h"
50#include "cs_time_plot.h"
51#include "cs_time_step.h"
52#include "cs_sdm.h"
53
54/*----------------------------------------------------------------------------*/
55
57
58/*============================================================================
59 * Macro definitions
60 *============================================================================*/
61
62/*============================================================================
63 * Type definitions
64 *============================================================================*/
65
71typedef enum {
72
83
94
96
97
107typedef struct {
108
116
129
140
159 cs_boundary_type_t *bf_type; /* Size: n_fc */
160 cs_real_t *pressure_bc_val; /* Size: n_fc */
161
167
168/*----------------------------------------------------------------------------*/
181/*----------------------------------------------------------------------------*/
182
183typedef void
185 const cs_cell_mesh_t *cm,
186 const cs_cdofb_navsto_builder_t *nsb,
187 cs_cell_sys_t *csys);
188
189/*============================================================================
190 * Public function prototypes
191 *============================================================================*/
192
193/*----------------------------------------------------------------------------*/
199/*----------------------------------------------------------------------------*/
200
201void
203
204/*----------------------------------------------------------------------------*/
213/*----------------------------------------------------------------------------*/
214
217 const cs_cdo_connect_t *connect);
218
219/*----------------------------------------------------------------------------*/
225/*----------------------------------------------------------------------------*/
226
227void
229
230/*----------------------------------------------------------------------------*/
242/*----------------------------------------------------------------------------*/
243
244void
246 const cs_navsto_param_t *nsp,
247 const cs_cell_mesh_t *cm,
248 const cs_cell_sys_t *csys,
249 const cs_cdo_bc_face_t *pr_bc,
250 const cs_boundary_type_t *bf_type,
252
253/*----------------------------------------------------------------------------*/
265/*----------------------------------------------------------------------------*/
266
267void
269 const cs_cdo_quantities_t *quant,
270 const cs_real_t *face_vel,
271 cs_real_t *mass_flux);
272
273/*----------------------------------------------------------------------------*/
286/*----------------------------------------------------------------------------*/
287
288double
290 const cs_cdo_quantities_t *quant,
291 const cs_adjacency_t *c2f,
292 const cs_real_t *f_vals);
293
294/*----------------------------------------------------------------------------*/
306/*----------------------------------------------------------------------------*/
307
308void
310 const cs_cdo_connect_t *connect,
311 const cs_cdo_quantities_t *quant,
312 const cs_time_step_t *ts,
313 const cs_navsto_param_t *nsp,
314 const cs_real_t *p_cell,
315 cs_real_t *p_face);
316
317/*----------------------------------------------------------------------------*/
327/*----------------------------------------------------------------------------*/
328
329void
330cs_cdofb_navsto_add_grad_div(short int n_fc,
331 const cs_real_t zeta,
332 const cs_real_t div[],
333 cs_sdm_t *mat);
334
335/*----------------------------------------------------------------------------*/
344/*----------------------------------------------------------------------------*/
345
346void
348 const cs_cdo_quantities_t *quant,
349 const cs_time_step_t *ts,
350 cs_field_t *pr);
351
352/*----------------------------------------------------------------------------*/
362/*----------------------------------------------------------------------------*/
363
364void
366 const cs_cdo_connect_t *connect,
367 const cs_time_step_t *ts,
368 cs_real_t *pr_f);
369
370/*----------------------------------------------------------------------------*/
379/*----------------------------------------------------------------------------*/
380
381void
383 const cs_cdo_quantities_t *quant,
384 cs_real_t values[]);
385
386/*----------------------------------------------------------------------------*/
394/*----------------------------------------------------------------------------*/
395
396void
398 cs_real_t values[]);
399
400/*----------------------------------------------------------------------------*/
429/*----------------------------------------------------------------------------*/
430
431void
433 const cs_mesh_t *mesh,
434 const cs_cdo_quantities_t *quant,
435 const cs_cdo_connect_t *connect,
436 const cs_time_step_t *ts,
437 cs_time_plot_t *time_plotter,
438 const cs_adv_field_t *adv_field,
439 const cs_real_t *mass_flux,
440 const cs_real_t *p_cell,
441 const cs_real_t *u_cell,
442 const cs_real_t *u_face);
443
444/*----------------------------------------------------------------------------*/
459/*----------------------------------------------------------------------------*/
460
461void
463 const cs_equation_param_t *eqp,
464 const cs_cell_mesh_t *cm,
465 const cs_property_data_t *pty,
467 cs_cell_sys_t *csys);
468
469/*----------------------------------------------------------------------------*/
486/*----------------------------------------------------------------------------*/
487
488void
490 const cs_equation_param_t *eqp,
491 const cs_cell_mesh_t *cm,
492 const cs_property_data_t *pty,
494 cs_cell_sys_t *csys);
495
496/*----------------------------------------------------------------------------*/
513/*----------------------------------------------------------------------------*/
514
515void
517 const cs_equation_param_t *eqp,
518 const cs_cell_mesh_t *cm,
519 const cs_property_data_t *pty,
521 cs_cell_sys_t *csys);
522
523/*----------------------------------------------------------------------------*/
540/*----------------------------------------------------------------------------*/
541
542void
544 const cs_equation_param_t *eqp,
545 const cs_cell_mesh_t *cm,
546 const cs_property_data_t *pty,
548 cs_cell_sys_t *csys);
549
550/*----------------------------------------------------------------------------*/
566/*----------------------------------------------------------------------------*/
567
568void
569cs_cdofb_symmetry(short int fb,
570 const cs_equation_param_t *eqp,
571 const cs_cell_mesh_t *cm,
572 const cs_property_data_t *pty,
574 cs_cell_sys_t *csys);
575
576/*----------------------------------------------------------------------------*/
589/*----------------------------------------------------------------------------*/
590
591void
592cs_cdofb_fixed_wall(short int fb,
593 const cs_equation_param_t *eqp,
594 const cs_cell_mesh_t *cm,
595 const cs_property_data_t *pty,
597 cs_cell_sys_t *csys);
598
599/*----------------------------------------------------------------------------*/
612/*----------------------------------------------------------------------------*/
613
616 const cs_real_t *pre_iterate,
617 cs_real_t *cur_iterate,
619
620/*----------------------------------------------------------------------------*/
629/*----------------------------------------------------------------------------*/
630
631void
633 cs_cdofb_navsto_source_t **p_func);
634
635/*----------------------------------------------------------------------------*/
647/*----------------------------------------------------------------------------*/
648
649void
651 const cs_cell_mesh_t *cm,
652 const cs_cdofb_navsto_builder_t *nsb,
653 cs_cell_sys_t *csys);
654
655/*----------------------------------------------------------------------------*/
668/*----------------------------------------------------------------------------*/
669
670void
672 const cs_cell_mesh_t *cm,
673 const cs_cdofb_navsto_builder_t *nsb,
674 cs_cell_sys_t *csys);
675
676/*----------------------------------------------------------------------------*/
689/*----------------------------------------------------------------------------*/
690
691void
693 const cs_cell_mesh_t *cm,
694 const cs_cdofb_navsto_builder_t *nsb,
695 cs_cell_sys_t *csys);
696
697/*----------------------------------------------------------------------------*/
710/*----------------------------------------------------------------------------*/
711
712void
714 const cs_cell_mesh_t *cm,
715 const cs_cdofb_navsto_builder_t *nsb,
716 cs_cell_sys_t *csys);
717
718/*----------------------------------------------------------------------------*/
730/*----------------------------------------------------------------------------*/
731
732void
734 const cs_lnum_t *elt_ids,
735 bool dense_output,
736 void *input,
737 cs_real_t *retval);
738
739/*----------------------------------------------------------------------------*/
740
742
743#endif /* __CS_CDOFB_NAVSTO_H__ */
int cs_boundary_type_t
Definition: cs_boundary.h:69
void cs_cdofb_navsto_gravity_term(const cs_navsto_param_t *nsp, const cs_cell_mesh_t *cm, const cs_cdofb_navsto_builder_t *nsb, cs_cell_sys_t *csys)
Take into account the gravity effects. Compute and add the source term to the local RHS....
Definition: cs_cdofb_navsto.cpp:2037
void cs_cdofb_block_dirichlet_pena(short int f, const cs_equation_param_t *eqp, const cs_cell_mesh_t *cm, const cs_property_data_t *pty, cs_cell_builder_t *cb, cs_cell_sys_t *csys)
Take into account a Dirichlet BCs on the three velocity components. For instance for a velocity inlet...
Definition: cs_cdofb_navsto.cpp:1507
void cs_cdofb_navsto_extra_op(const cs_navsto_param_t *nsp, const cs_mesh_t *mesh, const cs_cdo_quantities_t *quant, const cs_cdo_connect_t *connect, const cs_time_step_t *ts, cs_time_plot_t *time_plotter, const cs_adv_field_t *adv_field, const cs_real_t *mass_flux, const cs_real_t *p_cell, const cs_real_t *u_cell, const cs_real_t *u_face)
Perform extra-operation related to Fb schemes when solving Navier-Stokes. Computation of the followin...
Definition: cs_cdofb_navsto.cpp:969
void cs_cdofb_navsto_add_grad_div(short int n_fc, const cs_real_t zeta, const cs_real_t div[], cs_sdm_t *mat)
Add the grad-div part to the local matrix (i.e. for the current cell)
Definition: cs_cdofb_navsto.cpp:592
void cs_cdofb_navsto_boussinesq_at_face(const cs_navsto_param_t *nsp, const cs_cell_mesh_t *cm, const cs_cdofb_navsto_builder_t *nsb, cs_cell_sys_t *csys)
Take into account the buoyancy force with the Boussinesq approx. Compute and add the source term to t...
Definition: cs_cdofb_navsto.cpp:2138
void cs_cdofb_navsto_compute_face_pressure(const cs_mesh_t *mesh, const cs_cdo_connect_t *connect, const cs_cdo_quantities_t *quant, const cs_time_step_t *ts, const cs_navsto_param_t *nsp, const cs_real_t *p_cell, cs_real_t *p_face)
Compute an estimation of the pressure at faces.
Definition: cs_cdofb_navsto.cpp:520
void cs_cdofb_block_dirichlet_wsym(short int fb, const cs_equation_param_t *eqp, const cs_cell_mesh_t *cm, const cs_property_data_t *pty, cs_cell_builder_t *cb, cs_cell_sys_t *csys)
Take into account a Dirichlet BCs on the three velocity components. For instance for a velocity inlet...
Definition: cs_cdofb_navsto.cpp:1659
void cs_cdofb_symmetry(short int fb, const cs_equation_param_t *eqp, const cs_cell_mesh_t *cm, const cs_property_data_t *pty, cs_cell_builder_t *cb, cs_cell_sys_t *csys)
Take into account a boundary defined as 'symmetry' (treated as a sliding BCs on the three velocity co...
Definition: cs_cdofb_navsto.cpp:1766
cs_sles_convergence_state_t cs_cdofb_navsto_nl_algo_cvg(cs_param_nl_algo_t nl_algo_type, const cs_real_t *pre_iterate, cs_real_t *cur_iterate, cs_iter_algo_t *algo)
Test if one has to do one more non-linear iteration. Test if performed on the relative norm on the in...
Definition: cs_cdofb_navsto.cpp:1923
void cs_cdofb_navsto_init_pressure(const cs_navsto_param_t *nsp, const cs_cdo_quantities_t *quant, const cs_time_step_t *ts, cs_field_t *pr)
Initialize the pressure values.
Definition: cs_cdofb_navsto.cpp:666
void cs_cdofb_block_dirichlet_alge(short int f, const cs_equation_param_t *eqp, const cs_cell_mesh_t *cm, const cs_property_data_t *pty, cs_cell_builder_t *cb, cs_cell_sys_t *csys)
Take into account a Dirichlet BCs on the three velocity components. For instance for a velocity inlet...
Definition: cs_cdofb_navsto.cpp:1405
void cs_cdofb_fixed_wall(short int fb, const cs_equation_param_t *eqp, const cs_cell_mesh_t *cm, const cs_property_data_t *pty, cs_cell_builder_t *cb, cs_cell_sys_t *csys)
Take into account a wall BCs by a weak enforcement using Nitsche technique plus a symmetric treatment...
Definition: cs_cdofb_navsto.cpp:1879
cs_cdofb_navsto_boussinesq_type_t
Type of algorithm to compute the Boussinesq approximation.
Definition: cs_cdofb_navsto.h:71
@ CS_CDOFB_NAVSTO_BOUSSINESQ_CELL_DOF
Boussinesq approximation relyong on a cell contribution.
Definition: cs_cdofb_navsto.h:82
@ CS_CDOFB_NAVSTO_BOUSSINESQ_FACE_DOF
Boussinesq approximation relyong on face contributions.
Definition: cs_cdofb_navsto.h:93
void cs_cdofb_navsto_rescale_pressure_to_ref(const cs_navsto_param_t *nsp, const cs_cdo_quantities_t *quant, cs_real_t values[])
Update the pressure field in order to get a field with a mean-value equal to the reference value.
Definition: cs_cdofb_navsto.cpp:866
void cs_cdofb_navsto_free_builder(cs_cdofb_navsto_builder_t *nsb)
Destroy the given cs_cdofb_navsto_builder_t structure.
Definition: cs_cdofb_navsto.cpp:246
double cs_cdofb_navsto_cell_divergence(const cs_lnum_t c_id, const cs_cdo_quantities_t *quant, const cs_adjacency_t *c2f, const cs_real_t *f_vals)
Compute the divergence in a cell of a vector-valued array defined at faces (values are defined both a...
Definition: cs_cdofb_navsto.cpp:481
void cs_cdofb_navsto_set_boussinesq_algo(cs_cdofb_navsto_boussinesq_type_t type)
Set the way to compute the Boussinesq approximation.
Definition: cs_cdofb_navsto.cpp:198
void cs_cdofb_navsto_boussinesq_at_cell(const cs_navsto_param_t *nsp, const cs_cell_mesh_t *cm, const cs_cdofb_navsto_builder_t *nsb, cs_cell_sys_t *csys)
Take into account the buoyancy force with the Boussinesq approx. Compute and add the source term to t...
Definition: cs_cdofb_navsto.cpp:2078
void cs_cdofb_navsto_stream_source_term(cs_lnum_t n_elts, const cs_lnum_t *elt_ids, bool dense_output, void *input, cs_real_t *retval)
Get the source term for computing the stream function. This relies on the prototype associated to the...
Definition: cs_cdofb_navsto.cpp:2195
void cs_cdofb_navsto_init_face_pressure(const cs_navsto_param_t *nsp, const cs_cdo_connect_t *connect, const cs_time_step_t *ts, cs_real_t *pr_f)
Initialize the pressure values when the pressure is defined at faces.
Definition: cs_cdofb_navsto.cpp:762
void cs_cdofb_navsto_define_builder(cs_real_t t_eval, const cs_navsto_param_t *nsp, const cs_cell_mesh_t *cm, const cs_cell_sys_t *csys, const cs_cdo_bc_face_t *pr_bc, const cs_boundary_type_t *bf_type, cs_cdofb_navsto_builder_t *nsb)
Set the members of the cs_cdofb_navsto_builder_t structure.
Definition: cs_cdofb_navsto.cpp:270
void cs_cdofb_navsto_set_gravity_func(const cs_navsto_param_t *nsp, cs_cdofb_navsto_source_t **p_func)
Set the function pointer computing the source term in the momentum equation related to the gravity ef...
Definition: cs_cdofb_navsto.cpp:1993
cs_cdofb_navsto_builder_t cs_cdofb_navsto_create_builder(const cs_navsto_param_t *nsp, const cs_cdo_connect_t *connect)
Create and allocate a local NavSto builder when Fb schemes are used.
Definition: cs_cdofb_navsto.cpp:215
void cs_cdofb_block_dirichlet_weak(short int fb, const cs_equation_param_t *eqp, const cs_cell_mesh_t *cm, const cs_property_data_t *pty, cs_cell_builder_t *cb, cs_cell_sys_t *csys)
Take into account a Dirichlet BCs on the three velocity components. For instance for a velocity inlet...
Definition: cs_cdofb_navsto.cpp:1573
void cs_cdofb_navsto_set_zero_mean_pressure(const cs_cdo_quantities_t *quant, cs_real_t values[])
Update the pressure field in order to get a field with a zero-mean average.
Definition: cs_cdofb_navsto.cpp:904
void() cs_cdofb_navsto_source_t(const cs_navsto_param_t *nsp, const cs_cell_mesh_t *cm, const cs_cdofb_navsto_builder_t *nsb, cs_cell_sys_t *csys)
Compute and add a source term to the local RHS. This is a special treatment to enable source involvin...
Definition: cs_cdofb_navsto.h:184
void cs_cdofb_navsto_mass_flux(const cs_navsto_param_t *nsp, const cs_cdo_quantities_t *quant, const cs_real_t *face_vel, cs_real_t *mass_flux)
Compute the mass flux playing the role of the advection field in the Navier-Stokes equations One cons...
Definition: cs_cdofb_navsto.cpp:439
void cs_cdofb_navsto_boussinesq_by_part(const cs_navsto_param_t *nsp, const cs_cell_mesh_t *cm, const cs_cdofb_navsto_builder_t *nsb, cs_cell_sys_t *csys)
Take into account the buoyancy force with the Boussinesq approx. Compute and add the source term to t...
#define BEGIN_C_DECLS
Definition: cs_defs.h:542
double cs_real_t
Floating-point value.
Definition: cs_defs.h:342
#define END_C_DECLS
Definition: cs_defs.h:543
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:335
cs_param_nl_algo_t
Class of non-linear iterative algorithm.
Definition: cs_param_types.h:606
cs_sles_convergence_state_t
Definition: cs_sles.h:56
struct _cs_time_plot_t cs_time_plot_t
Definition: cs_time_plot.h:48
char * algo
Definition: field_names.h:100
Definition: mesh.f90:26
Definition: cs_mesh_adjacencies.h:68
Definition: cs_advection_field.h:151
Definition: cs_cdo_bc.h:106
Definition: cs_cdo_connect.h:61
Definition: cs_cdo_quantities.h:139
Structure storing additional arrays related to the building of the Navier-Stokes system.
Definition: cs_cdofb_navsto.h:107
cs_real_t mass_rhs
Definition: cs_cdofb_navsto.h:139
cs_real_t * div_op
Definition: cs_cdofb_navsto.h:128
cs_real_t rho_c
Definition: cs_cdofb_navsto.h:115
cs_boundary_type_t * bf_type
Definition: cs_cdofb_navsto.h:159
cs_real_t * pressure_bc_val
Definition: cs_cdofb_navsto.h:160
Set of local and temporary buffers.
Definition: cs_cdo_local.h:60
Set of local quantities and connectivities related to a mesh cell.
Definition: cs_cdo_local.h:203
Set of arrays and local (small) dense matrices related to a mesh cell This is a key structure for bui...
Definition: cs_cdo_local.h:147
Set of parameters to handle an unsteady convection-diffusion-reaction equation with term sources.
Definition: cs_equation_param.h:192
Field descriptor.
Definition: cs_field.h:131
Structure to handle the convergence of an iterative algorithm.
Definition: cs_iter_algo.h:290
Definition: cs_mesh.h:85
Structure storing the parameters related to the resolution of the Navier-Stokes system.
Definition: cs_navsto_param.h:263
Structure storing the evaluation of a property and its related data.
Definition: cs_property.h:223
time step descriptor
Definition: cs_time_step.h:64