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
56 BEGIN_C_DECLS
57
58 /*============================================================================
59 * Macro definitions
60 *============================================================================*/
61
62 /*============================================================================
63 * Type definitions
64 *============================================================================*/
65
66 /*!
67 * \enum cs_cdofb_navsto_boussinesq_type_t
68 * \brief Type of algorithm to compute the Boussinesq approximation
69 */
70
71 typedef enum {
72
73 /*!
74 * \brief Boussinesq approximation relyong on a cell contribution
75 *
76 * This algorithm uses cell DoFs for the Boussinesq part corresponding to
77 * rho0 * beta * (var[c_id] - var0) * g[] while the constant part equal to
78 * rho0 * g[] is built in order to be in balance with the pressure gradient
79 * (face DoF contributions).
80 */
81
82 CS_CDOFB_NAVSTO_BOUSSINESQ_CELL_DOF,
83
84 /*!
85 * \brief Boussinesq approximation relyong on face contributions
86 *
87 * This algorithm uses only face DoFs for the Boussinesq approximation.
88 * For the constant part (rho0 * g[]) as well as the variable part
89 * rho0 * beta * (var[c_id] - var0) * g[]
90 * The aim is to be in balance with the pressure gradient
91 */
92
93 CS_CDOFB_NAVSTO_BOUSSINESQ_FACE_DOF
94
95 } cs_cdofb_navsto_boussinesq_type_t;
96
97
98 /*! \struct cs_cdofb_navsto_builder_t
99 *
100 * \brief Structure storing additional arrays related to the building of the
101 * Navier-Stokes system.
102 *
103 * This structure is associated to a cell-wise building in case of CDO
104 * face-based scheme.
105 */
106
107 typedef struct {
108
109 /*!
110 * \var rho_c
111 * Value of the mass density for the current cell
112 *
113 */
114
115 cs_real_t rho_c;
116
117 /*!
118 * @}
119 * @name Operator(s)
120 * @{
121 *
122 * \var div_op
123 * Cell-wise divergence operator (in fact div_op = -|c| div).
124 * This operator is stored in an array of size 3*n_fc
125 *
126 */
127
128 cs_real_t *div_op;
129
130 /*!
131 * @}
132 * @name Boundary conditions
133 * @{
134 *
135 * \var bf_type
136 * Type of boundary to apply to each face.
137 *
138 * Array of size n_fc. Zero is set for interior faces.
139 *
140 *
141 * \var pressure_bc_val
142 * Store the value of the pressure on boundary faces.
143 *
144 * Array of size n_fc. Only useful if a Dirichlet boundary condition is set
145 * on a boundary face.
146 */
147
148 cs_boundary_type_t *bf_type; /* Size: n_fc */
149 cs_real_t *pressure_bc_val; /* Size: n_fc */
150
151 } cs_cdofb_navsto_builder_t;
152
153 /*----------------------------------------------------------------------------*/
154 /*!
155 * \brief Compute and add a source term to the local RHS.
156 * This is a special treatment to enable source involving face DoFs and
157 * potentially the local discrete divergence/gradient operators.
158 * In the standard case, only the cell DoFs are involved.
159 * Examples are gravity term or Bousinesq term(s)
160 *
161 * \param[in] nsp set of parameters to handle the Navier-Stokes system
162 * \param[in] cm pointer to a cs_cell_mesh_t structure
163 * \param[in] nsb pointer to a builder structure for the NavSto system
164 * \param[in, out] csys pointer to a cs_cell_sys_t structure
165 */
166 /*----------------------------------------------------------------------------*/
167
168 typedef void
169 (cs_cdofb_navsto_source_t)(const cs_navsto_param_t *nsp,
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 /*----------------------------------------------------------------------------*/
179 /*!
180 * \brief Compute the divergence vector associated to the current cell.
181 * WARNING: mind that, differently form the original definition, the
182 * result here is not divided by the cell volume
183 *
184 * \param[in] cm pointer to a \ref cs_cell_mesh_t structure
185 * \param[in, out] div array related to the divergence operator
186 */
187 /*----------------------------------------------------------------------------*/
188
189 static inline void
cs_cdofb_navsto_divergence_vect(const cs_cell_mesh_t * cm,cs_real_t div[])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 /*----------------------------------------------------------------------------*/
216 /*!
217 * \brief Set the way to compute the Boussinesq approximation
218 *
219 * \param[in] type type of algorithm to use
220 */
221 /*----------------------------------------------------------------------------*/
222
223 void
224 cs_cdofb_navsto_set_boussinesq_algo(cs_cdofb_navsto_boussinesq_type_t type);
225
226 /*----------------------------------------------------------------------------*/
227 /*!
228 * \brief Create and allocate a local NavSto builder when Fb schemes are used
229 *
230 * \param[in] nsp set of parameters to define the NavSto system
231 * \param[in] connect pointer to a cs_cdo_connect_t structure
232 *
233 * \return a cs_cdofb_navsto_builder_t structure
234 */
235 /*----------------------------------------------------------------------------*/
236
237 cs_cdofb_navsto_builder_t
238 cs_cdofb_navsto_create_builder(const cs_navsto_param_t *nsp,
239 const cs_cdo_connect_t *connect);
240
241 /*----------------------------------------------------------------------------*/
242 /*!
243 * \brief Destroy the given cs_cdofb_navsto_builder_t structure
244 *
245 * \param[in, out] nsb pointer to the cs_cdofb_navsto_builder_t to free
246 */
247 /*----------------------------------------------------------------------------*/
248
249 void
250 cs_cdofb_navsto_free_builder(cs_cdofb_navsto_builder_t *nsb);
251
252 /*----------------------------------------------------------------------------*/
253 /*!
254 * \brief Set the members of the cs_cdofb_navsto_builder_t structure
255 *
256 * \param[in] t_eval time at which one evaluates the pressure BC
257 * \param[in] nsp set of parameters to define the NavSto system
258 * \param[in] cm cellwise view of the mesh
259 * \param[in] csys cellwise view of the algebraic system
260 * \param[in] pr_bc set of definitions for the presuure BCs
261 * \param[in] bf_type type of boundaries for all boundary faces
262 * \param[in, out] nsb builder to update
263 */
264 /*----------------------------------------------------------------------------*/
265
266 void
267 cs_cdofb_navsto_define_builder(cs_real_t t_eval,
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,
273 cs_cdofb_navsto_builder_t *nsb);
274
275 /*----------------------------------------------------------------------------*/
276 /*!
277 * \brief Compute the mass flux playing the role of the advection field in
278 * the Navier-Stokes equations
279 * One considers the mass flux across primal faces which relies on the
280 * velocity vector defined on each face.
281 *
282 * \param[in] nsp set of parameters to define the NavSto system
283 * \param[in] quant set of additional geometrical quantities
284 * \param[in] face_vel velocity vectors for each face
285 * \param[in, out] mass_flux array of mass flux values to update (allocated)
286 */
287 /*----------------------------------------------------------------------------*/
288
289 void
290 cs_cdofb_navsto_mass_flux(const cs_navsto_param_t *nsp,
291 const cs_cdo_quantities_t *quant,
292 const cs_real_t *face_vel,
293 cs_real_t *mass_flux);
294
295 /*----------------------------------------------------------------------------*/
296 /*!
297 * \brief Compute the divergence in a cell of a vector-valued array defined at
298 * faces (values are defined both at interior and border faces).
299 * Variant based on the usage of \ref cs_cdo_quantities_t structure.
300 *
301 * \param[in] c_id cell id
302 * \param[in] quant pointer to a \ref cs_cdo_quantities_t
303 * \param[in] c2f pointer to cell-to-face \ref cs_adjacency_t
304 * \param[in] f_vals values of the face DoFs
305 *
306 * \return the divergence for the corresponding cell
307 */
308 /*----------------------------------------------------------------------------*/
309
310 double
311 cs_cdofb_navsto_cell_divergence(const cs_lnum_t c_id,
312 const cs_cdo_quantities_t *quant,
313 const cs_adjacency_t *c2f,
314 const cs_real_t *f_vals);
315
316 /*----------------------------------------------------------------------------*/
317 /*!
318 * \brief Compute an estimation of the pressure at faces
319 *
320 * \param[in] mesh pointer to a cs_mesh_t structure
321 * \param[in] connect pointer to a cs_cdo_connect_t structure
322 * \param[in] quant pointer to a cs_cdo_quantities_t structure
323 * \param[in] time_step pointer to a cs_time_step_t structure
324 * \param[in] nsp pointer to a \ref cs_navsto_param_t struct.
325 * \param[in] p_cell value of the pressure inside each cell
326 * \param[in, out] p_face value of the pressure at each face
327 */
328 /*----------------------------------------------------------------------------*/
329
330 void
331 cs_cdofb_navsto_compute_face_pressure(const cs_mesh_t *mesh,
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 /*----------------------------------------------------------------------------*/
340 /*!
341 * \brief Add the grad-div part to the local matrix (i.e. for the current
342 * cell)
343 *
344 * \param[in] n_fc local number of faces for the current cell
345 * \param[in] zeta scalar coefficient for the grad-div operator
346 * \param[in] div divergence
347 * \param[in, out] mat local system matrix to update
348 */
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 /*----------------------------------------------------------------------------*/
358 /*!
359 * \brief Initialize the pressure values
360 *
361 * \param[in] nsp pointer to a \ref cs_navsto_param_t structure
362 * \param[in] quant pointer to a \ref cs_cdo_quantities_t structure
363 * \param[in] ts pointer to a \ref cs_time_step_t structure
364 * \param[in, out] pr pointer to the pressure \ref cs_field_t structure
365 */
366 /*----------------------------------------------------------------------------*/
367
368 void
369 cs_cdofb_navsto_init_pressure(const cs_navsto_param_t *nsp,
370 const cs_cdo_quantities_t *quant,
371 const cs_time_step_t *ts,
372 cs_field_t *pr);
373
374 /*----------------------------------------------------------------------------*/
375 /*!
376 * \brief Initialize the pressure values when the pressure is defined at
377 * faces
378 *
379 * \param[in] nsp pointer to a \ref cs_navsto_param_t structure
380 * \param[in] connect pointer to a \ref cs_cdo_connect_t structure
381 * \param[in] ts pointer to a \ref cs_time_step_t structure
382 * \param[in, out] pr_f pointer to the pressure values at faces
383 */
384 /*----------------------------------------------------------------------------*/
385
386 void
387 cs_cdofb_navsto_init_face_pressure(const cs_navsto_param_t *nsp,
388 const cs_cdo_connect_t *connect,
389 const cs_time_step_t *ts,
390 cs_real_t *pr_f);
391
392 /*----------------------------------------------------------------------------*/
393 /*!
394 * \brief Update the pressure field in order to get a field with a mean-value
395 * equal to the reference value
396 *
397 * \param[in] nsp pointer to a cs_navsto_param_t structure
398 * \param[in] quant pointer to a cs_cdo_quantities_t structure
399 * \param[in, out] values pressure field values
400 */
401 /*----------------------------------------------------------------------------*/
402
403 void
404 cs_cdofb_navsto_rescale_pressure_to_ref(const cs_navsto_param_t *nsp,
405 const cs_cdo_quantities_t *quant,
406 cs_real_t values[]);
407
408 /*----------------------------------------------------------------------------*/
409 /*!
410 * \brief Update the pressure field in order to get a field with a zero-mean
411 * average
412 *
413 * \param[in] quant pointer to a cs_cdo_quantities_t structure
414 * \param[in, out] values pressure field values
415 */
416 /*----------------------------------------------------------------------------*/
417
418 void
419 cs_cdofb_navsto_set_zero_mean_pressure(const cs_cdo_quantities_t *quant,
420 cs_real_t values[]);
421
422 /*----------------------------------------------------------------------------*/
423 /*!
424 * \brief Perform extra-operation related to Fb schemes when solving
425 * Navier-Stokes. Computation of the following quantities according to
426 * post-processing flags beeing activated.
427 * - The mass flux accross the boundaries.
428 * - The global mass in the computational domain
429 * - The norm of the velocity divergence
430 * - the cellwise mass flux balance
431 * - the kinetic energy
432 * - the velocity gradient
433 * - the pressure gradient
434 * - the vorticity
435 * - the helicity
436 * - the enstrophy
437 * - the stream function
438 *
439 * \param[in] nsp pointer to a \ref cs_navsto_param_t struct.
440 * \param[in] mesh pointer to a cs_mesh_t structure
441 * \param[in] quant pointer to a \ref cs_cdo_quantities_t struct.
442 * \param[in] connect pointer to a \ref cs_cdo_connect_t struct.
443 * \param[in] ts pointer to a \ref cs_time_step_t struct.
444 * \param[in,out] time_plotter pointer to a \ref cs_time_plot_t struct.
445 * \param[in] adv_field pointer to a \ref cs_adv_field_t struct.
446 * \param[in] mass_flux scalar-valued mass flux for each face
447 * \param[in] p_cell scalar-valued pressure in each cell
448 * \param[in] u_cell vector-valued velocity in each cell
449 * \param[in] u_face vector-valued velocity on each face
450 */
451 /*----------------------------------------------------------------------------*/
452
453 void
454 cs_cdofb_navsto_extra_op(const cs_navsto_param_t *nsp,
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 /*----------------------------------------------------------------------------*/
467 /*!
468 * \brief Take into account a Dirichlet BCs on the three velocity components.
469 * For instance for a velocity inlet boundary or a wall
470 * Handle the velocity-block in the global algebraic system in case of
471 * an algebraic technique.
472 * This prototype matches the function pointer cs_cdo_apply_boundary_t
473 *
474 * \param[in] f face id in the cell mesh numbering
475 * \param[in] eqp pointer to a \ref cs_equation_param_t struct.
476 * \param[in] cm pointer to a \ref cs_cell_mesh_t structure
477 * \param[in] pty pointer to a \ref cs_property_data_t structure
478 * \param[in, out] cb pointer to a \ref cs_cell_builder_t structure
479 * \param[in, out] csys structure storing the cellwise system
480 */
481 /*----------------------------------------------------------------------------*/
482
483 void
484 cs_cdofb_block_dirichlet_alge(short int f,
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 /*----------------------------------------------------------------------------*/
492 /*!
493 * \brief Take into account a Dirichlet BCs on the three velocity components.
494 * For instance for a velocity inlet boundary or a wall
495 * Handle the velocity-block in the global algebraic system in case of
496 * a penalization technique (with a large coefficient).
497 * One assumes that static condensation has been performed and that
498 * the velocity-block has size 3*n_fc
499 * This prototype matches the function pointer cs_cdo_apply_boundary_t
500 *
501 * \param[in] f face id in the cell mesh numbering
502 * \param[in] eqp pointer to a \ref cs_equation_param_t struct.
503 * \param[in] cm pointer to a \ref cs_cell_mesh_t structure
504 * \param[in] pty pointer to a \ref cs_property_data_t structure
505 * \param[in, out] cb pointer to a \ref cs_cell_builder_t structure
506 * \param[in, out] csys structure storing the cellwise system
507 */
508 /*----------------------------------------------------------------------------*/
509
510 void
511 cs_cdofb_block_dirichlet_pena(short int f,
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 /*----------------------------------------------------------------------------*/
519 /*!
520 * \brief Take into account a Dirichlet BCs on the three velocity components.
521 * For instance for a velocity inlet boundary or a wall
522 * Handle the velocity-block in the global algebraic system in case of
523 * a weak penalization technique (Nitsche).
524 * One assumes that static condensation has not been performed yet and
525 * that the velocity-block has size 3*(n_fc + 1)
526 * This prototype matches the function pointer cs_cdo_apply_boundary_t
527 *
528 * \param[in] fb face id in the cell mesh numbering
529 * \param[in] eqp pointer to a \ref cs_equation_param_t struct.
530 * \param[in] cm pointer to a \ref cs_cell_mesh_t structure
531 * \param[in] pty pointer to a \ref cs_property_data_t structure
532 * \param[in, out] cb pointer to a \ref cs_cell_builder_t structure
533 * \param[in, out] csys structure storing the cellwise system
534 */
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 /*----------------------------------------------------------------------------*/
546 /*!
547 * \brief Take into account a Dirichlet BCs on the three velocity components.
548 * For instance for a velocity inlet boundary or a wall
549 * Handle the velocity-block in the global algebraic system in case of
550 * a weak penalization technique (symmetrized Nitsche).
551 * One assumes that static condensation has not been performed yet and
552 * that the velocity-block has size 3*(n_fc + 1)
553 * This prototype matches the function pointer cs_cdo_apply_boundary_t
554 *
555 * \param[in] fb face id in the cell mesh numbering
556 * \param[in] eqp pointer to a \ref cs_equation_param_t struct.
557 * \param[in] cm pointer to a \ref cs_cell_mesh_t structure
558 * \param[in] pty pointer to a \ref cs_property_data_t structure
559 * \param[in, out] cb pointer to a \ref cs_cell_builder_t structure
560 * \param[in, out] csys structure storing the cellwise system
561 */
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 /*----------------------------------------------------------------------------*/
573 /*!
574 * \brief Take into account a boundary defined as 'symmetry' (treated as a
575 * sliding BCs on the three velocity components.)
576 * A weak penalization technique (symmetrized Nitsche) is used. One
577 * assumes that static condensation has not been performed yet and that
578 * the velocity-block has (n_fc + 1) blocks of size 3x3.
579 * This prototype matches the function pointer cs_cdo_apply_boundary_t
580 *
581 * \param[in] fb face id in the cell mesh numbering
582 * \param[in] eqp pointer to a \ref cs_equation_param_t struct.
583 * \param[in] cm pointer to a \ref cs_cell_mesh_t structure
584 * \param[in] pty pointer to a \ref cs_property_data_t structure
585 * \param[in, out] cb pointer to a \ref cs_cell_builder_t structure
586 * \param[in, out] csys structure storing the cellwise system
587 */
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 /*----------------------------------------------------------------------------*/
599 /*!
600 * \brief Take into account a wall BCs by a weak enforcement using Nitsche
601 * technique plus a symmetric treatment.
602 * This prototype matches the function pointer cs_cdo_apply_boundary_t
603 *
604 * \param[in] fb face id in the cell mesh numbering
605 * \param[in] eqp pointer to a \ref cs_equation_param_t struct.
606 * \param[in] cm pointer to a \ref cs_cell_mesh_t structure
607 * \param[in] pty pointer to a \ref cs_property_data_t structure
608 * \param[in, out] cb pointer to a \ref cs_cell_builder_t structure
609 * \param[in, out] csys structure storing the cellwise system
610 */
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 /*----------------------------------------------------------------------------*/
622 /*!
623 * \brief Test if one has to do one more non-linear iteration.
624 * Test if performed on the relative norm on the increment between
625 * two iterations
626 *
627 * \param[in] nl_algo_type type of non-linear algorithm
628 * \param[in] pre_iterate previous state of the mass flux iterate
629 * \param[in] cur_iterate current state of the mass flux iterate
630 * \param[in, out] algo pointer to a cs_iter_algo_t structure
631 *
632 * \return the convergence state
633 */
634 /*----------------------------------------------------------------------------*/
635
636 cs_sles_convergence_state_t
637 cs_cdofb_navsto_nl_algo_cvg(cs_param_nl_algo_t nl_algo_type,
638 const cs_real_t *pre_iterate,
639 cs_real_t *cur_iterate,
640 cs_iter_algo_t *algo);
641
642 /*----------------------------------------------------------------------------*/
643 /*!
644 * \brief Set the function pointer computing the source term in the momentum
645 * equation related to the gravity effect (hydrostatic pressure or the
646 * Boussinesq approximation)
647 *
648 * \param[in] nsp set of parameters for the Navier-Stokes system
649 * \param[out] p_func way to compute the gravity effect
650 */
651 /*----------------------------------------------------------------------------*/
652
653 void
654 cs_cdofb_navsto_set_gravity_func(const cs_navsto_param_t *nsp,
655 cs_cdofb_navsto_source_t **p_func);
656
657 /*----------------------------------------------------------------------------*/
658 /*!
659 * \brief Take into account the gravity effects.
660 * Compute and add the source term to the local RHS.
661 * This is a special treatment since of face DoFs are involved
662 * contrary to the standard case where only the cell DoFs is involved.
663 *
664 * \param[in] nsp set of parameters to handle the Navier-Stokes system
665 * \param[in] cm pointer to a cs_cell_mesh_t structure
666 * \param[in] nsb pointer to a builder structure for the NavSto system
667 * \param[in, out] csys pointer to a cs_cell_sys_t structure
668 */
669 /*----------------------------------------------------------------------------*/
670
671 void
672 cs_cdofb_navsto_gravity_term(const cs_navsto_param_t *nsp,
673 const cs_cell_mesh_t *cm,
674 const cs_cdofb_navsto_builder_t *nsb,
675 cs_cell_sys_t *csys);
676
677 /*----------------------------------------------------------------------------*/
678 /*!
679 * \brief Take into account the buoyancy force with the Boussinesq approx.
680 * Compute and add the source term to the local RHS.
681 * This is the standard case where the face DoFs are used for the
682 * constant part rho0 . g[] and only the cell DoFs are involved for the
683 * remaining part (the Boussinesq approximation).
684 *
685 * \param[in] nsp set of parameters to handle the Navier-Stokes system
686 * \param[in] cm pointer to a cs_cell_mesh_t structure
687 * \param[in] nsb pointer to a builder structure for the NavSto system
688 * \param[in, out] csys pointer to a cs_cell_sys_t structure
689 */
690 /*----------------------------------------------------------------------------*/
691
692 void
693 cs_cdofb_navsto_boussinesq_at_cell(const cs_navsto_param_t *nsp,
694 const cs_cell_mesh_t *cm,
695 const cs_cdofb_navsto_builder_t *nsb,
696 cs_cell_sys_t *csys);
697
698 /*----------------------------------------------------------------------------*/
699 /*!
700 * \brief Take into account the buoyancy force with the Boussinesq approx.
701 * Compute and add the source term to the local RHS.
702 * This way to compute the Boussinesq approximation relies only on DoFs
703 * at faces. This should enable to keep a stable (no velocity) in case
704 * of a stratified configuration.
705 *
706 * \param[in] nsp set of parameters to handle the Navier-Stokes system
707 * \param[in] cm pointer to a cs_cell_mesh_t structure
708 * \param[in] nsb pointer to a builder structure for the NavSto system
709 * \param[in, out] csys pointer to a cs_cell_sys_t structure
710 */
711 /*----------------------------------------------------------------------------*/
712
713 void
714 cs_cdofb_navsto_boussinesq_at_face(const cs_navsto_param_t *nsp,
715 const cs_cell_mesh_t *cm,
716 const cs_cdofb_navsto_builder_t *nsb,
717 cs_cell_sys_t *csys);
718
719 /*----------------------------------------------------------------------------*/
720 /*!
721 * \brief Take into account the buoyancy force with the Boussinesq approx.
722 * Compute and add the source term to the local RHS.
723 * This way to compute the Boussinesq approximation relies only on DoFs
724 * at faces. This should enable to keep a stable (no velocity) in case
725 * of a stratified configuration.
726 *
727 * \param[in] nsp set of parameters to handle the Navier-Stokes system
728 * \param[in] cm pointer to a cs_cell_mesh_t structure
729 * \param[in] nsb pointer to a builder structure for the NavSto system
730 * \param[in, out] csys pointer to a cs_cell_sys_t structure
731 */
732 /*----------------------------------------------------------------------------*/
733
734 void
735 cs_cdofb_navsto_boussinesq_by_part(const cs_navsto_param_t *nsp,
736 const cs_cell_mesh_t *cm,
737 const cs_cdofb_navsto_builder_t *nsb,
738 cs_cell_sys_t *csys);
739
740 /*----------------------------------------------------------------------------*/
741 /*!
742 * \brief Get the source term for computing the stream function.
743 * This relies on the prototype associated to the generic function
744 * pointer \ref cs_dof_function_t
745 *
746 * \param[in] n_elts number of elements to consider
747 * \param[in] elt_ids list of elements ids
748 * \param[in] dense_output perform an indirection in retval or not
749 * \param[in] input NULL or pointer to a structure cast on-the-fly
750 * \param[in, out] retval result of the function. Must be allocated.
751 */
752 /*----------------------------------------------------------------------------*/
753
754 void
755 cs_cdofb_navsto_stream_source_term(cs_lnum_t n_elts,
756 const cs_lnum_t *elt_ids,
757 bool dense_output,
758 void *input,
759 cs_real_t *retval);
760
761 /*----------------------------------------------------------------------------*/
762
763 END_C_DECLS
764
765 #endif /* __CS_CDOFB_NAVSTO_H__ */
766