H2Lib  3.0
Data Structures | Typedefs | Functions
singquad2d

This module is responsible for standard and singular quadrature scheme that will be applied by bem3d and derived modules such as laplacebem3d. More...

Data Structures

struct  _singquad2d
 This struct collects all type of quadrature formulas needed by the computation of matrix entries within BEM. More...
 

Typedefs

typedef struct _singquad2d singquad2d
 
typedef singquad2dpsingquad2d
 
typedef const singquad2dpcsingquad2d
 

Functions

psingquad2d build_singquad2d (pcsurface3d gr, uint q, uint q2)
 Creates a new singquad2d object containing all necessary quadrature rules for singular integrals arising in BEM applications in 3 dimensional space. More...
 
void del_singquad2d (psingquad2d sq)
 Destructor for singquad2d objects. More...
 
void weight_basisfunc_ll_singquad2d (real *x, real *y, real *w, uint nq)
 Weighting a quadrature rule for a double integral with both linear basis functions. More...
 
void weight_basisfunc_cl_singquad2d (real *x, real *y, real *w, uint nq)
 Weighting a quadrature rule for a double integral with a combination of piecewise constant and linear basis functions. More...
 
void weight_basisfunc_lc_singquad2d (real *x, real *y, real *w, uint nq)
 Weighting a quadrature rule for a double integral with a combination of piecewise linear and constant basis functions. More...
 
void weight_basisfunc_l_singquad2d (real *x, real *y, real *w, uint nq)
 Weighting a quadrature rule for a single integral with linear basis functions. More...
 
uint fast_select_quadrature (uint(*geo_t)[3], uint t, uint s)
 Determine the number of common vertices of a pair of triangles. More...
 
uint select_quadrature_singquad2d (pcsingquad2d sq, const uint *tv, const uint *sv, uint *tp, uint *sp, real **x, real **y, real **w, uint *n, real *base)
 This function is designed to select the correct quadrature rule for a current pair of triangles. More...
 

Detailed Description

This module is responsible for standard and singular quadrature scheme that will be applied by bem3d and derived modules such as laplacebem3d.

Typedef Documentation

typedef const singquad2d* pcsingquad2d

Pointer to a constant singquad2d object.

Pointer to a singquad2d object.

typedef struct _singquad2d singquad2d

singquad2d is just an abbreviation for the struct _singquad2d. It is necessary for the computation of singular integral arising in BEM applications in 3 dimensional space.

Function Documentation

psingquad2d build_singquad2d ( pcsurface3d  gr,
uint  q,
uint  q2 
)

Creates a new singquad2d object containing all necessary quadrature rules for singular integrals arising in BEM applications in 3 dimensional space.

Parameters
grCurrently used geometry.
qOrder of gaussian quadrature rule used for construction of single triangle quadrature rule and in case of distant domains.
q2Order of gaussian quadrature rule used for construction of singular quadrature cases.
Returns
returns a new _singquad2d" "singquad2d" object.
void del_singquad2d ( psingquad2d  sq)

Destructor for singquad2d objects.

Delete a singquad2d object.

Parameters
sqsingquad2d object to be deleted.
uint fast_select_quadrature ( uint(*)  geo_t[3],
uint  t,
uint  s 
)

Determine the number of common vertices of a pair of triangles.

Parameters
geo_tArray of all triangles in the geometry.
tIndex for the first triangle.
sIndex for the second triangle.
Returns
Number of common vertices for t and s.
uint select_quadrature_singquad2d ( pcsingquad2d  sq,
const uint tv,
const uint sv,
uint tp,
uint sp,
real **  x,
real **  y,
real **  w,
uint n,
real base 
)

This function is designed to select the correct quadrature rule for a current pair of triangles.

Parameters
sqA singquad2d object containing all necessary quadrature rules.
tvAn array defining the 3 vertices of triangle $ t $.
svAn array defining the 3 vertices of triangle $ s $.
tpReturning a permutation array of the vertices for $ t $.
spReturning a permutation array of the vertices for $ s $.
xReturning the quadrature points for the triangle $ t $.
yReturning the quadrature points for the triangle $ s $.
wReturning the quadrature weights.
nReturning the total number of quadrature points.
baseReturning a constant offset.
Returns
Returns the number of common vertices for triangle $ t $ and $ s $, which defines the current quadrature case.
void weight_basisfunc_cl_singquad2d ( real x,
real y,
real w,
uint  nq 
)

Weighting a quadrature rule for a double integral with a combination of piecewise constant and linear basis functions.

Parameters
xPoints within the first triangle using piecewise constant basis functions.
yPoints within the second triangle using linear basis functions.
wQuadrature weights. Actually this is a 4 times nq two dimensional array of weights. The first 3 columns will be used for the 3 combinations of constant-linear basis functions weights and the last column will be remain for the standard constant-constant case.
nqNumber of quadrature points and weights stored within x, y, w.
Attention
The column w[3] has to be filled with valid quadrature weights before calling this function.
void weight_basisfunc_l_singquad2d ( real x,
real y,
real w,
uint  nq 
)

Weighting a quadrature rule for a single integral with linear basis functions.

Parameters
xX-component of Points within the triangle.
yY-component of Points within the triangle.
wQuadrature weights. Actually this is a 4 times nq two dimensional array of weights. The first 3 columns will be used for the 3 combinations of linear basis functions weights and the last column will be remain for the standard constant case.
nqNumber of quadrature points and weights stored within x, y, w.
Attention
The column w[3] has to be filled with valid quadrature weights before calling this function.
void weight_basisfunc_lc_singquad2d ( real x,
real y,
real w,
uint  nq 
)

Weighting a quadrature rule for a double integral with a combination of piecewise linear and constant basis functions.

Parameters
xPoints within the first triangle using piecewise linear basis functions.
yPoints within the second triangle using piecewise constant basis functions.
wQuadrature weights. Actually this is a 4 times nq two dimensional array of weights. The first 3 columns will be used for the 3 combinations of linear-constant basis functions weights and the last column will be remain for the standard constant-constant case.
nqNumber of quadrature points and weights stored within x, y, w.
Attention
The column w[3] has to be filled with valid quadrature weights before calling this function.
void weight_basisfunc_ll_singquad2d ( real x,
real y,
real w,
uint  nq 
)

Weighting a quadrature rule for a double integral with both linear basis functions.

Parameters
xPoints within the first triangle.
yPoints within the second triangle.
wQuadrature weights. Actually this is a 10 times nq two dimensional array of weights. The first 9 columns will be used for the 9 combinations of linear-linear basis functions weights and the last column will be remain for the standard constant-constant case.
nqNumber of quadrature points and weights stored within x, y, w.
Attention
The column w[9] has to be filled with valid quadrature weights before calling this function.