programmer's documentation
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
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-2017 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 #include "cs_math.h"
35 
36 /*----------------------------------------------------------------------------*/
37 
39 
40 /*============================================================================
41  * Macro definitions
42  *============================================================================*/
43 
44 /*============================================================================
45  * Type definitions
46  *============================================================================*/
47 
48 typedef enum {
49 
51  CS_QUADRATURE_BARY, /* Value at the barycenter * meas */
52  CS_QUADRATURE_BARY_SUBDIV, /* Value at the barycenter * meas on a sub-mesh */
53  CS_QUADRATURE_HIGHER, /* Unique weight but several Gauss points */
54  CS_QUADRATURE_HIGHEST, /* Specific weight for each Gauss points */
56 
58 
59 /*----------------------------------------------------------------------------*/
70 /*----------------------------------------------------------------------------*/
71 
72 typedef void
74  const cs_real_3_t v2,
75  double len,
76  cs_real_3_t gpts[],
77  double *weights);
78 
79 /*----------------------------------------------------------------------------*/
91 /*----------------------------------------------------------------------------*/
92 
93 typedef void
95  const cs_real_3_t v2,
96  const cs_real_3_t v3,
97  double area,
98  cs_real_3_t gpts[],
99  double *weights);
100 
101 /*----------------------------------------------------------------------------*/
113 /*----------------------------------------------------------------------------*/
114 
115 typedef void
117  const cs_real_3_t xe,
118  const cs_real_3_t xf,
119  const cs_real_3_t xc,
120  double vol,
121  cs_real_3_t gpts[],
122  double weights[]);
123 
124 /*============================================================================
125  * Public function prototypes
126  *============================================================================*/
127 
128 /*----------------------------------------------------------------------------*/
132 /*----------------------------------------------------------------------------*/
133 
134 void
135 cs_quadrature_setup(void);
136 
137 /*----------------------------------------------------------------------------*/
148 /*----------------------------------------------------------------------------*/
149 
150 static inline void
152  const cs_real_3_t v2,
153  double len,
154  cs_real_3_t gpts[],
155  double *w)
156 {
157  gpts[0][0] = 0.5*(v1[0] + v2[0]);
158  gpts[0][1] = 0.5*(v1[1] + v2[1]);
159  gpts[0][2] = 0.5*(v1[2] + v2[2]);
160  w[0] = len;
161 }
162 
163 /*----------------------------------------------------------------------------*/
174 /*----------------------------------------------------------------------------*/
175 
176 void
178  const cs_real_3_t v2,
179  double len,
180  cs_real_3_t gpts[],
181  double *w);
182 
183 /*----------------------------------------------------------------------------*/
194 /*----------------------------------------------------------------------------*/
195 
196 void
198  const cs_real_3_t v2,
199  double len,
200  cs_real_3_t gpts[],
201  double w[]);
202 
203 /*----------------------------------------------------------------------------*/
215 /*----------------------------------------------------------------------------*/
216 
217 static inline void
219  const cs_real_3_t v2,
220  const cs_real_3_t v3,
221  double area,
222  cs_real_3_t gpts[],
223  double *w)
224 {
225  gpts[0][0] = cs_math_onethird * (v1[0] + v2[0] + v3[0]);
226  gpts[0][1] = cs_math_onethird * (v1[1] + v2[1] + v3[1]);
227  gpts[0][2] = cs_math_onethird * (v1[2] + v2[2] + v3[2]);
228  w[0] = area;
229 }
230 
231 /*----------------------------------------------------------------------------*/
243 /*----------------------------------------------------------------------------*/
244 
245 void
247  const cs_real_3_t v2,
248  const cs_real_3_t v3,
249  double area,
250  cs_real_3_t gpts[],
251  double *w);
252 
253 /*----------------------------------------------------------------------------*/
265 /*----------------------------------------------------------------------------*/
266 
267 void
269  const cs_real_3_t v2,
270  const cs_real_3_t v3,
271  double area,
272  cs_real_3_t gpts[],
273  double w[]);
274 
275 /*----------------------------------------------------------------------------*/
287 /*----------------------------------------------------------------------------*/
288 
289 void
291  const cs_real_3_t v2,
292  const cs_real_3_t v3,
293  double area,
294  cs_real_3_t gpts[],
295  double w[]);
296 
297 /*----------------------------------------------------------------------------*/
310 /*----------------------------------------------------------------------------*/
311 
312 static inline void
314  const cs_real_3_t v2,
315  const cs_real_3_t v3,
316  const cs_real_3_t v4,
317  double vol,
318  cs_real_3_t gpts[],
319  double weight[])
320 {
321  gpts[0][0] = 0.25 * (v1[0] + v2[0] + v3[0] + v4[0]);
322  gpts[0][1] = 0.25 * (v1[1] + v2[1] + v3[1] + v4[1]);
323  gpts[0][2] = 0.25 * (v1[2] + v2[2] + v3[2] + v4[2]);
324  weight[0] = vol;
325 }
326 
327 /*----------------------------------------------------------------------------*/
340 /*----------------------------------------------------------------------------*/
341 
342 void
344  const cs_real_3_t xe,
345  const cs_real_3_t xf,
346  const cs_real_3_t xc,
347  double vol,
348  cs_real_3_t gpts[],
349  double weights[]);
350 
351 /*----------------------------------------------------------------------------*/
364 /*----------------------------------------------------------------------------*/
365 
366 void
368  const cs_real_3_t xe,
369  const cs_real_3_t xf,
370  const cs_real_3_t xc,
371  double vol,
372  cs_real_3_t gpts[],
373  double weights[]);
374 
375 /*----------------------------------------------------------------------------*/
388 /*----------------------------------------------------------------------------*/
389 
390 void
392  const cs_real_3_t xe,
393  const cs_real_3_t xf,
394  const cs_real_3_t xc,
395  double vol,
396  cs_real_3_t gpts[],
397  double weights[]);
398 
399 /*----------------------------------------------------------------------------*/
407 /*----------------------------------------------------------------------------*/
408 
409 const char *
411 
412 /*----------------------------------------------------------------------------*/
413 
415 
416 #endif /* __CS_QUADRATURE_H__ */
#define BEGIN_C_DECLS
Definition: cs_defs.h:453
Definition: cs_quadrature.h:52
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:186
static void cs_quadrature_edge_1pt(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.h:151
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:154
Definition: cs_quadrature.h:55
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
Definition: cs_quadrature.h:54
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:222
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
const char * cs_quadrature_get_type_name(const cs_quadrature_type_t type)
Return th name associated to a type of quadrature.
Definition: cs_quadrature.c:480
cs_quadrature_type_t
Definition: cs_quadrature.h:48
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:257
void cs_quadrature_setup(void)
Compute constant weights for all quadratures used.
Definition: cs_quadrature.c:109
void cs_quadrature_tet_15pts(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 5th order polynomials (order 6).
Definition: cs_quadrature.c:421
cs_real_t cs_real_3_t[3]
vector of 3 floating-point values
Definition: cs_defs.h:309
void( cs_quadrature_edge_t)(const cs_real_3_t v1, const cs_real_3_t v2, double len, cs_real_3_t gpts[], double *weights)
Generic function pointer to compute the quadrature points for an edge from v1 -> v2.
Definition: cs_quadrature.h:73
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:294
static void cs_quadrature_tria_1pt(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 (1 point) Exact for polynomial function up to order 1 (baryc...
Definition: cs_quadrature.h:218
const cs_real_t cs_math_onethird
#define END_C_DECLS
Definition: cs_defs.h:454
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 weights[])
Compute the quadrature in a tetrehedra. Exact for 2nd order polynomials (order 3).
Definition: cs_quadrature.c:337
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:377
Definition: cs_quadrature.h:53
static void cs_quadrature_tet_1pt(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 weight[])
Compute the quadrature in a tetrehedra. Exact for 1st order polynomials (order 2).
Definition: cs_quadrature.h:313
Definition: cs_quadrature.h:50
Definition: cs_quadrature.h:51