6.2
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  * Routines 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-2020 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_math.h"
46 #include "cs_matrix.h"
47 #include "cs_mesh.h"
48 #include "cs_navsto_param.h"
49 #include "cs_time_step.h"
50 #include "cs_sdm.h"
51 
52 /*----------------------------------------------------------------------------*/
53 
55 
56 /*============================================================================
57  * Macro definitions
58  *============================================================================*/
59 
60 /*============================================================================
61  * Type definitions
62  *============================================================================*/
63 
64 /* Structure storing additional arrays related to the building of the system in
65  case of CDO Face-based scheme. This structure is associated to a cell-wise
66  building */
67 
68 typedef struct {
69 
70  /* Value of the mass density for the current cell */
72 
73  /* Operator */
74  cs_real_t *div_op; /* Size: 3*n_fc
75  div_op = -|c|div */
76 
77  /* Boundary conditions */
78  cs_boundary_type_t *bf_type; /* Size: n_fc */
79  cs_real_t *pressure_bc_val; /* Size: n_fc */
80 
82 
83 /*============================================================================
84  * Static inline public function prototypes
85  *============================================================================*/
86 
87 /*----------------------------------------------------------------------------*/
96 /*----------------------------------------------------------------------------*/
97 
98 static inline void
99 cs_cdofb_navsto_divergence_vect(const cs_cell_mesh_t *cm,
100  cs_real_t div[])
101 {
102  /* D(\hat{u}) = \frac{1}{|c|} \sum_{f_c} \iota_{fc} u_f.f
103  * But, when integrating [[ p, q ]]_{P_c} = |c| p_c q_c
104  * Thus, the volume in the divergence drops
105  */
106 
107  for (short int f = 0; f < cm->n_fc; f++) {
108 
109  const cs_quant_t pfq = cm->face[f];
110  const cs_real_t i_f = cm->f_sgn[f] * pfq.meas;
111 
112  cs_real_t *_div_f = div + 3*f;
113  _div_f[0] = i_f * pfq.unitv[0];
114  _div_f[1] = i_f * pfq.unitv[1];
115  _div_f[2] = i_f * pfq.unitv[2];
116 
117  } /* Loop on cell faces */
118 }
119 
120 /*============================================================================
121  * Public function prototypes
122  *============================================================================*/
123 
124 /*----------------------------------------------------------------------------*/
133 /*----------------------------------------------------------------------------*/
134 
137  const cs_cdo_connect_t *connect);
138 
139 /*----------------------------------------------------------------------------*/
145 /*----------------------------------------------------------------------------*/
146 
147 void
149 
150 /*----------------------------------------------------------------------------*/
162 /*----------------------------------------------------------------------------*/
163 
164 void
166  const cs_navsto_param_t *nsp,
167  const cs_cell_mesh_t *cm,
168  const cs_cell_sys_t *csys,
169  const cs_cdo_bc_face_t *pr_bc,
170  const cs_boundary_type_t *bf_type,
172 
173 /*----------------------------------------------------------------------------*/
185 /*----------------------------------------------------------------------------*/
186 
187 void
189  const cs_cdo_quantities_t *quant,
190  const cs_real_t *face_vel,
191  cs_adv_field_t *adv);
192 
193 /*----------------------------------------------------------------------------*/
205 /*----------------------------------------------------------------------------*/
206 
207 cs_real_t
209  const cs_cdo_quantities_t *quant,
210  const cs_adjacency_t *c2f,
211  const cs_real_t *f_vals);
212 
213 /*----------------------------------------------------------------------------*/
223 /*----------------------------------------------------------------------------*/
224 
225 void
226 cs_cdofb_navsto_add_grad_div(short int n_fc,
227  const cs_real_t zeta,
228  const cs_real_t div[],
229  cs_sdm_t *mat);
230 
231 /*----------------------------------------------------------------------------*/
240 /*----------------------------------------------------------------------------*/
241 
242 void
244  const cs_cdo_quantities_t *quant,
245  const cs_time_step_t *ts,
246  cs_field_t *pr);
247 
248 /*----------------------------------------------------------------------------*/
258 /*----------------------------------------------------------------------------*/
259 
260 void
262  const cs_cdo_connect_t *connect,
263  const cs_time_step_t *ts,
264  cs_real_t *pr_f);
265 
266 /*----------------------------------------------------------------------------*/
275 /*----------------------------------------------------------------------------*/
276 
277 void
279  const cs_cdo_quantities_t *quant,
280  cs_real_t values[]);
281 
282 /*----------------------------------------------------------------------------*/
290 /*----------------------------------------------------------------------------*/
291 
292 void
294  cs_real_t values[]);
295 
296 /*----------------------------------------------------------------------------*/
318 /*----------------------------------------------------------------------------*/
319 
320 void
322  const cs_mesh_t *mesh,
323  const cs_cdo_quantities_t *quant,
324  const cs_cdo_connect_t *connect,
325  const cs_time_step_t *ts,
326  const cs_adv_field_t *adv_field,
327  const cs_real_t *mass_flux,
328  const cs_real_t *u_cell,
329  const cs_real_t *u_face);
330 
331 /*----------------------------------------------------------------------------*/
346 /*----------------------------------------------------------------------------*/
347 
348 void
350  const cs_equation_param_t *eqp,
351  const cs_cell_mesh_t *cm,
352  const cs_property_data_t *pty,
353  cs_cell_builder_t *cb,
354  cs_cell_sys_t *csys);
355 
356 /*----------------------------------------------------------------------------*/
373 /*----------------------------------------------------------------------------*/
374 
375 void
377  const cs_equation_param_t *eqp,
378  const cs_cell_mesh_t *cm,
379  const cs_property_data_t *pty,
380  cs_cell_builder_t *cb,
381  cs_cell_sys_t *csys);
382 
383 /*----------------------------------------------------------------------------*/
400 /*----------------------------------------------------------------------------*/
401 
402 void
404  const cs_equation_param_t *eqp,
405  const cs_cell_mesh_t *cm,
406  const cs_property_data_t *pty,
407  cs_cell_builder_t *cb,
408  cs_cell_sys_t *csys);
409 
410 /*----------------------------------------------------------------------------*/
427 /*----------------------------------------------------------------------------*/
428 
429 void
431  const cs_equation_param_t *eqp,
432  const cs_cell_mesh_t *cm,
433  const cs_property_data_t *pty,
434  cs_cell_builder_t *cb,
435  cs_cell_sys_t *csys);
436 
437 /*----------------------------------------------------------------------------*/
453 /*----------------------------------------------------------------------------*/
454 
455 void
456 cs_cdofb_symmetry(short int f,
457  const cs_equation_param_t *eqp,
458  const cs_cell_mesh_t *cm,
459  const cs_property_data_t *pty,
460  cs_cell_builder_t *cb,
461  cs_cell_sys_t *csys);
462 
463 /*----------------------------------------------------------------------------*/
477 /*----------------------------------------------------------------------------*/
478 
479 void
480 cs_cdofb_fixed_wall(short int f,
481  const cs_equation_param_t *eqp,
482  const cs_cell_mesh_t *cm,
483  const cs_property_data_t *pty,
484  cs_cell_builder_t *cb,
485  cs_cell_sys_t *csys);
486 
487 /*----------------------------------------------------------------------------*/
499 /*----------------------------------------------------------------------------*/
500 
501 void
503  const cs_lnum_t *elt_ids,
504  bool compact,
505  void *input,
506  cs_real_t *retval);
507 
508 /*----------------------------------------------------------------------------*/
520 /*----------------------------------------------------------------------------*/
521 
522 void
524  const cs_lnum_t *elt_ids,
525  bool compact,
526  void *input,
527  cs_real_t *retval);
528 
529 /*----------------------------------------------------------------------------*/
530 
532 
533 #endif /* __CS_CDOFB_NAVSTO_H__ */
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.c:240
time step descriptor
Definition: cs_time_step.h:64
void cs_cdofb_block_dirichlet_wsym(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.c:1372
int cs_boundary_type_t
Definition: cs_boundary.h:65
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.c:186
Definition: cs_advection_field.h:149
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.c:731
Field descriptor.
Definition: cs_field.h:125
void cs_cdofb_symmetry(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 symmetric boundary (treated as a sliding BCs on the three velocity components...
Definition: cs_cdofb_navsto.c:1469
Set of parameters to handle an unsteady convection-diffusion-reaction equation with term sources...
Definition: cs_equation_param.h:201
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.c:463
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, const cs_adv_field_t *adv_field, const cs_real_t *mass_flux, const cs_real_t *u_cell, const cs_real_t *u_face)
Perform extra-operation related to Fb schemes when solving Navier-Stokes.
Definition: cs_cdofb_navsto.c:829
cs_real_t * pressure_bc_val
Definition: cs_cdofb_navsto.h:79
Structure storing the evaluation of a property and its related data.
Definition: cs_property.h:174
Definition: cs_mesh_adjacencies.h:67
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_adv_field_t *adv)
Compute the mass flux playing the role of the advection field in the Navier-Stokes equations One cons...
Definition: cs_cdofb_navsto.c:384
#define BEGIN_C_DECLS
Definition: cs_defs.h:495
double meas
Definition: cs_cdo_quantities.h:118
cs_real_t rho_c
Definition: cs_cdofb_navsto.h:71
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.c:1224
Set of local quantities and connectivities related to a mesh cell This is a key structure for all cel...
Definition: cs_cdo_local.h:158
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.c:216
Definition: cs_cdo_connect.h:76
Structure storing the parameters related to the resolution of the Navier-Stokes system.
Definition: cs_navsto_param.h:462
cs_real_t * div_op
Definition: cs_cdofb_navsto.h:74
double cs_real_t
Floating-point value.
Definition: cs_defs.h:307
Definition: cs_cdo_quantities.h:124
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:107
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.c:533
short int * f_sgn
Definition: cs_cdo_local.h:191
Definition: cs_mesh.h:74
cs_real_t 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 of a cell using the cs_cdo_quantities_t structure.
Definition: cs_cdofb_navsto.c:426
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.c:1124
short int n_fc
Definition: cs_cdo_local.h:188
Definition: cs_cdo_bc.h:88
void cs_cdofb_fixed_wall(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 wall BCs by a weak enforcement using Nitsche technique plus a symmetric treatment...
Definition: cs_cdofb_navsto.c:1573
Set of local and temporary buffers useful for building the algebraic system with a cellwise process...
Definition: cs_cdo_local.h:69
cs_quant_t * face
Definition: cs_cdo_local.h:194
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.c:629
Definition: cs_cdo_quantities.h:116
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:301
double unitv[3]
Definition: cs_cdo_quantities.h:119
#define END_C_DECLS
Definition: cs_defs.h:496
void cs_cdofb_navsto_stream_source_term(cs_lnum_t n_elts, const cs_lnum_t *elt_ids, bool compact, 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.c:1656
void cs_cdofb_block_dirichlet_weak(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.c:1291
cs_boundary_type_t * bf_type
Definition: cs_cdofb_navsto.h:78
void cs_cdofb_navsto_boussinesq_source_term(cs_lnum_t n_elts, const cs_lnum_t *elt_ids, bool compact, void *input, cs_real_t *retval)
Get the source term for computing the Boussinesq approximation This relies on the prototype associate...
Definition: cs_cdofb_navsto.c:1616
double precision, dimension(ncharm), save zeta
Definition: cpincl.f90:99
Definition: cs_cdofb_navsto.h:68
Definition: mesh.f90:26
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.c:770