8.1
general documentation
cs_hho_builder.h
Go to the documentation of this file.
1 #ifndef __CS_HHO_BUILDER_H__
2 #define __CS_HHO_BUILDER_H__
3 
4 /*============================================================================
5  * Low-level functions and structures used to build the algebraic system with
6  * a cellwise process when Hybrid High Order schemes are set for the space
7  * discretization
8  *============================================================================*/
9 
10 /*
11  This file is part of code_saturne, a general-purpose CFD tool.
12 
13  Copyright (C) 1998-2023 EDF S.A.
14 
15  This program is free software; you can redistribute it and/or modify it under
16  the terms of the GNU General Public License as published by the Free Software
17  Foundation; either version 2 of the License, or (at your option) any later
18  version.
19 
20  This program is distributed in the hope that it will be useful, but WITHOUT
21  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
22  FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
23  details.
24 
25  You should have received a copy of the GNU General Public License along with
26  this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
27  Street, Fifth Floor, Boston, MA 02110-1301, USA.
28 */
29 
30 /*----------------------------------------------------------------------------
31  * Local headers
32  *----------------------------------------------------------------------------*/
33 
34 #include "cs_base.h"
35 #include "cs_basis_func.h"
36 #include "cs_cdo_connect.h"
37 #include "cs_property.h"
38 #include "cs_sdm.h"
39 #include "cs_xdef.h"
40 
41 /*----------------------------------------------------------------------------*/
42 
44 
45 /*============================================================================
46  * Macro definitions
47  *============================================================================*/
48 
49 /*============================================================================
50  * Type definitions
51  *============================================================================*/
52 
53 /* Cellwise builder for HHO discretization */
54 typedef struct {
55 
56  /* Current and maximal number of face basis allocated. This number is equal
57  to the max. number of faces for a cell */
58  short int n_face_basis;
59  short int n_max_face_basis;
60 
61  cs_basis_func_t **face_basis; /* P_(d-1)^k polynomial basis */
62  cs_basis_func_t *cell_basis; /* P_d^k polynomial basis */
63  cs_basis_func_t *grad_basis; /* P_d^(k+1) \ P_d^0 polynomial basis */
64 
65  cs_sdm_t *grad_reco_op; /* Gradient operator; Rectangular matrix */
66 
67  /* Temporary matrices defined by blocks */
68  cs_sdm_t *tmp; /* Temporary block matrix (fs x ts) */
69  cs_sdm_t *bf_t; /* Transposed of Bf (used in stabilization) */
70  cs_sdm_t *jstab; /* Stabilization part related to a face */
71  cs_sdm_t *hdg; /* Another temporary matrix */
72 
74 
75 /*============================================================================
76  * Public function prototypes
77  *============================================================================*/
78 
79 /*----------------------------------------------------------------------------*/
88 /*----------------------------------------------------------------------------*/
89 
91 cs_hho_builder_create(int order,
92  int n_fc);
93 
94 /*----------------------------------------------------------------------------*/
100 /*----------------------------------------------------------------------------*/
101 
102 void
104 
105 /*----------------------------------------------------------------------------*/
116 /*----------------------------------------------------------------------------*/
117 
118 void
120  cs_cell_builder_t *cb,
121  cs_hho_builder_t *hhob);
122 
123 /*----------------------------------------------------------------------------*/
131 /*----------------------------------------------------------------------------*/
132 
133 static inline void
135  cs_cell_builder_t *cb,
136  cs_hho_builder_t *hhob)
137 {
138  if (hhob == NULL)
139  return;
140 
141  /* Sanity checks */
142  assert(cm != NULL);
143 
144  /* Setup cell basis functions */
145  hhob->cell_basis->setup(hhob->cell_basis, cm, 0, cm->xc, cb);
146  hhob->n_face_basis = 0;
147 }
148 
149 /*----------------------------------------------------------------------------*/
162 /*----------------------------------------------------------------------------*/
163 
164 void
166  const cs_property_data_t *diff_pty,
167  cs_cell_builder_t *cb,
168  cs_hho_builder_t *hhob);
169 
170 /*----------------------------------------------------------------------------*/
180 /*----------------------------------------------------------------------------*/
181 
182 void
184  const cs_property_data_t *diff_pty,
185  cs_cell_builder_t *cb,
186  cs_hho_builder_t *hhob);
187 
188 /*----------------------------------------------------------------------------*/
202 /*----------------------------------------------------------------------------*/
203 
204 void
206  const cs_cell_mesh_t *cm,
207  cs_real_t t_eval,
208  cs_cell_builder_t *cb,
209  cs_hho_builder_t *hhob,
210  cs_real_t red[]);
211 
212 /*----------------------------------------------------------------------------*/
228 /*----------------------------------------------------------------------------*/
229 
230 void
232  const cs_cell_mesh_t *cm,
233  cs_real_t t_eval,
234  cs_cell_builder_t *cb,
235  cs_hho_builder_t *hhob,
236  cs_real_t red[]);
237 
238 /*----------------------------------------------------------------------------*/
251 /*----------------------------------------------------------------------------*/
252 
253 void
255  short int f,
256  const cs_cell_mesh_t *cm,
257  cs_real_t t_eval,
258  cs_cell_builder_t *cb,
259  cs_hho_builder_t *hhob,
260  cs_real_t res[]);
261 
262 /*----------------------------------------------------------------------------*/
275 /*----------------------------------------------------------------------------*/
276 
277 void
279  short int f,
280  const cs_cell_mesh_t *cm,
281  cs_real_t t_eval,
282  cs_cell_builder_t *cb,
283  cs_hho_builder_t *hhob,
284  cs_real_t res[]);
285 
286 /*----------------------------------------------------------------------------*/
287 
289 
290 #endif /* __CS_HHO_BUILDER_H__ */
#define BEGIN_C_DECLS
Definition: cs_defs.h:514
double cs_real_t
Floating-point value.
Definition: cs_defs.h:319
#define END_C_DECLS
Definition: cs_defs.h:515
cs_hho_builder_t * cs_hho_builder_create(int order, int n_fc)
Allocate a cs_hho_builder_t structure.
Definition: cs_hho_builder.c:751
void cs_hho_builder_free(cs_hho_builder_t **p_builder)
Free a cs_hho_builder_t structure.
Definition: cs_hho_builder.c:821
void cs_hho_builder_compute_dirichlet_v(const cs_xdef_t *def, short int f, const cs_cell_mesh_t *cm, cs_real_t t_eval, cs_cell_builder_t *cb, cs_hho_builder_t *hhob, cs_real_t res[])
Compute the projection of the Dirichlet boundary conditions onto the polynomial spaces on faces....
Definition: cs_hho_builder.c:1964
void cs_hho_builder_compute_grad_reco(const cs_cell_mesh_t *cm, const cs_property_data_t *diff_pty, cs_cell_builder_t *cb, cs_hho_builder_t *hhob)
Compute the gradient operator stemming from the relation stiffness * grad_op = rhs where stiffness is...
Definition: cs_hho_builder.c:912
static void cs_hho_builder_cellbasis_setup(const cs_cell_mesh_t *cm, cs_cell_builder_t *cb, cs_hho_builder_t *hhob)
Set-up the basis functions related to a cell only.
Definition: cs_hho_builder.h:134
void cs_hho_builder_compute_dirichlet(const cs_xdef_t *def, short int f, const cs_cell_mesh_t *cm, cs_real_t t_eval, cs_cell_builder_t *cb, cs_hho_builder_t *hhob, cs_real_t res[])
Compute the projection of the Dirichlet boundary conditions onto the polynomial spaces on faces.
Definition: cs_hho_builder.c:1836
void cs_hho_builder_diffusion(const cs_cell_mesh_t *cm, const cs_property_data_t *diff_pty, cs_cell_builder_t *cb, cs_hho_builder_t *hhob)
Compute the diffusion operator. The gradient reconstruction operator has to be built just before this...
Definition: cs_hho_builder.c:1153
void cs_hho_builder_reduction_from_analytic(const cs_xdef_t *def, const cs_cell_mesh_t *cm, cs_real_t t_eval, cs_cell_builder_t *cb, cs_hho_builder_t *hhob, cs_real_t red[])
Compute the reduction onto the polynomial spaces (cell and faces) of a function defined by an analyti...
Definition: cs_hho_builder.c:1461
void cs_hho_builder_reduction_from_analytic_v(const cs_xdef_t *def, const cs_cell_mesh_t *cm, cs_real_t t_eval, cs_cell_builder_t *cb, cs_hho_builder_t *hhob, cs_real_t red[])
Compute the reduction onto the polynomial spaces (cell and faces) of a function defined by an analyti...
Definition: cs_hho_builder.c:1633
void cs_hho_builder_cellwise_setup(const cs_cell_mesh_t *cm, cs_cell_builder_t *cb, cs_hho_builder_t *hhob)
Set-up the basis functions related to a cell, its gradient and to the faces of this cell....
Definition: cs_hho_builder.c:861
Definition: cs_basis_func.h:181
cs_basis_func_setup_t * setup
Definition: cs_basis_func.h:199
Set of local and temporary buffers.
Definition: cs_cdo_local.h:60
Set of local quantities and connectivities related to a mesh cell.
Definition: cs_cdo_local.h:203
cs_real_3_t xc
Definition: cs_cdo_local.h:217
Definition: cs_hho_builder.h:54
cs_sdm_t * tmp
Definition: cs_hho_builder.h:68
cs_basis_func_t ** face_basis
Definition: cs_hho_builder.h:61
cs_basis_func_t * grad_basis
Definition: cs_hho_builder.h:63
short int n_max_face_basis
Definition: cs_hho_builder.h:59
short int n_face_basis
Definition: cs_hho_builder.h:58
cs_sdm_t * hdg
Definition: cs_hho_builder.h:71
cs_basis_func_t * cell_basis
Definition: cs_hho_builder.h:62
cs_sdm_t * grad_reco_op
Definition: cs_hho_builder.h:65
cs_sdm_t * jstab
Definition: cs_hho_builder.h:70
cs_sdm_t * bf_t
Definition: cs_hho_builder.h:69
Structure storing the evaluation of a property and its related data.
Definition: cs_property.h:223
Structure storing medata for defining a quantity in a very flexible way.
Definition: cs_xdef.h:160