1 #ifndef __CS_EQUATION_COMMON_H__
2 #define __CS_EQUATION_COMMON_H__
3
4 /*============================================================================
5 * Functions to handle common equation features for building algebraic system
6 * in CDO schemes
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 * Local headers
31 *----------------------------------------------------------------------------*/
32
33 #include "cs_cdo_bc.h"
34 #include "cs_cdo_connect.h"
35 #include "cs_cdo_local.h"
36 #include "cs_cdo_quantities.h"
37 #include "cs_enforcement.h"
38 #include "cs_equation_param.h"
39 #include "cs_flag.h"
40 #include "cs_matrix.h"
41 #include "cs_time_step.h"
42 #include "cs_timer.h"
43 #include "cs_sles.h"
44 #include "cs_source_term.h"
45
46 /*----------------------------------------------------------------------------*/
47
48 BEGIN_C_DECLS
49
50 /*============================================================================
51 * Macro definitions
52 *============================================================================*/
53
54 /* Strategy of synchronization of values shared across several cells
55 * This applies to vertices, edges and faces
56 *
57 * CS_EQUATION_SYNC_ZERO_VALUE
58 * If zero is a possible value then set this value, otherwise one takes
59 * the mean-value
60 *
61 * CS_EQUATION_SYNC_MEAN_VALUE
62 * Compute the mean-value across values to set
63 *
64 */
65
66 #define CS_EQUATION_SYNC_ZERO_VALUE 1
67 #define CS_EQUATION_SYNC_MEAN_VALUE 2
68
69 /*============================================================================
70 * Type definitions
71 *============================================================================*/
72
73 typedef struct _equation_builder_t cs_equation_builder_t;
74
75 /*----------------------------------------------------------------------------*/
76 /*!
77 * \brief Generic function prototype for a hook during the cellwise building
78 * of the linear system
79 * Enable an advanced user to get a fine control of the discretization
80 *
81 * \param[in] eqp pointer to a cs_equation_param_t structure
82 * \param[in] eqb pointer to a cs_equation_builder_t structure
83 * \param[in] eqc context to cast for this discretization
84 * \param[in] cm pointer to a cellwise view of the mesh
85 * \param[in, out] mass_hodge pointer to a cs_hodge_t structure (mass matrix)
86 * \param[in, out] diff_hodge pointer to a cs_hodge_t structure (diffusion)
87 * \param[in, out] csys pointer to a cellwise view of the system
88 * \param[in, out] cb pointer to a cellwise builder
89 */
90 /*----------------------------------------------------------------------------*/
91
92 typedef void
93 (cs_equation_user_hook_t)(const cs_equation_param_t *eqp,
94 const cs_equation_builder_t *eqb,
95 const void *eq_context,
96 const cs_cell_mesh_t *cm,
97 cs_hodge_t *mass_hodge,
98 cs_hodge_t *diff_hodge,
99 cs_cell_sys_t *csys,
100 cs_cell_builder_t *cb);
101
102 /*! \struct cs_equation_builder_t
103 * \brief Store common elements used when building an algebraic system
104 * related to an equation
105 */
106
107 struct _equation_builder_t {
108
109 bool init_step; /*!< true if this is the initialization step */
110
111 /*!
112 * @name Flags to know what to build and how to build such terms
113 * @{
114 */
115
116 cs_eflag_t msh_flag; /*!< Information related to what to build in a
117 * \ref cs_cell_mesh_t structure for a generic
118 * cell */
119 cs_eflag_t bd_msh_flag; /*!< Information related to what to build in a
120 * \ref cs_cell_mesh_t structure for a cell close
121 * to the boundary */
122 cs_eflag_t st_msh_flag; /*!< Information related to what to build in a
123 * \ref cs_cell_mesh_t structure when only the
124 * source term has to be built */
125 cs_flag_t sys_flag; /*!< Information related to the sytem */
126
127 /*!
128 * @}
129 * @name Metadata related to associated physical properties
130 * @{
131 */
132
133 bool diff_pty_uniform; /*!< Is diffusion property uniform ? */
134 bool curlcurl_pty_uniform; /*!< Is curl-curl property uniform ? */
135 bool graddiv_pty_uniform; /*!< Is grad-div property uniform ? */
136 bool time_pty_uniform; /*!< Is time property uniform ? */
137 bool reac_pty_uniform[CS_CDO_N_MAX_REACTIONS]; /*!< Is each reaction
138 * property uniform ? */
139
140 /*!
141 * @}
142 * @name Source terms
143 * @{
144 */
145
146 cs_mask_t *source_mask; /*!< NULL if no source term or one source term
147 * is defined. Allocated to n_cells in order to
148 * know in each cell which source term has to be
149 * computed */
150
151 /*! \var compute_source
152 * Pointer to functions which compute the value of the source term
153 */
154
155 cs_source_term_cellwise_t *compute_source[CS_N_MAX_SOURCE_TERMS];
156
157 /*!
158 * @}
159 * @name Enforcement of degrees of freedom (DoFs)
160 * @{
161 */
162
163 cs_real_t *enforced_values;
164
165 /*!
166 * @}
167 * @name User hook
168 * @{
169 *
170 * \var user_hook_context
171 * Pointer to a shared structure (the lifecycle of this structure is not
172 * managed by the current cs_equation_builder_t structure)
173 *
174 * \var user_hook_function
175 * Function pointer associated to a predefined prototype
176 * This function enables a user to modify the cellwise system (matrix and
177 * rhs) before applying the time scheme, the static condensation if needed or
178 * the strong/penalized enforcement of boundary conditions.
179 * This is useful to add a term in the equation like an advanced source term
180 * without the need to allocate an array and with an access to the local
181 * structure such as the local cell mesh, the cell builder and the high-level
182 * structures related to an equation
183 */
184
185 void *user_hook_context;
186 cs_equation_user_hook_t *user_hook_function;
187
188 /*!
189 * @}
190 * @name Boundary conditions
191 * @{
192 *
193 * \var face_bc
194 * face_bc should not change during the simulation.
195 * The case of a definition of the BCs which changes of type during the
196 * simulation is possible but not implemented.
197 * You just have to call the initialization step each time the type of BCs
198 * is modified to define an updated \ref cs_cdo_bc_face_t structure.
199 */
200
201 cs_cdo_bc_face_t *face_bc; /*!< Information about boundary conditions
202 applied to faces */
203
204 cs_real_t *dir_values; /*!< Array storing the Dirichlet values at
205 DoFs */
206
207 /*!
208 * @}
209 * @name Performance monitoring
210 * @{
211 *
212 * Monitoring the efficiency of the algorithm used to manipulate/build
213 * an equation.
214 */
215
216 cs_timer_counter_t tcb; /*!< Cumulated elapsed time for building the
217 * current system */
218 cs_timer_counter_t tcs; /*!< Cumulated elapsed time for solving the
219 * current system */
220 cs_timer_counter_t tce; /*!< Cumulated elapsed time for computing
221 * all extra operations (post, balance,
222 * fluxes...) */
223
224 /*! @} */
225
226 };
227
228 /*
229 * Structure used to store information generated during the analysis
230 * of the balance of each term of an equation
231 */
232
233 typedef struct {
234
235 /* where balance is computed: primal vertices or primal cells */
236
237 cs_flag_t location;
238 cs_lnum_t size;
239 cs_real_t *balance;
240
241 /* Balance for each main term */
242
243 cs_real_t *unsteady_term;
244 cs_real_t *reaction_term;
245 cs_real_t *diffusion_term;
246 cs_real_t *advection_term;
247 cs_real_t *source_term;
248 cs_real_t *boundary_term;
249
250 } cs_equation_balance_t;
251
252 /*============================================================================
253 * Inline public function prototypes
254 *============================================================================*/
255
256 /*----------------------------------------------------------------------------*/
257 /*!
258 * \brief Retrieve the flag to give for building a cs_cell_mesh_t structure
259 *
260 * \param[in] cell_flag flag related to the current cell
261 * \param[in] eqb pointer to a cs_equation_builder_t structure
262 *
263 * \return the flag to set for the current cell
264 */
265 /*----------------------------------------------------------------------------*/
266
267 static inline cs_eflag_t
cs_equation_cell_mesh_flag(cs_flag_t cell_flag,const cs_equation_builder_t * eqb)268 cs_equation_cell_mesh_flag(cs_flag_t cell_flag,
269 const cs_equation_builder_t *eqb)
270 {
271 cs_eflag_t _flag = eqb->msh_flag | eqb->st_msh_flag;
272
273 if (cell_flag & CS_FLAG_BOUNDARY_CELL_BY_FACE)
274 _flag |= eqb->bd_msh_flag;
275
276 return _flag;
277 }
278
279 /*============================================================================
280 * Public function prototypes
281 *============================================================================*/
282
283 /*----------------------------------------------------------------------------*/
284 /*!
285 * \brief Allocate a pointer to a buffer of size at least the n_cells for
286 * managing temporary usage of memory when dealing with equations
287 * The size of the temporary buffer can be bigger according to the
288 * numerical settings
289 * Set also shared pointers from the main domain members
290 *
291 * \param[in] connect pointer to a cs_cdo_connect_t structure
292 * \param[in] quant pointer to additional mesh quantities struct.
293 * \param[in] time_step pointer to a time step structure
294 * \param[in] eb_flag metadata for Edge-based schemes
295 * \param[in] fb_flag metadata for Face-based schemes
296 * \param[in] vb_flag metadata for Vertex-based schemes
297 * \param[in] vcb_flag metadata for Vertex+Cell-basde schemes
298 * \param[in] hho_flag metadata for HHO schemes
299 */
300 /*----------------------------------------------------------------------------*/
301
302 void
303 cs_equation_common_init(const cs_cdo_connect_t *connect,
304 const cs_cdo_quantities_t *quant,
305 const cs_time_step_t *time_step,
306 cs_flag_t eb_flag,
307 cs_flag_t fb_flag,
308 cs_flag_t vb_flag,
309 cs_flag_t vcb_flag,
310 cs_flag_t hho_flag);
311
312 /*----------------------------------------------------------------------------*/
313 /*!
314 * \brief Allocate a pointer to a buffer of size at least the 2*n_cells for
315 * managing temporary usage of memory when dealing with equations
316 * Call specific structure allocation related to a numerical scheme
317 * according to the scheme flag
318 * The size of the temporary buffer can be bigger according to the
319 * numerical settings
320 */
321 /*----------------------------------------------------------------------------*/
322
323 void
324 cs_equation_common_finalize(void);
325
326 /*----------------------------------------------------------------------------*/
327 /*!
328 * \brief Allocate a new structure to handle the building of algebraic system
329 * related to an cs_equation_t structure
330 *
331 * \param[in] eqp pointer to a cs_equation_param_t structure
332 * \param[in] mesh pointer to a cs_mesh_t structure
333 *
334 * \return a pointer to a new allocated cs_equation_builder_t structure
335 */
336 /*----------------------------------------------------------------------------*/
337
338 cs_equation_builder_t *
339 cs_equation_builder_init(const cs_equation_param_t *eqp,
340 const cs_mesh_t *mesh);
341
342 /*----------------------------------------------------------------------------*/
343 /*!
344 * \brief Free a cs_equation_builder_t structure
345 *
346 * \param[in, out] p_builder pointer of pointer to the cs_equation_builder_t
347 * structure to free
348 */
349 /*----------------------------------------------------------------------------*/
350
351 void
352 cs_equation_builder_free(cs_equation_builder_t **p_builder);
353
354 /*----------------------------------------------------------------------------*/
355 /*!
356 * \brief Free some members of a cs_equation_builder_t structure
357 *
358 * \param[in, out] eqb pointer to the cs_equation_builder_t structure
359 */
360 /*----------------------------------------------------------------------------*/
361
362 void
363 cs_equation_builder_reset(cs_equation_builder_t *eqb);
364
365 /*----------------------------------------------------------------------------*/
366 /*!
367 * \brief Compute the value of the renormalization coefficient for the
368 * the residual norm of the linear system
369 * A pre-computation arising from cellwise contribution may have
370 * been done before this call according to the requested type of
371 * renormalization.
372 *
373 * \param[in] type type of renormalization
374 * \param[in] rhs_size size of the rhs array
375 * \param[in] rhs array related to the right-hand side
376 * \param[in, out] normalization value of the residual normalization
377 */
378 /*----------------------------------------------------------------------------*/
379
380 void
381 cs_equation_sync_rhs_normalization(cs_param_resnorm_type_t type,
382 cs_lnum_t rhs_size,
383 const cs_real_t rhs[],
384 double *normalization);
385
386 /*----------------------------------------------------------------------------*/
387 /*!
388 * \brief Prepare a linear system and synchronize buffers to handle parallelism.
389 * Transfer a mesh-based description of arrays x0 and rhs into an
390 * algebraic description for the linear system in x and b.
391 *
392 * \param[in] stride stride to apply to the range set operations
393 * \param[in] x_size size of the vector unknowns (scatter view)
394 * \param[in] matrix pointer to a cs_matrix_t structure
395 * \param[in] rset pointer to a range set structure
396 * \param[in] rhs_redux do or not a parallel sum reduction on the RHS
397 * \param[in, out] x array of unknowns (in: initial guess)
398 * \param[in, out] b right-hand side
399 */
400 /*----------------------------------------------------------------------------*/
401
402 void
403 cs_equation_prepare_system(int stride,
404 cs_lnum_t x_size,
405 const cs_matrix_t *matrix,
406 const cs_range_set_t *rset,
407 bool rhs_redux,
408 cs_real_t *x,
409 cs_real_t *b);
410
411 /*----------------------------------------------------------------------------*/
412 /*!
413 * \brief Solve a linear system arising with scalar-valued cell-based DoFs
414 * No rotation is taken into account when synchronizing halo.
415 *
416 * \param[in] n_dofs local number of DoFs
417 * \param[in] slesp pointer to a cs_param_sles_t structure
418 * \param[in] matrix pointer to a cs_matrix_t structure
419 * \param[in] normalization value used for the residual normalization
420 * \param[in, out] sles pointer to a cs_sles_t structure
421 * \param[in, out] x solution of the linear system (in: initial guess)
422 * \param[in, out] b right-hand side (scatter/gather if needed)
423 *
424 * \return the number of iterations of the linear solver
425 */
426 /*----------------------------------------------------------------------------*/
427
428 int
429 cs_equation_solve_scalar_cell_system(cs_lnum_t n_dofs,
430 const cs_param_sles_t *slesp,
431 const cs_matrix_t *matrix,
432 cs_real_t normalization,
433 cs_sles_t *sles,
434 cs_real_t *x,
435 cs_real_t *b);
436
437 /*----------------------------------------------------------------------------*/
438 /*!
439 * \brief Solve a linear system arising from CDO schemes with scalar-valued
440 * degrees of freedom
441 *
442 * \param[in] n_scatter_dofs local number of DoFs (may be != n_gather_elts)
443 * \param[in] slesp pointer to a cs_param_sles_t structure
444 * \param[in] matrix pointer to a cs_matrix_t structure
445 * \param[in] rs pointer to a cs_range_set_t structure
446 * \param[in] normalization value used for the residual normalization
447 * \param[in] rhs_redux do or not a parallel sum reduction on the RHS
448 * \param[in, out] sles pointer to a cs_sles_t structure
449 * \param[in, out] x solution of the linear system (in: initial guess)
450 * \param[in, out] b right-hand side (scatter/gather if needed)
451 *
452 * \return the number of iterations of the linear solver
453 */
454 /*----------------------------------------------------------------------------*/
455
456 int
457 cs_equation_solve_scalar_system(cs_lnum_t n_scatter_dofs,
458 const cs_param_sles_t *slesp,
459 const cs_matrix_t *matrix,
460 const cs_range_set_t *rset,
461 cs_real_t normalization,
462 bool rhs_redux,
463 cs_sles_t *sles,
464 cs_real_t *x,
465 cs_real_t *b);
466
467 /*----------------------------------------------------------------------------*/
468 /*!
469 * \brief Print a message in the performance output file related to the
470 * monitoring of equation
471 *
472 * \param[in] eqname pointer to the name of the current equation
473 * \param[in] eqb pointer to a cs_equation_builder_t structure
474 */
475 /*----------------------------------------------------------------------------*/
476
477 void
478 cs_equation_write_monitoring(const char *eqname,
479 const cs_equation_builder_t *eqb);
480
481 /*----------------------------------------------------------------------------*/
482 /*!
483 * \brief Initialize all reaction properties. This function is shared across
484 * all CDO schemes. The \ref cs_cell_builder_t structure stores the
485 * computed property values.
486 *
487 * \param[in] eqp pointer to a cs_equation_param_t structure
488 * \param[in] eqb pointer to a cs_equation_builder_t structure
489 * \param[in] t_eval time at which one performs the evaluation
490 * \param[in, out] cb pointer to a cs_cell_builder_t structure
491 */
492 /*----------------------------------------------------------------------------*/
493
494 void
495 cs_equation_init_reaction_properties(const cs_equation_param_t *eqp,
496 const cs_equation_builder_t *eqb,
497 cs_real_t t_eval,
498 cs_cell_builder_t *cb);
499
500 /*----------------------------------------------------------------------------*/
501 /*!
502 * \brief Initialize all reaction properties. This function is shared across
503 * all CDO schemes. The \ref cs_cell_builder_t structure stores the
504 * computed property values.
505 * If the property is uniform, a first call to the function
506 * \ref cs_equation_init_reaction_properties or to the function
507 * \ref cs_equation_init_properties has to be done before the
508 * loop on cells
509 *
510 * \param[in] eqp pointer to a cs_equation_param_t structure
511 * \param[in] eqb pointer to a cs_equation_builder_t structure
512 * \param[in] cm pointer to a \ref cs_cell_mesh_t structure
513 * \param[in, out] cb pointer to a \ref cs_cell_builder_t structure
514 */
515 /*----------------------------------------------------------------------------*/
516
517 void
518 cs_equation_set_reaction_properties_cw(const cs_equation_param_t *eqp,
519 const cs_equation_builder_t *eqb,
520 const cs_cell_mesh_t *cm,
521 cs_cell_builder_t *cb);
522
523 /*----------------------------------------------------------------------------*/
524 /*!
525 * \brief Initialize all properties potentially useful to build the algebraic
526 * system. This function is shared across all CDO schemes.
527 * The \ref cs_cell_builder_t structure stores property values related
528 * to the reaction term, unsteady term and grad-div term.
529 *
530 * \param[in] eqp pointer to a cs_equation_param_t structure
531 * \param[in] eqb pointer to a cs_equation_builder_t structure
532 * \param[in, out] diffusion_hodge pointer to the diffusion hodge structure
533 * \param[in, out] cb pointer to a cs_cell_builder_t structure
534 */
535 /*----------------------------------------------------------------------------*/
536
537 void
538 cs_equation_init_properties(const cs_equation_param_t *eqp,
539 const cs_equation_builder_t *eqb,
540 cs_hodge_t *diffusion_hodge,
541 cs_cell_builder_t *cb);
542
543 /*----------------------------------------------------------------------------*/
544 /*!
545 * \brief Build the list of degrees of freedom (DoFs) related to an internal
546 * enforcement. If set to NULL, the array dof_ids (storing the
547 * indirection) is allocated to n_x.
548 *
549 * \param[in] n_x number of entities where DoFs are defined
550 * \param[in] c2x cell --> x connectivity
551 * \param[in] ifs pointer to a fvm_interface_set_t structure
552 * \param[in] eqp set of parameters related to an equation
553 * \param[in, out] p_dof_ids double pointer on DoF ids subject to enforcement
554 */
555 /*----------------------------------------------------------------------------*/
556
557 void
558 cs_equation_build_dof_enforcement(cs_lnum_t n_x,
559 const cs_adjacency_t *c2x,
560 const cs_interface_set_t *ifs,
561 const cs_equation_param_t *eqp,
562 cs_lnum_t *p_dof_ids[]);
563
564 /*----------------------------------------------------------------------------*/
565 /*!
566 * \brief Take into account the enforcement of internal DoFs. Apply an
567 * algebraic manipulation. Update members of the cs_cell_sys_t
568 * structure related to the internal enforcement.
569 *
570 * | | | | | | | | | |
571 * | Aii | Aie | | Aii | 0 | |bi| |bi -Aid.x_enf|
572 * |------------| --> |------------| and |--| --> |-------------|
573 * | | | | | | | | | |
574 * | Aei | Aee | | 0 | Id | |be| | x_enf |
575 *
576 * where x_enf is the value of the enforcement for the selected internal DoFs
577 *
578 * \param[in] eqb pointer to a cs_equation_builder_t structure
579 * \param[in, out] cb pointer to a cs_cell_builder_t structure
580 * \param[in, out] csys structure storing the cell-wise system
581 */
582 /*----------------------------------------------------------------------------*/
583
584 void
585 cs_equation_enforced_internal_dofs(const cs_equation_builder_t *eqb,
586 cs_cell_builder_t *cb,
587 cs_cell_sys_t *csys);
588
589 /*----------------------------------------------------------------------------*/
590 /*!
591 * \brief Take into account the enforcement of internal DoFs. Case of matrices
592 * defined by blocks. Apply an algebraic manipulation. Update members
593 * of the cs_cell_sys_t structure related to the internal enforcement.
594 *
595 * | | | | | | | | | |
596 * | Aii | Aie | | Aii | 0 | |bi| |bi -Aid.x_enf|
597 * |------------| --> |------------| and |--| --> |-------------|
598 * | | | | | | | | | |
599 * | Aei | Aee | | 0 | Id | |be| | x_enf |
600 *
601 * where x_enf is the value of the enforcement for the selected internal DoFs
602 *
603 * \param[in] eqb pointer to a cs_equation_builder_t structure
604 * \param[in, out] cb pointer to a cs_cell_builder_t structure
605 * \param[in, out] csys structure storing the cell-wise system
606 */
607 /*----------------------------------------------------------------------------*/
608
609 void
610 cs_equation_enforced_internal_block_dofs(const cs_equation_builder_t *eqb,
611 cs_cell_builder_t *cb,
612 cs_cell_sys_t *csys);
613
614 /*----------------------------------------------------------------------------*/
615 /*!
616 * \brief Retrieve a pointer to a buffer of size at least the 2*n_cells
617 * The size of the temporary buffer can be bigger according to the
618 * numerical settings
619 *
620 * \return a pointer to an array of double
621 */
622 /*----------------------------------------------------------------------------*/
623
624 cs_real_t *
625 cs_equation_get_tmpbuf(void);
626
627 /*----------------------------------------------------------------------------*/
628 /*!
629 * \brief Get the allocation size of the temporary buffer
630 *
631 * \return the size of the temporary buffer
632 */
633 /*----------------------------------------------------------------------------*/
634
635 size_t
636 cs_equation_get_tmpbuf_size(void);
637
638 /*----------------------------------------------------------------------------*/
639 /*!
640 * \brief Allocate a cs_equation_balance_t structure
641 *
642 * \param[in] location where the balance is performed
643 * \param[in] size size of arrays in the structure
644 *
645 * \return a pointer to the new allocated structure
646 */
647 /*----------------------------------------------------------------------------*/
648
649 cs_equation_balance_t *
650 cs_equation_balance_create(cs_flag_t location,
651 cs_lnum_t size);
652
653 /*----------------------------------------------------------------------------*/
654 /*!
655 * \brief Reset a cs_equation_balance_t structure
656 *
657 * \param[in, out] b pointer to a cs_equation_balance_t to reset
658 */
659 /*----------------------------------------------------------------------------*/
660
661 void
662 cs_equation_balance_reset(cs_equation_balance_t *b);
663
664 /*----------------------------------------------------------------------------*/
665 /*!
666 * \brief Synchronize balance terms if this is a parallel computation
667 *
668 * \param[in] connect pointer to a cs_cdo_connect_t structure
669 * \param[in, out] b pointer to a cs_equation_balance_t to rsync
670 */
671 /*----------------------------------------------------------------------------*/
672
673 void
674 cs_equation_balance_sync(const cs_cdo_connect_t *connect,
675 cs_equation_balance_t *b);
676
677 /*----------------------------------------------------------------------------*/
678 /*!
679 * \brief Free a cs_equation_balance_t structure
680 *
681 * \param[in, out] p_balance pointer to the pointer to free
682 */
683 /*----------------------------------------------------------------------------*/
684
685 void
686 cs_equation_balance_destroy(cs_equation_balance_t **p_balance);
687
688 /*----------------------------------------------------------------------------*/
689 /*!
690 * \brief Synchronize the volumetric definitions to consider at each vertex
691 *
692 * \param[in] connect pointer to a cs_cdo_connect_t structure
693 * \param[in] n_defs number of definitions
694 * \param[in] defs number of times the values has been updated
695 * \param[in, out] def2v_idx index array to define
696 * \param[in, out] def2v_ids array of ids to define
697 */
698 /*----------------------------------------------------------------------------*/
699
700 void
701 cs_equation_sync_vol_def_at_vertices(const cs_cdo_connect_t *connect,
702 int n_defs,
703 cs_xdef_t **defs,
704 cs_lnum_t def2v_idx[],
705 cs_lnum_t def2v_ids[]);
706
707 /*----------------------------------------------------------------------------*/
708 /*!
709 * \brief Synchronize the volumetric definitions to consider at each edge
710 *
711 * \param[in] connect pointer to a cs_cdo_connect_t structure
712 * \param[in] n_defs number of definitions
713 * \param[in] defs number of times the values has been updated
714 * \param[in, out] def2e_idx index array to define
715 * \param[in, out] def2e_ids array of ids to define
716 */
717 /*----------------------------------------------------------------------------*/
718
719 void
720 cs_equation_sync_vol_def_at_edges(const cs_cdo_connect_t *connect,
721 int n_defs,
722 cs_xdef_t **defs,
723 cs_lnum_t def2e_idx[],
724 cs_lnum_t def2e_ids[]);
725
726 /*----------------------------------------------------------------------------*/
727 /*!
728 * \brief Synchronize the volumetric definitions to consider at each face
729 *
730 * \param[in] connect pointer to a cs_cdo_connect_t structure
731 * \param[in] n_defs number of definitions
732 * \param[in] defs number of times the values has been updated
733 * \param[in, out] def2f_idx index array to define
734 * \param[in, out] def2f_ids array of ids to define
735 */
736 /*----------------------------------------------------------------------------*/
737
738 void
739 cs_equation_sync_vol_def_at_faces(const cs_cdo_connect_t *connect,
740 int n_defs,
741 cs_xdef_t **defs,
742 cs_lnum_t def2f_idx[],
743 cs_lnum_t def2f_ids[]);
744
745 /*----------------------------------------------------------------------------*/
746 /*!
747 * \brief Compute the mean-value across ranks at each vertex
748 *
749 * \param[in] connect pointer to a cs_cdo_connect_t structure
750 * \param[in] dim number of entries for each vertex
751 * \param[in] counter number of occurences on this rank
752 * \param[in, out] values array to update
753 */
754 /*----------------------------------------------------------------------------*/
755
756 void
757 cs_equation_sync_vertex_mean_values(const cs_cdo_connect_t *connect,
758 int dim,
759 int *counter,
760 cs_real_t *values);
761
762 /*----------------------------------------------------------------------------*/
763
764 END_C_DECLS
765
766 #endif /* __CS_EQUATION_COMMON_H__ */
767