programmer's documentation
cs_quadrature.h
Go to the documentation of this file.
1 #ifndef __CS_QUADRATURE_H__
2 #define __CS_QUADRATURE_H__
3 
4 /*============================================================================
5  * Routines to handle quadrature rules
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  * Local headers
30  *----------------------------------------------------------------------------*/
31 
32 #include "cs_base.h"
33 #include "cs_defs.h"
34 
35 /*----------------------------------------------------------------------------*/
36 
38 
39 /*============================================================================
40  * Macro definitions
41  *============================================================================*/
42 
43 /*============================================================================
44  * Type definitions
45  *============================================================================*/
46 
47 typedef enum {
48 
50  CS_QUADRATURE_BARY, /* Value at the barycenter * meas */
51  CS_QUADRATURE_BARY_SUBDIV, /* Value at the barycenter * meas on a sub-mesh */
52  CS_QUADRATURE_HIGHER, /* Unique weight but several Gauss points */
53  CS_QUADRATURE_HIGHEST, /* Specific weight for each Gauss points */
55 
57 
58 /*============================================================================
59  * Public function prototypes
60  *============================================================================*/
61 
62 /*----------------------------------------------------------------------------*/
66 /*----------------------------------------------------------------------------*/
67 
68 void
70 
71 /*----------------------------------------------------------------------------*/
82 /*----------------------------------------------------------------------------*/
83 
84 void
86  const cs_real_3_t v2,
87  double len,
88  cs_real_3_t gpts[],
89  double *w);
90 
91 /*----------------------------------------------------------------------------*/
102 /*----------------------------------------------------------------------------*/
103 
104 void
106  const cs_real_3_t v2,
107  double len,
108  cs_real_3_t gpts[],
109  double w[]);
110 
111 /*----------------------------------------------------------------------------*/
123 /*----------------------------------------------------------------------------*/
124 
125 void
127  const cs_real_3_t v2,
128  const cs_real_3_t v3,
129  double area,
130  cs_real_3_t gpts[],
131  double *w);
132 
133 /*----------------------------------------------------------------------------*/
145 /*----------------------------------------------------------------------------*/
146 
147 void
149  const cs_real_3_t v2,
150  const cs_real_3_t v3,
151  double area,
152  cs_real_3_t gpts[],
153  double w[]);
154 
155 /*----------------------------------------------------------------------------*/
167 /*----------------------------------------------------------------------------*/
168 
169 void
171  const cs_real_3_t v2,
172  const cs_real_3_t v3,
173  double area,
174  cs_real_3_t gpts[],
175  double w[]);
176 
177 /*----------------------------------------------------------------------------*/
190 /*----------------------------------------------------------------------------*/
191 
192 void
194  const cs_real_3_t xe,
195  const cs_real_3_t xf,
196  const cs_real_3_t xc,
197  double vol,
198  cs_real_3_t gpts[],
199  double *w);
200 
201 /*----------------------------------------------------------------------------*/
214 /*----------------------------------------------------------------------------*/
215 
216 void
218  const cs_real_3_t xe,
219  const cs_real_3_t xf,
220  const cs_real_3_t xc,
221  double vol,
222  cs_real_3_t gpts[],
223  double weights[]);
224 
225 /*----------------------------------------------------------------------------*/
233 /*----------------------------------------------------------------------------*/
234 
235 const char *
237 
238 /*----------------------------------------------------------------------------*/
239 
241 
242 #endif /* __CS_QUADRATURE_H__ */
Definition: cs_quadrature.h:53
Definition: cs_quadrature.h:51
#define BEGIN_C_DECLS
Definition: cs_defs.h:451
void cs_quadrature_edge_3pts(const cs_real_3_t v1, const cs_real_3_t v2, double len, cs_real_3_t gpts[], double w[])
Compute quadrature points for an edge from v1 -> v2 (3 points) Exact for polynomial function up to or...
Definition: cs_quadrature.c:166
void cs_quadrature_edge_2pts(const cs_real_3_t v1, const cs_real_3_t v2, double len, cs_real_3_t gpts[], double *w)
Compute quadrature points for an edge from v1 -> v2 (2 points) Exact for polynomial function up to or...
Definition: cs_quadrature.c:134
Definition: cs_quadrature.h:54
Definition: cs_quadrature.h:52
void cs_quadrature_tria_3pts(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 *w)
Compute quadrature points for a triangle (3 points) Exact for polynomial function up to order 2...
Definition: cs_quadrature.c:202
Definition: cs_quadrature.h:50
void cs_quadrature_tria_4pts(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 w[])
Compute quadrature points for a triangle (4 points) Exact for polynomial function up to order 3...
Definition: cs_quadrature.c:237
void cs_quadrature_setup(void)
Compute constant weights for all quadratures used.
Definition: cs_quadrature.c:99
cs_real_t cs_real_3_t[3]
vector of 3 floating-point values
Definition: cs_defs.h:309
void cs_quadrature_tria_7pts(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 w[])
Compute quadrature points for a triangle (7 points) Exact for polynomial function up to order 5...
Definition: cs_quadrature.c:274
cs_quadra_type_t
Definition: cs_quadrature.h:47
Definition: cs_quadrature.h:49
#define END_C_DECLS
Definition: cs_defs.h:452
const char * cs_quadrature_get_type_name(const cs_quadra_type_t type)
Return th name associated to a type of quadrature.
Definition: cs_quadrature.c:397
void cs_quadrature_tet_4pts(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 *w)
Compute the quadrature in a tetrehedra. Exact for 2nd order polynomials (order 3).
Definition: cs_quadrature.c:317
void cs_quadrature_tet_5pts(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[])
Compute the quadrature in a tetrehedra. Exact for 3rd order polynomials (order 4).
Definition: cs_quadrature.c:357