programmer's documentation
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
Typedefs | Enumerations | Functions
cs_quadrature.h File Reference
#include "cs_base.h"
#include "cs_defs.h"
#include "cs_math.h"
Include dependency graph for cs_quadrature.h:

Go to the source code of this file.

Typedefs

typedef 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. More...
 
typedef 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. More...
 
typedef 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. More...
 

Enumerations

enum  cs_quadrature_type_t {
  CS_QUADRATURE_NONE, CS_QUADRATURE_BARY, CS_QUADRATURE_BARY_SUBDIV, CS_QUADRATURE_HIGHER,
  CS_QUADRATURE_HIGHEST, CS_QUADRATURE_N_TYPES
}
 

Functions

void cs_quadrature_setup (void)
 Compute constant weights for all quadratures used. More...
 
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 order 1. More...
 
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 order 3. More...
 
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 order 5. More...
 
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 (barycentric approx.) More...
 
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. More...
 
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. More...
 
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. More...
 
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). More...
 
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). More...
 
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). More...
 
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). More...
 
const char * cs_quadrature_get_type_name (const cs_quadrature_type_t type)
 Return th name associated to a type of quadrature. More...
 

Typedef Documentation

typedef 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.

Parameters
[in]v1first vertex
[in]v2second vertex
[in]lenlength of edge [v1, v2]
[in,out]gptsgauss points
[in,out]wweight (same weight for the two points)
typedef 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.

Parameters
[in]xvfirst vertex
[in]xesecond vertex
[in]xfthird vertex
[in]xcfourth vertex
[in]volvolume of tetrahedron {xv, xe, xf, xc}
[in,out]gptsGauss points
[in,out]weightsweigths related to each Gauss point
typedef 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.

Parameters
[in]v1first vertex
[in]v2second vertex
[in]v3third vertex
[in]areaarea of triangle {v1, v2, v3}
[in,out]gptsGauss points
[in,out]weightsweights values

Enumeration Type Documentation

Enumerator
CS_QUADRATURE_NONE 
CS_QUADRATURE_BARY 
CS_QUADRATURE_BARY_SUBDIV 
CS_QUADRATURE_HIGHER 
CS_QUADRATURE_HIGHEST 
CS_QUADRATURE_N_TYPES 

Function Documentation

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 
)
inlinestatic

Compute quadrature points for an edge from v1 -> v2 (2 points) Exact for polynomial function up to order 1.

Parameters
[in]v1first vertex
[in]v2second vertex
[in]lenlength of edge [v1, v2]
[in,out]gptsgauss points
[in,out]wweight (same weight for the two points)
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 order 3.

Parameters
[in]v1first vertex
[in]v2second vertex
[in]lenlength of edge [v1, v2]
[in,out]gptsgauss points
[in,out]wweight (same weight for the two points)
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 order 5.

Parameters
[in]v1first vertex
[in]v2second vertex
[in]lenlength of edge [v1, v2]
[in,out]gptsgauss points
[in,out]wweights
[in]v1first vertex
[in]v2second vertex
[in]lenlength of edge [v1, v2]
[in,out]gptsgauss points
[in,out]wweights
const char* cs_quadrature_get_type_name ( const cs_quadrature_type_t  type)

Return th name associated to a type of quadrature.

Parameters
[in]typecs_quadrature_type_t
Returns
the name associated to a given type of quadrature
void cs_quadrature_setup ( void  )

Compute constant weights for all quadratures used.

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).

Parameters
[in]xvfirst vertex
[in]xesecond vertex
[in]xfthird vertex
[in]xcfourth vertex
[in]volvolume of tetrahedron {xv, xe, xf, xc}
[in,out]gpts15 Gauss points (size = 3*15)
[in,out]weights15 weigths related to each Gauss point
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[] 
)
inlinestatic

Compute the quadrature in a tetrehedra. Exact for 1st order polynomials (order 2).

Parameters
[in]v1first vertex
[in]v2second vertex
[in]v3third vertex
[in]v4fourth vertex
[in]volvolume of tetrahedron
[in,out]gpts1 Gauss point
[in,out]weightsweight related to this point (= volume)
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).

Parameters
[in]xvfirst vertex
[in]xesecond vertex
[in]xfthird vertex
[in]xcfourth vertex
[in]volvolume of tetrahedron {xv, xe, xf, xc}
[in,out]gpts4 Gauss points (size = 3*4)
[in,out]weightsweight (same value for all points)
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).

Parameters
[in]xvfirst vertex
[in]xesecond vertex
[in]xfthird vertex
[in]xcfourth vertex
[in]volvolume of tetrahedron {xv, xe, xf, xc}
[in,out]gpts5 Gauss points (size = 3*5)
[in,out]weights5 weigths related to each Gauss point
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 
)
inlinestatic

Compute quadrature points for a triangle (1 point) Exact for polynomial function up to order 1 (barycentric approx.)

Parameters
[in]v1first vertex
[in]v2second vertex
[in]v3third vertex
[in]areaarea of triangle {v1, v2, v3}
[in,out]gptsgauss points
[in,out]wweight
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.

Parameters
[in]v1first vertex
[in]v2second vertex
[in]v3third vertex
[in]areaarea of triangle {v1, v2, v3}
[in,out]gptsgauss points
[in,out]wweight (same weight for the three points)
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.

Parameters
[in]v1first vertex
[in]v2second vertex
[in]v3third vertex
[in]areaarea of triangle {v1, v2, v3}
[in,out]gptsgauss points
[in,out]wweights
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.

Parameters
[in]v1first vertex
[in]v2second vertex
[in]v3third vertex
[in]areaarea of triangle {v1, v2, v3}
[in,out]gptsgauss points
[in,out]wweights