7.1
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-2021 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 
71 typedef enum {
72 
83 
94 
96 
97 
107 typedef struct {
108 
116 
129 
148  cs_boundary_type_t *bf_type; /* Size: n_fc */
149  cs_real_t *pressure_bc_val; /* Size: n_fc */
150 
152 
153 /*----------------------------------------------------------------------------*/
166 /*----------------------------------------------------------------------------*/
167 
168 typedef void
170  const cs_cell_mesh_t *cm,
171  const cs_cdofb_navsto_builder_t *nsb,
172  cs_cell_sys_t *csys);
173 
174 /*============================================================================
175  * Static inline public function prototypes
176  *============================================================================*/
177 
178 /*----------------------------------------------------------------------------*/
187 /*----------------------------------------------------------------------------*/
188 
189 static inline void
190 cs_cdofb_navsto_divergence_vect(const cs_cell_mesh_t *cm,
191  cs_real_t div[])
192 {
193  /* D(\hat{u}) = \frac{1}{|c|} \sum_{f_c} \iota_{fc} u_f.f
194  * But, when integrating [[ p, q ]]_{P_c} = |c| p_c q_c
195  * Thus, the volume in the divergence drops
196  */
197 
198  for (short int f = 0; f < cm->n_fc; f++) {
199 
200  const cs_quant_t pfq = cm->face[f];
201  const cs_real_t i_f = cm->f_sgn[f] * pfq.meas;
202 
203  cs_real_t *_div_f = div + 3*f;
204  _div_f[0] = i_f * pfq.unitv[0];
205  _div_f[1] = i_f * pfq.unitv[1];
206  _div_f[2] = i_f * pfq.unitv[2];
207 
208  } /* Loop on cell faces */
209 }
210 
211 /*============================================================================
212  * Public function prototypes
213  *============================================================================*/
214 
215 /*----------------------------------------------------------------------------*/
221 /*----------------------------------------------------------------------------*/
222 
223 void
225 
226 /*----------------------------------------------------------------------------*/
235 /*----------------------------------------------------------------------------*/
236 
239  const cs_cdo_connect_t *connect);
240 
241 /*----------------------------------------------------------------------------*/
247 /*----------------------------------------------------------------------------*/
248 
249 void
251 
252 /*----------------------------------------------------------------------------*/
264 /*----------------------------------------------------------------------------*/
265 
266 void
268  const cs_navsto_param_t *nsp,
269  const cs_cell_mesh_t *cm,
270  const cs_cell_sys_t *csys,
271  const cs_cdo_bc_face_t *pr_bc,
272  const cs_boundary_type_t *bf_type,
274 
275 /*----------------------------------------------------------------------------*/
287 /*----------------------------------------------------------------------------*/
288 
289 void
291  const cs_cdo_quantities_t *quant,
292  const cs_real_t *face_vel,
293  cs_real_t *mass_flux);
294 
295 /*----------------------------------------------------------------------------*/
308 /*----------------------------------------------------------------------------*/
309 
310 double
312  const cs_cdo_quantities_t *quant,
313  const cs_adjacency_t *c2f,
314  const cs_real_t *f_vals);
315 
316 /*----------------------------------------------------------------------------*/
328 /*----------------------------------------------------------------------------*/
329 
330 void
332  const cs_cdo_connect_t *connect,
333  const cs_cdo_quantities_t *quant,
334  const cs_time_step_t *ts,
335  const cs_navsto_param_t *nsp,
336  const cs_real_t *p_cell,
337  cs_real_t *p_face);
338 
339 /*----------------------------------------------------------------------------*/
349 /*----------------------------------------------------------------------------*/
350 
351 void
352 cs_cdofb_navsto_add_grad_div(short int n_fc,
353  const cs_real_t zeta,
354  const cs_real_t div[],
355  cs_sdm_t *mat);
356 
357 /*----------------------------------------------------------------------------*/
366 /*----------------------------------------------------------------------------*/
367 
368 void
370  const cs_cdo_quantities_t *quant,
371  const cs_time_step_t *ts,
372  cs_field_t *pr);
373 
374 /*----------------------------------------------------------------------------*/
384 /*----------------------------------------------------------------------------*/
385 
386 void
388  const cs_cdo_connect_t *connect,
389  const cs_time_step_t *ts,
390  cs_real_t *pr_f);
391 
392 /*----------------------------------------------------------------------------*/
401 /*----------------------------------------------------------------------------*/
402 
403 void
405  const cs_cdo_quantities_t *quant,
406  cs_real_t values[]);
407 
408 /*----------------------------------------------------------------------------*/
416 /*----------------------------------------------------------------------------*/
417 
418 void
420  cs_real_t values[]);
421 
422 /*----------------------------------------------------------------------------*/
451 /*----------------------------------------------------------------------------*/
452 
453 void
455  const cs_mesh_t *mesh,
456  const cs_cdo_quantities_t *quant,
457  const cs_cdo_connect_t *connect,
458  const cs_time_step_t *ts,
459  cs_time_plot_t *time_plotter,
460  const cs_adv_field_t *adv_field,
461  const cs_real_t *mass_flux,
462  const cs_real_t *p_cell,
463  const cs_real_t *u_cell,
464  const cs_real_t *u_face);
465 
466 /*----------------------------------------------------------------------------*/
481 /*----------------------------------------------------------------------------*/
482 
483 void
485  const cs_equation_param_t *eqp,
486  const cs_cell_mesh_t *cm,
487  const cs_property_data_t *pty,
488  cs_cell_builder_t *cb,
489  cs_cell_sys_t *csys);
490 
491 /*----------------------------------------------------------------------------*/
508 /*----------------------------------------------------------------------------*/
509 
510 void
512  const cs_equation_param_t *eqp,
513  const cs_cell_mesh_t *cm,
514  const cs_property_data_t *pty,
515  cs_cell_builder_t *cb,
516  cs_cell_sys_t *csys);
517 
518 /*----------------------------------------------------------------------------*/
535 /*----------------------------------------------------------------------------*/
536 
537 void
538 cs_cdofb_block_dirichlet_weak(short int fb,
539  const cs_equation_param_t *eqp,
540  const cs_cell_mesh_t *cm,
541  const cs_property_data_t *pty,
542  cs_cell_builder_t *cb,
543  cs_cell_sys_t *csys);
544 
545 /*----------------------------------------------------------------------------*/
562 /*----------------------------------------------------------------------------*/
563 
564 void
565 cs_cdofb_block_dirichlet_wsym(short int fb,
566  const cs_equation_param_t *eqp,
567  const cs_cell_mesh_t *cm,
568  const cs_property_data_t *pty,
569  cs_cell_builder_t *cb,
570  cs_cell_sys_t *csys);
571 
572 /*----------------------------------------------------------------------------*/
588 /*----------------------------------------------------------------------------*/
589 
590 void
591 cs_cdofb_symmetry(short int fb,
592  const cs_equation_param_t *eqp,
593  const cs_cell_mesh_t *cm,
594  const cs_property_data_t *pty,
595  cs_cell_builder_t *cb,
596  cs_cell_sys_t *csys);
597 
598 /*----------------------------------------------------------------------------*/
611 /*----------------------------------------------------------------------------*/
612 
613 void
614 cs_cdofb_fixed_wall(short int fb,
615  const cs_equation_param_t *eqp,
616  const cs_cell_mesh_t *cm,
617  const cs_property_data_t *pty,
618  cs_cell_builder_t *cb,
619  cs_cell_sys_t *csys);
620 
621 /*----------------------------------------------------------------------------*/
634 /*----------------------------------------------------------------------------*/
635 
638  const cs_real_t *pre_iterate,
639  cs_real_t *cur_iterate,
640  cs_iter_algo_t *algo);
641 
642 /*----------------------------------------------------------------------------*/
651 /*----------------------------------------------------------------------------*/
652 
653 void
655  cs_cdofb_navsto_source_t **p_func);
656 
657 /*----------------------------------------------------------------------------*/
669 /*----------------------------------------------------------------------------*/
670 
671 void
673  const cs_cell_mesh_t *cm,
674  const cs_cdofb_navsto_builder_t *nsb,
675  cs_cell_sys_t *csys);
676 
677 /*----------------------------------------------------------------------------*/
690 /*----------------------------------------------------------------------------*/
691 
692 void
694  const cs_cell_mesh_t *cm,
695  const cs_cdofb_navsto_builder_t *nsb,
696  cs_cell_sys_t *csys);
697 
698 /*----------------------------------------------------------------------------*/
711 /*----------------------------------------------------------------------------*/
712 
713 void
715  const cs_cell_mesh_t *cm,
716  const cs_cdofb_navsto_builder_t *nsb,
717  cs_cell_sys_t *csys);
718 
719 /*----------------------------------------------------------------------------*/
732 /*----------------------------------------------------------------------------*/
733 
734 void
736  const cs_cell_mesh_t *cm,
737  const cs_cdofb_navsto_builder_t *nsb,
738  cs_cell_sys_t *csys);
739 
740 /*----------------------------------------------------------------------------*/
752 /*----------------------------------------------------------------------------*/
753 
754 void
756  const cs_lnum_t *elt_ids,
757  bool dense_output,
758  void *input,
759  cs_real_t *retval);
760 
761 /*----------------------------------------------------------------------------*/
762 
764 
765 #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:306
time step descriptor
Definition: cs_time_step.h:64
int cs_boundary_type_t
Definition: cs_boundary.h:69
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:252
Definition: cs_advection_field.h:149
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.c:1591
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:875
Field descriptor.
Definition: cs_field.h:125
Set of parameters to handle an unsteady convection-diffusion-reaction equation with term sources...
Definition: cs_equation_param.h:177
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:609
Structure storing the evaluation of a property and its related data.
Definition: cs_property.h:178
Definition: cs_mesh_adjacencies.h:68
cs_param_nl_algo_t
Class of non-linear iterative algorithm.
Definition: cs_param_types.h:511
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.c:1674
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 &#39;symmetry&#39; (treated as a sliding BCs on the three velocity co...
Definition: cs_cdofb_navsto.c:1775
#define BEGIN_C_DECLS
Definition: cs_defs.h:510
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:169
double meas
Definition: cs_cdo_quantities.h:123
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.c:978
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:1524
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.c:2164
Set of local quantities and connectivities related to a mesh cell.
Definition: cs_cdo_local.h:202
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:282
Definition: cs_cdo_connect.h:79
cs_real_t * div_op
Definition: cs_cdofb_navsto.h:128
Structure storing the parameters related to the resolution of the Navier-Stokes system.
Definition: cs_navsto_param.h:605
double cs_real_t
Floating-point value.
Definition: cs_defs.h:322
Definition: cs_cdo_quantities.h:129
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...
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:146
struct _cs_time_plot_t cs_time_plot_t
Definition: cs_time_plot.h:48
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:679
Set of common parameters to manage an iterative algorithm.
Definition: cs_iter_algo.h:101
short int * f_sgn
Definition: cs_cdo_local.h:235
Definition: cs_mesh.h:84
cs_sles_convergence_state_t
Convergence status indicator.
Definition: cs_sles.h:56
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:1424
Boussinesq approximation relyong on face contributions.
Definition: cs_cdofb_navsto.h:93
cs_cdofb_navsto_boussinesq_type_t
Type of algorithm to compute the Boussinesq approximation.
Definition: cs_cdofb_navsto.h:71
short int n_fc
Definition: cs_cdo_local.h:232
Definition: cs_cdo_bc.h:88
Set of local and temporary buffers.
Definition: cs_cdo_local.h:60
cs_quant_t * face
Definition: cs_cdo_local.h:238
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:773
Definition: cs_cdo_quantities.h:121
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:316
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.c:2063
double unitv[3]
Definition: cs_cdo_quantities.h:124
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.c:235
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.c:456
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.c:498
#define END_C_DECLS
Definition: cs_defs.h:511
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.c:2221
cs_real_t rho_c
Definition: cs_cdofb_navsto.h:115
cs_real_t * pressure_bc_val
Definition: cs_cdofb_navsto.h:149
Boussinesq approximation relyong on a cell contribution.
Definition: cs_cdofb_navsto.h:82
double precision, dimension(ncharm), save zeta
Definition: cpincl.f90:99
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.c:537
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.c:1926
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.c:2019
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.c:2104
Structure storing additional arrays related to the building of the Navier-Stokes system.
Definition: cs_cdofb_navsto.h:107
cs_boundary_type_t * bf_type
Definition: cs_cdofb_navsto.h:148
Definition: mesh.f90:26
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.c:1882
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:913