programmer's documentation
cs_basis_func.h
Go to the documentation of this file.
1 #ifndef __CS_BASIS_FUNC_H__
2 #define __CS_BASIS_FUNC_H__
3 
4 /*============================================================================
5  * Build a set of basis functions for cells and faces and cell gradients
6  *============================================================================*/
7 
8 /*
9  This file is part of Code_Saturne, a general-purpose CFD tool.
10 
11  Copyright (C) 1998-2018 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 
30 /*----------------------------------------------------------------------------
31  * Local headers
32  *----------------------------------------------------------------------------*/
33 
34 #include "cs_cdo.h"
35 #include "cs_cdo_local.h"
36 #include "cs_quadrature.h"
37 
38 /*----------------------------------------------------------------------------*/
39 
41 
42 /*============================================================================
43  * Macro definitions
44  *============================================================================*/
45 
46 /* Default basis function are rescale and modified according to the
47  geometrical inertia tensor of the related entity unless a flag "monomial"
48  is set */
49 
50 #define CS_BASIS_FUNC_MONOMIAL (1 << 0) /* 1: Use mononial basis functions */
51 #define CS_BASIS_FUNC_GRADIENT (1 << 1) /* 2: Basis functions for gradient */
52 
53 /*============================================================================
54  * Type definitions
55  *============================================================================*/
56 
57 /*----------------------------------------------------------------------------*/
66 /*----------------------------------------------------------------------------*/
67 
68 typedef void
70  const cs_real_t coords[3],
71  cs_real_t *eval);
72 
73 /*----------------------------------------------------------------------------*/
84 /*----------------------------------------------------------------------------*/
85 
86 typedef void
88  const cs_real_t coords[3],
89  short int start,
90  short int end,
91  cs_real_t *eval);
92 
93 /*----------------------------------------------------------------------------*/
103 /*----------------------------------------------------------------------------*/
104 
105 typedef void
107  const cs_cell_mesh_t *cm,
108  const short int id,
109  const cs_real_t center[3],
110  cs_cell_builder_t *cb);
111 
112 /*----------------------------------------------------------------------------*/
121 /*----------------------------------------------------------------------------*/
122 
123 typedef void
125  const cs_cell_mesh_t *cm,
126  const short int id);
127 
128 /*----------------------------------------------------------------------------*/
135 /*----------------------------------------------------------------------------*/
136 
137 typedef void
139 
140 /*----------------------------------------------------------------------------*/
152 /*----------------------------------------------------------------------------*/
153 
154 typedef void
155 (cs_basis_func_project_t) (const void *pbf,
156  const cs_real_t *array,
157  cs_real_t *dof);
158 
159 /*----------------------------------------------------------------------------*/
166 /*----------------------------------------------------------------------------*/
167 
168 typedef void
169 (cs_basis_func_dump_proj_t) (const void *pbf);
170 
171 
172 /* Structure storing information to build and use a set of basis functions
173  *
174  * - When dealing with volumetric basis fucntions (eg cell- or grad-basis),
175  * the functions are ordered as follows:
176  * 1, x, y, z, xx, xy, xz, yy, yz, zz
177  * - When dealing with face-basis fucntions the functions are ordered as
178  * follows:
179  * 1, x, y, xx, xy, yy
180  */
181 
182 typedef struct {
183 
184  cs_flag_t flag; // metadata
185  short int poly_order; // max. polynomial order of the basis
186  short int dim; // 2D or 3D associated geometrical entities
187  int size; // number of elementary basis functions
188 
189  cs_real_t phi0; // constant basis function
190  cs_nvec3_t *axis; // over_diam * main axis of the entity
191  cs_real_3_t center; // center used for rescaling
192 
193  /* For poly_order > 1 we store the degree related to each monone for
194  each basis function */
196  short int *deg; /* size = n_deg_elts * dim */
197 
198  /* Function used for setting up the quantities used for defining
199  the set of basis functions */
201 
202  /* Function used for evaluating basis functions on a set of points */
205 
206  /* Projector is a mass matrix related to the basis function */
207  cs_sdm_t *projector;
210 
211  /* Modified Cholesky factorization. To perform the inversion agains a set
212  of right-hand sides */
216 
217  /* Quadrature function to use for computing integral over o volume (tet) or
218  on a surface (tria) */
223 
225 
226 /*============================================================================
227  * Public function prototypes
228  *============================================================================*/
229 
230 /*----------------------------------------------------------------------------*/
240 /*----------------------------------------------------------------------------*/
241 
244  short int k,
245  short int dim);
246 
247 /*----------------------------------------------------------------------------*/
259 /*----------------------------------------------------------------------------*/
260 
263 
264 /*----------------------------------------------------------------------------*/
272 /*----------------------------------------------------------------------------*/
273 
274 void
276  cs_basis_func_t *rcv);
277 
278 /*----------------------------------------------------------------------------*/
286 /*----------------------------------------------------------------------------*/
287 
290 
291 /*----------------------------------------------------------------------------*/
297 /*----------------------------------------------------------------------------*/
298 
299 void
301 
302 /*----------------------------------------------------------------------------*/
309 /*----------------------------------------------------------------------------*/
310 
311 void
313  cs_flag_t cell_flag);
314 
315 /*----------------------------------------------------------------------------*/
322 /*----------------------------------------------------------------------------*/
323 
324 void
326  cs_flag_t *cell_flag);
327 
328 /*----------------------------------------------------------------------------*/
329 
331 
332 #endif /* __CS_BASIS_FUNC_H__ */
Definition: cs_basis_func.h:182
Definition: cs_field_pointer.h:70
cs_flag_t flag
Definition: cs_basis_func.h:184
int n_deg_elts
Definition: cs_basis_func.h:195
cs_basis_func_eval_at_point_t * eval_at_point
Definition: cs_basis_func.h:204
cs_real_t phi0
Definition: cs_basis_func.h:189
cs_basis_func_t * cs_basis_func_free(cs_basis_func_t *pbf)
Free a cs_basis_func_t structure.
void() cs_basis_func_eval_all_at_point_t(const void *bf, const cs_real_t coords[3], cs_real_t *eval)
Generic prototype for evaluating all basis functions at a given point.
Definition: cs_basis_func.h:69
void() cs_basis_func_setup_t(void *pbf, const cs_cell_mesh_t *cm, const short int id, const cs_real_t center[3], cs_cell_builder_t *cb)
Generic prototype for setting up a set of basis functions.
Definition: cs_basis_func.h:106
cs_basis_func_t * cs_basis_func_grad_create(const cs_basis_func_t *ref)
Allocate a cs_basis_func_t structure which is associated to an existing set of basis functions...
int n_gpts_tetra
Definition: cs_basis_func.h:221
short int poly_order
Definition: cs_basis_func.h:185
#define BEGIN_C_DECLS
Definition: cs_defs.h:453
Definition: cs_cdo_local.h:138
cs_basis_func_t * cs_basis_func_create(cs_flag_t flag, short int k, short int dim)
Allocate a cs_basis_func_t structure.
void cs_basis_func_set_hho_flag(cs_flag_t face_flag, cs_flag_t cell_flag)
Set options for basis functions when using HHO schemes.
void cs_basis_func_copy_setup(const cs_basis_func_t *ref, cs_basis_func_t *rcv)
Copy the center and the different axis from the reference basis Up to now, only cell basis functions ...
cs_basis_func_compute_facto_t * compute_factorization
Definition: cs_basis_func.h:213
double cs_real_t
Floating-point value.
Definition: cs_defs.h:297
Definition: cs_cdo.h:130
cs_nvec3_t * axis
Definition: cs_basis_func.h:190
void() cs_quadrature_tet_t(const cs_real_3_t xv, const cs_real_3_t xe, const cs_real_3_t xf, const cs_real_3_t xc, double vol, cs_real_3_t gpts[], double weights[])
Generic function to compute the quadrature points in a tetrehedra.
Definition: cs_quadrature.h:116
short int dim
Definition: cs_basis_func.h:186
cs_basis_func_setup_t * setup
Definition: cs_basis_func.h:200
cs_real_t * facto
Definition: cs_basis_func.h:214
cs_quadrature_tet_t * quadrature_tetra
Definition: cs_basis_func.h:222
cs_real_3_t center
Definition: cs_basis_func.h:191
short int * deg
Definition: cs_basis_func.h:196
Definition: cs_cdo_local.h:71
cs_real_t cs_real_3_t[3]
vector of 3 floating-point values
Definition: cs_defs.h:309
void() cs_basis_func_eval_at_point_t(const void *bf, const cs_real_t coords[3], short int start, short int end, cs_real_t *eval)
Generic prototype for evaluating a set of basis functions at a given point.
Definition: cs_basis_func.h:87
void() cs_basis_func_project_t(const void *pbf, const cs_real_t *array, cs_real_t *dof)
Generic prototype for projecting an array on the polynomial basis function. This results from the app...
Definition: cs_basis_func.h:155
void() cs_basis_func_compute_proj_t(void *pbf, const cs_cell_mesh_t *cm, const short int id)
Generic prototype for computing the projector to the space spanned by the basis functions.
Definition: cs_basis_func.h:124
cs_basis_func_project_t * project
Definition: cs_basis_func.h:215
void() cs_basis_func_compute_facto_t(void *pbf)
Generic prototype for computing the Modified Choloesky factorization of the projection matrix (mass m...
Definition: cs_basis_func.h:138
int n_gpts_tria
Definition: cs_basis_func.h:219
int size
Definition: cs_basis_func.h:187
cs_basis_func_eval_all_at_point_t * eval_all_at_point
Definition: cs_basis_func.h:203
cs_basis_func_dump_proj_t * dump_projector
Definition: cs_basis_func.h:209
#define END_C_DECLS
Definition: cs_defs.h:454
unsigned short int cs_flag_t
Definition: cs_defs.h:299
void() cs_basis_func_dump_proj_t(const void *pbf)
Generic prototype for printing the projector matrix related to the given basis function.
Definition: cs_basis_func.h:169
void cs_basis_func_dump(const cs_basis_func_t *pbf)
Dump a cs_basis_func_t structure.
void() cs_quadrature_tria_t(const cs_real_3_t v1, const cs_real_3_t v2, const cs_real_3_t v3, double area, cs_real_3_t gpts[], double *weights)
Generic functoin pointer to compute the quadrature points for a triangle.
Definition: cs_quadrature.h:94
void cs_basis_func_get_hho_flag(cs_flag_t *face_flag, cs_flag_t *cell_flag)
Get options for basis functions when using HHO schemes.
cs_sdm_t * projector
Definition: cs_basis_func.h:207
cs_basis_func_compute_proj_t * compute_projector
Definition: cs_basis_func.h:208
cs_quadrature_tria_t * quadrature_tria
Definition: cs_basis_func.h:220