8.3
general 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-2024 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_local.h"
35#include "cs_quadrature.h"
36
37/*----------------------------------------------------------------------------*/
38
40
41/*============================================================================
42 * Macro definitions
43 *============================================================================*/
44
45/* Default basis function are rescale and modified according to the
46 geometrical inertia tensor of the related entity unless a flag "monomial"
47 is set */
48
49#define CS_BASIS_FUNC_MONOMIAL (1 << 0) /* 1: Use mononial basis functions */
50#define CS_BASIS_FUNC_GRADIENT (1 << 1) /* 2: Basis functions for gradient */
51
52/*============================================================================
53 * Type definitions
54 *============================================================================*/
55
56/*----------------------------------------------------------------------------*/
65/*----------------------------------------------------------------------------*/
66
67typedef void
68(cs_basis_func_eval_all_at_point_t) (const void *bf,
69 const cs_real_t coords[3],
70 cs_real_t *eval);
71
72/*----------------------------------------------------------------------------*/
83/*----------------------------------------------------------------------------*/
84
85typedef void
86(cs_basis_func_eval_at_point_t) (const void *bf,
87 const cs_real_t coords[3],
88 short int start,
89 short int end,
90 cs_real_t *eval);
91
92/*----------------------------------------------------------------------------*/
102/*----------------------------------------------------------------------------*/
103
104typedef void
105(cs_basis_func_setup_t) (void *pbf,
106 const cs_cell_mesh_t *cm,
107 const short int id,
108 const cs_real_t center[3],
110
111/*----------------------------------------------------------------------------*/
120/*----------------------------------------------------------------------------*/
121
122typedef void
124 const cs_cell_mesh_t *cm,
125 const short int id);
126
127/*----------------------------------------------------------------------------*/
134/*----------------------------------------------------------------------------*/
135
136typedef void
138
139/*----------------------------------------------------------------------------*/
151/*----------------------------------------------------------------------------*/
152
153typedef void
154(cs_basis_func_project_t) (const void *pbf,
155 const cs_real_t *array,
156 cs_real_t *dof);
157
158/*----------------------------------------------------------------------------*/
165/*----------------------------------------------------------------------------*/
166
167typedef void
168(cs_basis_func_dump_proj_t) (const void *pbf);
169
170
171/* Structure storing information to build and use a set of basis functions
172 *
173 * - When dealing with volumetric basis fucntions (eg cell- or grad-basis),
174 * the functions are ordered as follows:
175 * 1, x, y, z, xx, xy, xz, yy, yz, zz
176 * - When dealing with face-basis fucntions the functions are ordered as
177 * follows:
178 * 1, x, y, xx, xy, yy
179 */
180
181typedef struct {
182
183 cs_flag_t flag; // metadata
184 short int poly_order; // max. polynomial order of the basis
185 short int dim; // 2D or 3D associated geometrical entities
186 int size; // number of elementary basis functions
187
188 cs_real_t phi0; // constant basis function
189 cs_nvec3_t *axis; // over_diam * main axis of the entity
190 cs_real_3_t center; // center used for rescaling
191
192 /* For poly_order > 1 we store the degree related to each monone for
193 each basis function */
195 short int *deg; /* size = n_deg_elts * dim */
196
197 /* Function used for setting up the quantities used for defining
198 the set of basis functions */
200
201 /* Function used for evaluating basis functions on a set of points */
204
205 /* Projector is a mass matrix related to the basis function */
206 cs_sdm_t *projector;
209
210 /* Modified Cholesky factorization. To perform the inversion agains a set
211 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 * Static inline public function prototypes
228 *============================================================================*/
229
230/*----------------------------------------------------------------------------*/
238/*----------------------------------------------------------------------------*/
239
240static inline short int
242{
243 if (bf == NULL)
244 return -1;
245 else
246 return bf->poly_order;
247}
248
249/*============================================================================
250 * Public function prototypes
251 *============================================================================*/
252
253/*----------------------------------------------------------------------------*/
263/*----------------------------------------------------------------------------*/
264
267 short int k,
268 short int dim);
269
270/*----------------------------------------------------------------------------*/
282/*----------------------------------------------------------------------------*/
283
286
287/*----------------------------------------------------------------------------*/
295/*----------------------------------------------------------------------------*/
296
297void
299 cs_basis_func_t *rcv);
300
301/*----------------------------------------------------------------------------*/
309/*----------------------------------------------------------------------------*/
310
313
314/*----------------------------------------------------------------------------*/
321/*----------------------------------------------------------------------------*/
322
323void
325 cs_flag_t cell_flag);
326
327/*----------------------------------------------------------------------------*/
334/*----------------------------------------------------------------------------*/
335
336void
338 cs_flag_t *cell_flag);
339
340/*----------------------------------------------------------------------------*/
346/*----------------------------------------------------------------------------*/
347
348void
350
351/*----------------------------------------------------------------------------*/
361/*----------------------------------------------------------------------------*/
362
363void
364cs_basis_func_fprintf(FILE *fp,
365 const char *fname,
366 const cs_basis_func_t *pbf);
367
368/*----------------------------------------------------------------------------*/
369
371
372#endif /* __CS_BASIS_FUNC_H__ */
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.
Definition: cs_basis_func.cpp:2703
void cs_basis_func_dump(const cs_basis_func_t *pbf)
Dump a cs_basis_func_t structure.
Definition: cs_basis_func.cpp:2736
static short int cs_basis_func_get_poly_order(const cs_basis_func_t *bf)
Get the polynomial order of the given basis.
Definition: cs_basis_func.h:241
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:154
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.
Definition: cs_basis_func.cpp:2367
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....
Definition: cs_basis_func.cpp:2564
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:168
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 ...
Definition: cs_basis_func.cpp:2649
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.
Definition: cs_basis_func.cpp:2720
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:123
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:68
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:137
cs_basis_func_t * cs_basis_func_free(cs_basis_func_t *pbf)
Free a cs_basis_func_t structure.
Definition: cs_basis_func.cpp:2675
void cs_basis_func_fprintf(FILE *fp, const char *fname, const cs_basis_func_t *pbf)
Print a cs_basis_func_t structure Print into the file f if given otherwise open a new file named fnam...
Definition: cs_basis_func.cpp:2777
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:86
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:105
#define BEGIN_C_DECLS
Definition: cs_defs.h:542
double cs_real_t
Floating-point value.
Definition: cs_defs.h:342
cs_real_t cs_real_3_t[3]
vector of 3 floating-point values
Definition: cs_defs.h:359
#define END_C_DECLS
Definition: cs_defs.h:543
unsigned short int cs_flag_t
Definition: cs_defs.h:344
@ k
Definition: cs_field_pointer.h:72
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:130
void() cs_quadrature_tet_t(const cs_real_3_t v1, const cs_real_3_t v2, const cs_real_3_t v3, const cs_real_3_t v4, double vol, cs_real_3_t gpts[], double weights[])
Generic function to compute the quadrature points in a tetrehedra.
Definition: cs_quadrature.h:152
Definition: cs_basis_func.h:181
cs_nvec3_t * axis
Definition: cs_basis_func.h:189
cs_basis_func_compute_facto_t * compute_factorization
Definition: cs_basis_func.h:212
cs_basis_func_dump_proj_t * dump_projector
Definition: cs_basis_func.h:208
int size
Definition: cs_basis_func.h:186
short int * deg
Definition: cs_basis_func.h:195
int n_deg_elts
Definition: cs_basis_func.h:194
cs_basis_func_setup_t * setup
Definition: cs_basis_func.h:199
cs_basis_func_eval_all_at_point_t * eval_all_at_point
Definition: cs_basis_func.h:202
cs_basis_func_compute_proj_t * compute_projector
Definition: cs_basis_func.h:207
cs_basis_func_project_t * project
Definition: cs_basis_func.h:213
cs_real_t * facto
Definition: cs_basis_func.h:214
int n_gpts_tria
Definition: cs_basis_func.h:219
cs_flag_t flag
Definition: cs_basis_func.h:183
short int poly_order
Definition: cs_basis_func.h:184
cs_real_t phi0
Definition: cs_basis_func.h:188
short int dim
Definition: cs_basis_func.h:185
cs_basis_func_eval_at_point_t * eval_at_point
Definition: cs_basis_func.h:203
int n_gpts_tetra
Definition: cs_basis_func.h:221
cs_sdm_t * projector
Definition: cs_basis_func.h:206
cs_real_3_t center
Definition: cs_basis_func.h:190
cs_quadrature_tria_t * quadrature_tria
Definition: cs_basis_func.h:220
int facto_max_size
Definition: cs_basis_func.h:215
cs_quadrature_tet_t * quadrature_tetra
Definition: cs_basis_func.h:222
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
Definition: cs_defs.h:400