1 #ifndef __CS_HODGE_H__
2 #define __CS_HODGE_H__
3
4 /*============================================================================
5 * Manage discrete Hodge operators and closely related operators
6 *============================================================================*/
7
8 /*
9 This file is part of Code_Saturne, a general-purpose CFD tool.
10
11 Copyright (C) 1998-2021 EDF S.A.
12
13 This program is free software; you can redistribute it and/or modify it under
14 the terms of the GNU General Public License as published by the Free Software
15 Foundation; either version 2 of the License, or (at your option) any later
16 version.
17
18 This program is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
21 details.
22
23 You should have received a copy of the GNU General Public License along with
24 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
25 Street, Fifth Floor, Boston, MA 02110-1301, USA.
26 */
27
28 /*----------------------------------------------------------------------------
29 * Standard C library headers
30 *----------------------------------------------------------------------------*/
31
32 /*----------------------------------------------------------------------------
33 * Local headers
34 *----------------------------------------------------------------------------*/
35
36 #include "cs_cdo_connect.h"
37 #include "cs_cdo_local.h"
38 #include "cs_cdo_quantities.h"
39 #include "cs_param_cdo.h"
40 #include "cs_property.h"
41 #include "cs_sdm.h"
42
43 /*----------------------------------------------------------------------------*/
44
45 BEGIN_C_DECLS
46
47 /*============================================================================
48 * Macro and type definitions
49 *============================================================================*/
50
51 /*! \enum cs_hodge_type_t
52 * \brief Type of discrete Hodge operators
53 *
54 * \var CS_HODGE_TYPE_VPCD
55 * Hodge operator from primal vertices to dual cells. This operator is also
56 * equivalent to a mass matrix for vertex-based schemes
57 *
58 * \var CS_HODGE_TYPE_EPFD
59 * Hodge operator from primal edges to dual faces
60 *
61 * \var CS_HODGE_TYPE_FPED
62 * Hodge operator from primal faces to dual edges
63 *
64 * \var CS_HODGE_TYPE_EDFP
65 * Hodge operator from dual edges to primal faces
66 *
67 * \var CS_HODGE_TYPE_CPVD
68 * Hodge operator from primal cells to dual vertices
69 *
70 * \var CS_HODGE_TYPE_FB
71 * Hodge operator playing the role of mass matrix for face-based schemes.
72 * This operator deals with a hybrid space of degrees of freedom (faces and
73 * cells)
74 *
75 * \var CS_HODGE_TYPE_VC
76 * Hodge operator playing the role of mass matrix for vertex+cell-based schemes.
77 * This operator deals with a hybrid space of degrees of freedom (vertices and
78 * cells)
79 */
80
81 typedef enum {
82
83 CS_HODGE_TYPE_VPCD,
84 CS_HODGE_TYPE_EPFD,
85 CS_HODGE_TYPE_FPED,
86 CS_HODGE_TYPE_EDFP,
87 CS_HODGE_TYPE_CPVD,
88 CS_HODGE_TYPE_FB,
89 CS_HODGE_TYPE_VC,
90
91 CS_HODGE_N_TYPES
92
93 } cs_hodge_type_t;
94
95 /*! \enum cs_hodge_algo_t
96 * \brief Type of algorithm to build a discrete Hodge operator
97 *
98 * \var CS_HODGE_ALGO_VORONOI
99 * This algorithm assumes that an orthogonality condition holds between
100 * geometrical entities (e.g. the face normal/dual edge or the primal edge/dual
101 * face normal). This algorithm leads to a diagonal discrete Hodge operator
102 * (and is related to a two-point flux approximation). Be aware that using this
103 * technique on a non-orthogonal mesh leads to a consistency error which does
104 * not decrease when refining the mesh.
105 *
106 * \var CS_HODGE_ALGO_WBS
107 * WBS means "Whitney Barycentric Subdivision".
108 * This algorithm relies on a subdivision into tetrahedra of each polyhedral
109 * cell and then applies algorithms on this subdivision shape functions close
110 * to what exists in the Finite Element litterature
111 *
112 * \var CS_HODGE_ALGO_COST
113 * Orthogonal decomposition between the consistency (CO) and the stabilization
114 * (ST) parts. Several discretizations share this way to build discrete Hodge
115 * operators like SUSHI, DGA and GCR (these discretizations only differ from
116 * the scaling coefficient in front of the stabilization term)
117 *
118 * \var CS_HODGE_ALGO_OCS2
119 * This algorithm is close to \ref CS_HODGE_ALGO_COST but relies on a
120 * subdivision of each polyhedral cell corresponding to the refinement
121 * considered in \ref CS_HODGE_ALGO_COST for the building of the stabilization
122 * term
123 *
124 * \var CS_HODGE_ALGO_BUBBLE
125 * This algorithm also relies on an orthogonal decomposition between the
126 * consistency (CO) and the stabilization (ST) parts but the stabilization part
127 * relies on a bubble function (similar to what can be found in the Finite
128 * Element literature)
129 *
130 * \var CS_HODGE_ALGO_AUTO
131 * Automatic switch between the above-mentioned algorithms (not used for the
132 * moment).
133 */
134
135 typedef enum {
136
137 CS_HODGE_ALGO_VORONOI,
138 CS_HODGE_ALGO_WBS,
139 CS_HODGE_ALGO_COST,
140 CS_HODGE_ALGO_OCS2,
141 CS_HODGE_ALGO_BUBBLE,
142 CS_HODGE_ALGO_AUTO,
143
144 CS_HODGE_N_ALGOS
145
146 } cs_hodge_algo_t;
147
148 /*!
149 * \struct cs_hodge_param_t
150 * \brief Structure storing all metadata/parameters related to the usage of a
151 * discrete Hodge operator
152 */
153
154 typedef struct {
155
156 bool inv_pty;
157
158 /*!< \brief Inversion of the property evaluation or not
159 *
160 * Each Hodge operator is associated to a property. Since either the value or
161 * its reciprocal is considered, this parameter helps one to know in which
162 * situation we are. If this is set to true, then one needs to consider the
163 * reciprocal of the associated property */
164
165 cs_hodge_type_t type; /*!< Type of discrete Hodge operator */
166 cs_hodge_algo_t algo; /*!< Type of algorithm used to build the operator */
167
168 double coef;
169
170 /*!< \brief Scaling coefficient value
171 *
172 * Value of the coefficient scaling the stabilization part if the COST or
173 * OCS2 algorithm is used. Otherwise the value is set to 0 and ignored.
174 */
175
176 } cs_hodge_param_t;
177
178 /* DISCRETE HODGE OPERATORS */
179 /* ======================== */
180
181 /*!
182 * \struct cs_hodge_t
183 * \brief Structure associated to a discrete Hodge operator *
184 */
185
186 typedef struct {
187
188 const cs_hodge_param_t *param; /*!< Set of parameters (shared) */
189
190 cs_property_data_t *pty_data;
191
192 /*!< \brief Data evaluation
193 *
194 * Each Hodge operator is associated to a property. Pointer to a structure
195 * storing the evaluation of the associated property
196 */
197
198 cs_sdm_t *matrix; /*!< Matrix storing operator values */
199
200 } cs_hodge_t;
201
202 /*----------------------------------------------------------------------------*/
203 /*!
204 * \brief Build a discrete Hodge operator or a related operator (such as the
205 * stiffmess matrix) for a given cell.
206 * The discrete Hodge operator is stored in hodge->matrix whereas the
207 * associated operator is stored in cb->loc
208 *
209 * \param[in] cm pointer to a cs_cell_mesh_t struct.
210 * \param[in, out] hodge pointer to a cs_hodge_t structure
211 * \param[in, out] cb pointer to a cs_cell_builder_t structure
212 */
213 /*----------------------------------------------------------------------------*/
214
215 typedef void
216 (cs_hodge_compute_t)(const cs_cell_mesh_t *cm,
217 cs_hodge_t *hodge,
218 cs_cell_builder_t *cb);
219
220 /*============================================================================
221 * Static inline public function definitions
222 *============================================================================*/
223
224 /*----------------------------------------------------------------------------*/
225 /*!
226 * \brief Check if two sets of parameters related to how build a discrete
227 * Hodge operator are similar.
228 *
229 * \param[in] h1_info pointer to a first cs_hodge_param_t structure
230 * \param[in] h2_info pointer to a second cs_hodge_param_t structure
231 *
232 * \return true or false
233 */
234 /*----------------------------------------------------------------------------*/
235
236 static inline bool
cs_hodge_param_is_similar(const cs_hodge_param_t h1_info,const cs_hodge_param_t h2_info)237 cs_hodge_param_is_similar(const cs_hodge_param_t h1_info,
238 const cs_hodge_param_t h2_info)
239 {
240 if (h1_info.type != h2_info.type)
241 return false;
242 if (h1_info.algo != h2_info.algo)
243 return false;
244 if (h1_info.algo == CS_HODGE_ALGO_COST ||
245 h1_info.algo == CS_HODGE_ALGO_BUBBLE) {
246 if (fabs(h1_info.coef - h2_info.coef) > 0)
247 return false;
248 else
249 return true;
250 }
251 else
252 return true;
253 }
254
255 /*============================================================================
256 * Public function definitions
257 *============================================================================*/
258
259 /*----------------------------------------------------------------------------*/
260 /*!
261 * \brief Create and initialize a pointer to a cs_hodge_t structure
262 *
263 * \param[in] connect pointer to cs_cdo_connect_t structure
264 * \param[in] property pointer to a property structure
265 * \param[in] hp pointer to a cs_hodge_param_t structure
266 * \param[in] need_tensor true if one needs a tensor otherwise false
267 * \param[in] need_eigen true if one needs to compute eigen valuese
268 *
269 * \return a pointer to the new allocated cs_hodge_t structure
270 */
271 /*----------------------------------------------------------------------------*/
272
273 cs_hodge_t *
274 cs_hodge_create(const cs_cdo_connect_t *connect,
275 const cs_property_t *property,
276 const cs_hodge_param_t *hp,
277 bool need_tensor,
278 bool need_eigen);
279
280 /*----------------------------------------------------------------------------*/
281 /*!
282 * \brief Create and initialize an array of pointers to a cs_hodge_t
283 * structures. This array is of size the number of OpenMP threads.
284 * Only the one associated to the current thread is set.
285 *
286 * \param[in] connect pointer to cs_cdo_connect_t structure
287 * \param[in] property pointer to a property structure
288 * \param[in] hp pointer to a cs_hodge_param_t structure
289 * \param[in] need_tensor true if one needs a tensor otherwise false
290 * \param[in] need_eigen true if one needs to compute eigen valuese
291 *
292 * \return an array of pointers of cs_hodge_t structures
293 */
294 /*----------------------------------------------------------------------------*/
295
296 cs_hodge_t **
297 cs_hodge_init_context(const cs_cdo_connect_t *connect,
298 const cs_property_t *property,
299 const cs_hodge_param_t *hp,
300 bool need_tensor,
301 bool need_eigen);
302
303 /*----------------------------------------------------------------------------*/
304 /*!
305 * \brief Free a cs_hodge_t structure
306 *
307 * \param[in, out] p_hodge double pointer to a cs_hodge_t structure
308 */
309 /*----------------------------------------------------------------------------*/
310
311 void
312 cs_hodge_free(cs_hodge_t **p_hodge);
313
314 /*----------------------------------------------------------------------------*/
315 /*!
316 * \brief Free a set of cs_hodge_t structures
317 *
318 * \param[in, out] p_hodges triple pointer to cs_hodge_t structures
319 */
320 /*----------------------------------------------------------------------------*/
321
322 void
323 cs_hodge_free_context(cs_hodge_t ***p_hodges);
324
325 /*----------------------------------------------------------------------------*/
326 /*!
327 * \brief Retrieve a function pointer to compute a discrete Hodge operator
328 *
329 * \param[in] calling_func name of the calling function
330 * \param[in] hp a cs_hodge_param_t structure
331 *
332 * \return a pointer to the corresponding function
333 */
334 /*----------------------------------------------------------------------------*/
335
336 cs_hodge_compute_t *
337 cs_hodge_get_func(const char *calling_func,
338 const cs_hodge_param_t hp);
339
340 /*----------------------------------------------------------------------------*/
341 /*!
342 * \brief Check the consistency of the settings between terms related to a
343 * mass matrix and define the common algorithm to use.
344 * If a term should not be considered, set the algorithm to
345 * CS_HODGE_N_ALGOS
346 *
347 * \param[in] eqname name of the equation to check
348 * \param[in] reac_algo optional algo. used for the reaction term
349 * \param[in] time_algo optional algo. used for the unsteady term
350 * \param[in] srct_algo optional algo. used for the source term
351 *
352 * \return the common algorithm to use
353 */
354 /*----------------------------------------------------------------------------*/
355
356 cs_hodge_algo_t
357 cs_hodge_set_mass_algo(const char *eqname,
358 cs_hodge_algo_t reac_algo,
359 cs_hodge_algo_t time_algo,
360 cs_hodge_algo_t srct_algo);
361
362 /*----------------------------------------------------------------------------*/
363 /*!
364 * \brief Output the settings related to a cs_hodge_param_t structure
365 *
366 * \param[in] prefix optional string
367 * \param[in] property optional pointer to a property structure
368 * \param[in] hp a cs_hodge_param_t structure
369 */
370 /*----------------------------------------------------------------------------*/
371
372 void
373 cs_hodge_param_log(const char *prefix,
374 const cs_property_t *property,
375 const cs_hodge_param_t hp);
376
377 /*----------------------------------------------------------------------------*/
378 /*!
379 * \brief Copy the set of parameters associated to a discrete Hodge operator
380 * to another one
381 *
382 * \param[in] h_ref reference set of parameters
383 * \param[in, out] h_cpy set of parameters to update
384 */
385 /*----------------------------------------------------------------------------*/
386
387 void
388 cs_hodge_copy_parameters(const cs_hodge_param_t *h_ref,
389 cs_hodge_param_t *h_cpy);
390
391 /*----------------------------------------------------------------------------*/
392 /*!
393 * \brief Set the property value (scalar- or tensor-valued) related to a
394 * discrete Hodge operator inside a cell and if needed other related
395 * quantities
396 *
397 * \param[in] c_id id of the cell to deal with
398 * \param[in] t_eval time at which one performs the evaluation
399 * \param[in] c_flag flag related to this cell
400 * \param[in, out] hodge pointer to a cs_hodge_t structure
401 */
402 /*----------------------------------------------------------------------------*/
403
404 void
405 cs_hodge_set_property_value(const cs_lnum_t c_id,
406 const cs_real_t t_eval,
407 const cs_flag_t c_flag,
408 cs_hodge_t *hodge);
409
410 /*----------------------------------------------------------------------------*/
411 /*!
412 * \brief Set the property value (scalar- or tensor-valued) related to a
413 * discrete Hodge operator inside a cell and if needed ohter related
414 * quantities.
415 * Cell-wise variant (usage of cs_cell_mesh_t structure)
416 *
417 * \param[in] cm pointer to a cs_cell_mesh_t structure
418 * \param[in] t_eval time at which one performs the evaluation
419 * \param[in] c_flag flag related to this cell
420 * \param[in, out] hodge pointer to a cs_hodge_t structure
421 */
422 /*----------------------------------------------------------------------------*/
423
424 void
425 cs_hodge_set_property_value_cw(const cs_cell_mesh_t *cm,
426 const cs_real_t t_eval,
427 const cs_flag_t c_flag,
428 cs_hodge_t *hodge);
429
430 /*----------------------------------------------------------------------------*/
431 /*!
432 * \brief Build a local stiffness matrix using the generic COST algo.
433 * The computed matrix is stored in cb->loc and the related discrete
434 * hodge operator in hodge->matrix
435 * Case of CDO face-based schemes
436 *
437 * \param[in] cm pointer to a cs_cell_mesh_t structure
438 * \param[in, out] hodge pointer to a cs_hodge_t structure
439 * \param[in, out] cb pointer to a cs_cell_builder_t structure
440 */
441 /*----------------------------------------------------------------------------*/
442
443 void
444 cs_hodge_fb_cost_get_stiffness(const cs_cell_mesh_t *cm,
445 cs_hodge_t *hodge,
446 cs_cell_builder_t *cb);
447
448 /*----------------------------------------------------------------------------*/
449 /*!
450 * \brief Build a local stiffness matrix using the generic COST algo.
451 * with the usage of bubble stabilization.
452 * The computed matrix is stored in cb->loc and the related discrete
453 * hodge operator in hodge->matrix
454 * Case of CDO face-based schemes
455 *
456 * \param[in] cm pointer to a cs_cell_mesh_t structure
457 * \param[in, out] hodge pointer to a cs_hodge_t structure
458 * \param[in, out] cb pointer to a cs_cell_builder_t structure
459 */
460 /*----------------------------------------------------------------------------*/
461
462 void
463 cs_hodge_fb_bubble_get_stiffness(const cs_cell_mesh_t *cm,
464 cs_hodge_t *hodge,
465 cs_cell_builder_t *cb);
466
467 /*----------------------------------------------------------------------------*/
468 /*!
469 * \brief Build a local stiffness matrix using the Voronoi algorithm
470 * The computed matrix is stored in cb->loc and the related discrete
471 * hodge operator in hodge->matrix
472 * Case of CDO face-based schemes
473 *
474 * \param[in] cm pointer to a cs_cell_mesh_t structure
475 * \param[in, out] hodge pointer to a cs_hodge_t structure
476 * \param[in, out] cb pointer to a cs_cell_builder_t structure
477 */
478 /*----------------------------------------------------------------------------*/
479
480 void
481 cs_hodge_fb_voro_get_stiffness(const cs_cell_mesh_t *cm,
482 cs_hodge_t *hodge,
483 cs_cell_builder_t *cb);
484
485 /*----------------------------------------------------------------------------*/
486 /*!
487 * \brief Build a local stiffness matrix using the generic COST algo.
488 * The computed matrix is stored in cb->loc and the related discrete
489 * hodge operator in hodge->matrix
490 * Case of CDO vertex-based schemes and an isotropic property
491 *
492 * \param[in] cm pointer to a cs_cell_mesh_t structure
493 * \param[in, out] hodge pointer to a cs_hodge_t structure
494 * \param[in, out] cb pointer to a cs_cell_builder_t structure
495 */
496 /*----------------------------------------------------------------------------*/
497
498 void
499 cs_hodge_vb_cost_get_iso_stiffness(const cs_cell_mesh_t *cm,
500 cs_hodge_t *hodge,
501 cs_cell_builder_t *cb);
502
503 /*----------------------------------------------------------------------------*/
504 /*!
505 * \brief Build a local stiffness matrix using the generic COST algo.
506 * The computed matrix is stored in cb->loc and the related discrete
507 * hodge operator in hodge->matrix
508 * Case of CDO vertex-based schemes and an anistropic property
509 *
510 * \param[in] cm pointer to a cs_cell_mesh_t structure
511 * \param[in, out] hodge pointer to a cs_hodge_t structure
512 * \param[in, out] cb pointer to a cs_cell_builder_t structure
513 */
514 /*----------------------------------------------------------------------------*/
515
516 void
517 cs_hodge_vb_cost_get_aniso_stiffness(const cs_cell_mesh_t *cm,
518 cs_hodge_t *hodge,
519 cs_cell_builder_t *cb);
520
521 /*----------------------------------------------------------------------------*/
522 /*!
523 * \brief Build a local stiffness matrix using the generic Bubble algo.
524 * The computed matrix is stored in cb->loc and the related discrete
525 * hodge operator in hodge->matrix
526 * Case of CDO vertex-based schemes and isotropic material property
527 *
528 * \param[in] cm pointer to a cs_cell_mesh_t structure
529 * \param[in, out] hodge pointer to a cs_hodge_t structure
530 * \param[in, out] cb pointer to a cs_cell_builder_t structure
531 */
532 /*----------------------------------------------------------------------------*/
533
534 void
535 cs_hodge_vb_bubble_get_iso_stiffness(const cs_cell_mesh_t *cm,
536 cs_hodge_t *hodge,
537 cs_cell_builder_t *cb);
538
539 /*----------------------------------------------------------------------------*/
540 /*!
541 * \brief Build a local stiffness matrix using the generic Bubble algo.
542 * The computed matrix is stored in cb->loc and the related discrete
543 * hodge operator in hodge->matrix
544 * Case of CDO vertex-based schemes and anisotropic material property
545 *
546 * \param[in] cm pointer to a cs_cell_mesh_t structure
547 * \param[in, out] hodge pointer to a cs_hodge_t structure
548 * \param[in, out] cb pointer to a cs_cell_builder_t structure
549 */
550 /*----------------------------------------------------------------------------*/
551
552 void
553 cs_hodge_vb_bubble_get_aniso_stiffness(const cs_cell_mesh_t *cm,
554 cs_hodge_t *hodge,
555 cs_cell_builder_t *cb);
556
557 /*----------------------------------------------------------------------------*/
558 /*!
559 * \brief Build a local stiffness matrix using the Orthogonal
560 * Consistent/Sub-Stabilization decomposition (OCS2) with a
561 * subdivision of pvol_{e,c}.
562 * The computed matrix is stored in cb->loc and the related discrete
563 * hodge operator in hodge->matrix
564 * Case Vb schemes and an anisotropic material property
565 *
566 * \param[in] cm pointer to a cs_cell_mesh_t structure
567 * \param[in, out] hodge pointer to a cs_hodge_t structure
568 * \param[in, out] cb pointer to a cs_cell_builder_t structure
569 */
570 /*----------------------------------------------------------------------------*/
571
572 void
573 cs_hodge_vb_ocs2_get_aniso_stiffness(const cs_cell_mesh_t *cm,
574 cs_hodge_t *hodge,
575 cs_cell_builder_t *cb);
576
577 /*----------------------------------------------------------------------------*/
578 /*!
579 * \brief Build a local stiffness matrix using the generic COST algo.
580 * The computed matrix is stored in cb->loc and the related discrete
581 * hodge operator in hodge->matrix
582 * Case of CDO vertex-based schemes
583 *
584 * \param[in] cm pointer to a cs_cell_mesh_t structure
585 * \param[in, out] hodge pointer to a cs_hodge_t structure
586 * \param[in, out] cb pointer to a cs_cell_builder_t structure
587 */
588 /*----------------------------------------------------------------------------*/
589
590 void
591 cs_hodge_vb_cost_get_stiffness(const cs_cell_mesh_t *cm,
592 cs_hodge_t *hodge,
593 cs_cell_builder_t *cb);
594
595 /*----------------------------------------------------------------------------*/
596 /*!
597 * \brief Build a local stiffness matrix using the Voronoi algorithm
598 * The computed matrix is stored in cb->loc and the related discrete
599 * hodge operator in hodge->matrix
600 * Case of CDO vertex-based schemes
601 *
602 * \param[in] cm pointer to a cs_cell_mesh_t structure
603 * \param[in, out] hodge pointer to a cs_hodge_t structure
604 * \param[in, out] cb pointer to a cs_cell_builder_t structure
605 */
606 /*----------------------------------------------------------------------------*/
607
608 void
609 cs_hodge_vb_voro_get_stiffness(const cs_cell_mesh_t *cm,
610 cs_hodge_t *hodge,
611 cs_cell_builder_t *cb);
612
613 /*----------------------------------------------------------------------------*/
614 /*!
615 * \brief Build a local stiffness matrix using the generic WBS algo.
616 * WBS standing for Whitney Barycentric Subdivision (WBS)
617 * algo.
618 * The computed matrix is stored in cb->loc
619 *
620 * \param[in] cm pointer to a cs_cell_mesh_t structure
621 * \param[in, out] hodge pointer to a cs_hodge_t structure
622 * \param[in, out] cb pointer to a cs_cell_builder_t structure
623 */
624 /*----------------------------------------------------------------------------*/
625
626 void
627 cs_hodge_vb_wbs_get_stiffness(const cs_cell_mesh_t *cm,
628 cs_hodge_t *hodge,
629 cs_cell_builder_t *cb);
630
631 /*----------------------------------------------------------------------------*/
632 /*!
633 * \brief Build a local stiffness matrix using the generic WBS algo.
634 * WBS standing for Whitney Barycentric Subdivision (WBS) algo.
635 * The computed matrix is stored in cb->loc
636 *
637 * \param[in] cm pointer to a cs_cell_mesh_t structure
638 * \param[in, out] hodge pointer to a cs_hodge_t structure
639 * \param[in, out] cb pointer to a cs_cell_builder_t structure
640 */
641 /*----------------------------------------------------------------------------*/
642
643 void
644 cs_hodge_vcb_get_stiffness(const cs_cell_mesh_t *cm,
645 cs_hodge_t *hodge,
646 cs_cell_builder_t *cb);
647
648 /*----------------------------------------------------------------------------*/
649 /*!
650 * \brief Build a local Hodge operator on a given cell which is equivalent of
651 * a mass matrix. It relies on a CO+ST algo. and is specific to CDO-Fb
652 * schemes.
653 * The discrete Hodge operator is stored in hodge->matrix
654 *
655 * \param[in] cm pointer to a cs_cell_mesh_t structure
656 * \param[in, out] hodge pointer to a cs_hodge_t structure
657 * \param[in, out] cb pointer to a cs_cell_builder_t structure
658 */
659 /*----------------------------------------------------------------------------*/
660
661 void
662 cs_hodge_fb_get(const cs_cell_mesh_t *cm,
663 cs_hodge_t *hodge,
664 cs_cell_builder_t *cb);
665
666 /*----------------------------------------------------------------------------*/
667 /*!
668 * \brief Build a local Hodge operator for a given cell using the Voronoi
669 * algo. This leads to a diagonal operator.
670 * This function is specific for vertex+cell-based schemes
671 * The discrete Hodge operator is stored in hodge->matrix
672 *
673 * \param[in] cm pointer to a cs_cell_mesh_t structure
674 * \param[in, out] hodge pointer to a cs_hodge_t structure
675 * \param[in, out] cb pointer to a cs_cell_builder_t structure
676 */
677 /*----------------------------------------------------------------------------*/
678
679 void
680 cs_hodge_vcb_voro_get(const cs_cell_mesh_t *cm,
681 cs_hodge_t *hodge,
682 cs_cell_builder_t *cb);
683
684 /*----------------------------------------------------------------------------*/
685 /*!
686 * \brief Build a local Hodge operator for a given cell using the WBS algo.
687 * This function is specific for vertex+cell-based schemes
688 * The discrete Hodge operator is stored in hodge->matrix
689 *
690 * \param[in] cm pointer to a cs_cell_mesh_t structure
691 * \param[in, out] hodge pointer to a cs_hodge_t structure
692 * \param[in, out] cb pointer to a cs_cell_builder_t structure
693 */
694 /*----------------------------------------------------------------------------*/
695
696 void
697 cs_hodge_vcb_wbs_get(const cs_cell_mesh_t *cm,
698 cs_hodge_t *hodge,
699 cs_cell_builder_t *cb);
700
701 /*----------------------------------------------------------------------------*/
702 /*!
703 * \brief Build a local Hodge operator for a given cell using WBS algo.
704 * Hodge op. from primal vertices to dual cells.
705 * This function is specific for vertex-based schemes
706 * The discrete Hodge operator is stored in hodge->matrix
707 *
708 * \param[in] cm pointer to a cs_cell_mesh_t structure
709 * \param[in, out] hodge pointer to a cs_hodge_t structure
710 * \param[in, out] cb pointer to a cs_cell_builder_t structure
711 */
712 /*----------------------------------------------------------------------------*/
713
714 void
715 cs_hodge_vpcd_wbs_get(const cs_cell_mesh_t *cm,
716 cs_hodge_t *hodge,
717 cs_cell_builder_t *cb);
718
719 /*----------------------------------------------------------------------------*/
720 /*!
721 * \brief Build a local Hodge operator for a given cell using VORONOI algo.
722 * Hodge op. from primal vertices to dual cells.
723 * This function is specific for vertex-based schemes
724 * The discrete Hodge operator is stored in hodge->matrix
725 *
726 * \param[in] cm pointer to a cs_cell_mesh_t structure
727 * \param[in, out] hodge pointer to a cs_hodge_t structure
728 * \param[in, out] cb pointer to a cs_cell_builder_t structure
729 */
730 /*----------------------------------------------------------------------------*/
731
732 void
733 cs_hodge_vpcd_voro_get(const cs_cell_mesh_t *cm,
734 cs_hodge_t *hodge,
735 cs_cell_builder_t *cb);
736
737 /*----------------------------------------------------------------------------*/
738 /*!
739 * \brief Build a local Hodge operator for a given cell using VORONOI algo.
740 * Hodge op. from primal edges to dual faces.
741 * The discrete Hodge operator is stored in hodge->matrix
742 * This function is specific for vertex-based schemes
743 *
744 * \param[in] cm pointer to a cs_cell_mesh_t structure
745 * \param[in, out] hodge pointer to a cs_hodge_t structure
746 * \param[in, out] cb pointer to a cs_cell_builder_t structure
747 */
748 /*----------------------------------------------------------------------------*/
749
750 void
751 cs_hodge_epfd_voro_get(const cs_cell_mesh_t *cm,
752 cs_hodge_t *hodge,
753 cs_cell_builder_t *cb);
754
755 /*----------------------------------------------------------------------------*/
756 /*!
757 * \brief Build a local Hodge operator for a given cell using the COST algo.
758 * Hodge op. from primal edges to dual faces.
759 * The discrete Hodge operator is stored in hodge->matrix
760 * This function is specific for vertex-based schemes
761 *
762 * \param[in] cm pointer to a cs_cell_mesh_t structure
763 * \param[in, out] hodge pointer to a cs_hodge_t structure
764 * \param[in, out] cb pointer to a cs_cell_builder_t structure
765 */
766 /*----------------------------------------------------------------------------*/
767
768 void
769 cs_hodge_epfd_cost_get(const cs_cell_mesh_t *cm,
770 cs_hodge_t *hodge,
771 cs_cell_builder_t *cb);
772
773 /*----------------------------------------------------------------------------*/
774 /*!
775 * \brief Build a local Hodge operator for a given cell using the COST algo.
776 * with a bubble stabilization.
777 * The discrete Hodge operator is stored in hodge->matrix
778 * Hodge op. from primal edges to dual faces. This function is
779 * specific for vertex-based schemes
780 *
781 * \param[in] cm pointer to a cs_cell_mesh_t structure
782 * \param[in, out] hodge pointer to a cs_hodge_t structure
783 * \param[in, out] cb pointer to a cs_cell_builder_t structure
784 */
785 /*----------------------------------------------------------------------------*/
786
787 void
788 cs_hodge_epfd_bubble_get(const cs_cell_mesh_t *cm,
789 cs_hodge_t *hodge,
790 cs_cell_builder_t *cb);
791
792 /*----------------------------------------------------------------------------*/
793 /*!
794 * \brief Build a local Hodge operator for a given cell using the Orthogonal
795 * Consistent/Sub-Stabilization decomposition (OCS2) with a
796 * subdivision of pvol_{e,c}.
797 * The discrete Hodge operator is stored in hodge->matrix
798 * Hodge op. from primal edges to dual faces.
799 * This function is specific for vertex-based schemes
800 *
801 * \param[in] cm pointer to a cs_cell_mesh_t structure
802 * \param[in, out] hodge pointer to a cs_hodge_t structure
803 * \param[in, out] cb pointer to a cs_cell_builder_t structure
804 */
805 /*----------------------------------------------------------------------------*/
806
807 void
808 cs_hodge_epfd_ocs2_get(const cs_cell_mesh_t *cm,
809 cs_hodge_t *hodge,
810 cs_cell_builder_t *cb);
811
812 /*----------------------------------------------------------------------------*/
813 /*!
814 * \brief Build a local Hodge operator for a given cell using VORONOI algo.
815 * Hodge op. from primal faces to dual edges.
816 * The discrete Hodge operator is stored in hodge->matrix
817 * This function is related to cell-based schemes
818 *
819 * \param[in] cm pointer to a cs_cell_mesh_t structure
820 * \param[in, out] hodge pointer to a cs_hodge_t structure
821 * \param[in, out] cb pointer to a cs_cell_builder_t structure
822 */
823 /*----------------------------------------------------------------------------*/
824
825 void
826 cs_hodge_fped_voro_get(const cs_cell_mesh_t *cm,
827 cs_hodge_t *hodge,
828 cs_cell_builder_t *cb);
829
830 /*----------------------------------------------------------------------------*/
831 /*!
832 * \brief Build a local Hodge operator for a given cell using the COST algo.
833 * Hodge op. from primal faces to dual edges.
834 * The discrete Hodge operator is stored in hodge->matrix
835 * This function is related to cell-based schemes
836 *
837 * \param[in] cm pointer to a cs_cell_mesh_t structure
838 * \param[in, out] hodge pointer to a cs_hodge_t structure
839 * \param[in, out] cb pointer to a cs_cell_builder_t structure
840 */
841 /*----------------------------------------------------------------------------*/
842
843 void
844 cs_hodge_fped_cost_get(const cs_cell_mesh_t *cm,
845 cs_hodge_t *hodge,
846 cs_cell_builder_t *cb);
847
848 /*----------------------------------------------------------------------------*/
849 /*!
850 * \brief Build a local Hodge operator for a given cell using the Bubble algo.
851 * Hodge op. from primal faces to dual edges.
852 * The discrete Hodge operator is stored in hodge->matrix
853 * This function is related to cell-based schemes
854 *
855 * \param[in] cm pointer to a cs_cell_mesh_t structure
856 * \param[in, out] hodge pointer to a cs_hodge_t structure
857 * \param[in, out] cb pointer to a cs_cell_builder_t structure
858 */
859 /*----------------------------------------------------------------------------*/
860
861 void
862 cs_hodge_fped_bubble_get(const cs_cell_mesh_t *cm,
863 cs_hodge_t *hodge,
864 cs_cell_builder_t *cb);
865
866 /*----------------------------------------------------------------------------*/
867 /*!
868 * \brief Build a local Hodge operator for a given cell using VORONOI algo.
869 * Hodge op. from dual edges to primal faces.
870 * The discrete Hodge operator is stored in hodge->matrix
871 * This function is related to face-based schemes
872 *
873 * \param[in] cm pointer to a cs_cell_mesh_t structure
874 * \param[in, out] hodge pointer to a cs_hodge_t structure
875 * \param[in, out] cb pointer to a cs_cell_builder_t structure
876 */
877 /*----------------------------------------------------------------------------*/
878
879 void
880 cs_hodge_edfp_voro_get(const cs_cell_mesh_t *cm,
881 cs_hodge_t *hodge,
882 cs_cell_builder_t *cb);
883
884 /*----------------------------------------------------------------------------*/
885 /*!
886 * \brief Build a local Hodge operator for a given cell using the COST algo.
887 * Hodge op. from dual edges to primal faces.
888 * This function is related to face-based schemes
889 *
890 * \param[in] cm pointer to a cs_cell_mesh_t struct.
891 * \param[in, out] hodge pointer to a cs_hodge_t structure
892 * \param[in, out] cb pointer to a cs_cell_builder_t structure
893 */
894 /*----------------------------------------------------------------------------*/
895
896 void
897 cs_hodge_edfp_cost_get(const cs_cell_mesh_t *cm,
898 cs_hodge_t *hodge,
899 cs_cell_builder_t *cb);
900
901 /*----------------------------------------------------------------------------*/
902 /*!
903 * \brief Build a local Hodge operator for a given cell using the Bubble algo.
904 * Hodge op. from dual edges to primal faces.
905 * The discrete Hodge operator is stored in hodge->matrix
906 * This function is related to face-based schemes
907 *
908 * \param[in] cm pointer to a cs_cell_mesh_t structure
909 * \param[in, out] hodge pointer to a cs_hodge_t structure
910 * \param[in, out] cb pointer to a cs_cell_builder_t structure
911 */
912 /*----------------------------------------------------------------------------*/
913
914 void
915 cs_hodge_edfp_bubble_get(const cs_cell_mesh_t *cm,
916 cs_hodge_t *hodge,
917 cs_cell_builder_t *cb);
918
919 /*----------------------------------------------------------------------------*/
920 /*!
921 * \brief Build a local Hodge operator for a given cell using the COST algo.
922 * Hodge op. from dual edges to primal faces.
923 * This function is related to face-based schemes
924 *
925 * \param[in] cm pointer to a cs_cell_mesh_t struct.
926 * \param[in, out] hodge pointer to a cs_hodge_t structure
927 * \param[in, out] cb pointer to a cs_cell_builder_t structure
928 */
929 /*----------------------------------------------------------------------------*/
930
931 void
932 cs_hodge_edfp_cost_get_opt(const cs_cell_mesh_t *cm,
933 cs_hodge_t *hodge,
934 cs_cell_builder_t *cb);
935
936 /*----------------------------------------------------------------------------*/
937 /*!
938 * \brief Build a local Hodge operator for a given cell using the Bubble algo.
939 * Hodge op. from dual edges to primal faces.
940 * The discrete Hodge operator is stored in hodge->matrix
941 * This function is related to face-based schemes
942 *
943 * \param[in] cm pointer to a cs_cell_mesh_t structure
944 * \param[in, out] hodge pointer to a cs_hodge_t structure
945 * \param[in, out] cb pointer to a cs_cell_builder_t structure
946 */
947 /*----------------------------------------------------------------------------*/
948
949 void
950 cs_hodge_edfp_bubble_get(const cs_cell_mesh_t *cm,
951 cs_hodge_t *hodge,
952 cs_cell_builder_t *cb);
953
954 /*----------------------------------------------------------------------------*/
955 /*!
956 * \brief Compute cellwise a discrete hodge operator and multiple it with
957 * a vector
958 *
959 * \param[in] connect pointer to a cs_cdo_connect_t structure
960 * \param[in] quant pointer to a cs_cdo_quantities_t structure
961 * \param[in] hodgep cs_hodge_param_t structure
962 * \param[in] pty pointer to a cs_property_t structure or NULL
963 * \param[in] in_vals vector to multiply with the discrete Hodge op.
964 * \param[in] t_eval time at which one performs the evaluation
965 * \param[in, out] result array storing the resulting matrix-vector product
966 */
967 /*----------------------------------------------------------------------------*/
968
969 void
970 cs_hodge_matvec(const cs_cdo_connect_t *connect,
971 const cs_cdo_quantities_t *quant,
972 const cs_hodge_param_t hodgep,
973 const cs_property_t *pty,
974 const cs_real_t in_vals[],
975 cs_real_t t_eval,
976 cs_real_t result[]);
977
978 /*----------------------------------------------------------------------------*/
979 /*!
980 * \brief Compute cellwise a discrete hodge operator in order to define
981 * a circulation array from a flux array
982 *
983 * \param[in] connect pointer to a cs_cdo_connect_t structure
984 * \param[in] quant pointer to a cs_cdo_quantities_t structure
985 * \param[in] t_eval time at which one performs the evaluation
986 * \param[in] hodgep cs_hodge_param_t structure
987 * \param[in] pty pointer to a cs_property_t structure or NULL
988 * \param[in] flux vector to multiply with the discrete Hodge op.
989 * \param[in, out] circul array storing the resulting matrix-vector product
990 */
991 /*----------------------------------------------------------------------------*/
992
993 void
994 cs_hodge_circulation_from_flux(const cs_cdo_connect_t *connect,
995 const cs_cdo_quantities_t *quant,
996 cs_real_t t_eval,
997 const cs_hodge_param_t hodgep,
998 const cs_property_t *pty,
999 const cs_real_t flux[],
1000 cs_real_t circul[]);
1001
1002 /*----------------------------------------------------------------------------*/
1003 /*!
1004 * \brief Compute the hodge operator related to a face (i.e. a mass matrix
1005 * with unity property) using a Whitney Barycentric Subdivision (WBS)
1006 * algorithm
1007 *
1008 * \param[in] fm pointer to a cs_face_mesh_t structure
1009 * \param[in, out] hf pointer to a cs_sdm_t structure to define
1010 */
1011 /*----------------------------------------------------------------------------*/
1012
1013 void
1014 cs_hodge_compute_wbs_surfacic(const cs_face_mesh_t *fm,
1015 cs_sdm_t *hf);
1016
1017 /*----------------------------------------------------------------------------*/
1018
1019 END_C_DECLS
1020
1021 #endif /* __CS_HODGE_H__ */
1022