sfepy.discrete.quadratures module

quadrature_tables are organized as follows:

quadrature_tables = {
    '<geometry1>' : {
        order1 : QuadraturePoints(args1),
        order2 : QuadraturePoints(args2),
        ...
    },
    '<geometry2>' : {
        order1 : QuadraturePoints(args1),
        order2 : QuadraturePoints(args2),
        ...
    },
    ...
}

Note The order for quadratures on tensor product domains (‘2_4’, ‘3_8’ geometries) in case of composite Gauss quadratures (products of 1D quadratures) holds for each component separately, so the actual polynomial order may be much higher (up to order * dimension).

Naming conventions in problem description files:

`<family>_<order>_<dimension>`

Integral ‘family’ is just an arbitrary name given by user.

Low order quadrature coordinates and weights copied from The Finite Element Method Displayed by Gouri Dhatt and Gilbert Touzat, Wiley-Interscience Production, 1984.

The line integral (geometry ‘1_2’) coordinates and weights are from Abramowitz, M. and Stegun, I.A., Handbook of Mathematical Functions, Dover Publications, New York, 1972. The triangle (geometry ‘2_3’) coordinates and weights are from Dunavant, D.A., High Degree Efficient Symmetrical Gaussian Quadrature Rules for the Triangle, Int. J. Num. Meth. Eng., 21 (1985) pp 1129-1148 - only rules with points inside the reference triangle are used. The actual values were copied from PHAML (http://math.nist.gov/phaml/), see also Mitchell, W.F., PHAML User’s Guide, NISTIR 7374, 2006.

Quadrature rules for the quadrilateral (geometry ‘2_4’) and hexahedron (geometry ‘3_8’) of order higher than 5 are computed as the tensor product of the line (geometry ‘1_2’) rules.

Quadrature rules for the triangle (geometry ‘2_3’) and tetrahedron (geometry ‘3_4’) of order higher than 19 and 6, respectively follow A. Grundmann and H.M. Moeller, Invariant integration formulas for the n-simplex by combinatorial methods, SIAM J. Numer. Anal. 15 (1978), 282–290. The generating function was adapted from pytools/hegde codes (http://mathema.tician.de/software/hedge) by Andreas Kloeckner.

class sfepy.discrete.quadratures.QuadraturePoints(data, coors=None, weights=None, bounds=None, tp_fix=1.0, weight_fix=1.0, symmetric=False)[source]

Representation of a set of quadrature points.

Parameters:
data : array_like

The array of shape (n_point, dim + 1) of quadrature point coordinates (first dim columns) and weights (the last column).

coors : array_like, optional

Optionally, instead of using data, the coordinates and weights can be provided separately - data are then ignored.

weights : array_like, optional

Optionally, instead of using data, the coordinates and weights can be provided separately - data are then ignored.

bounds : (float, float), optional

The coordinates and weights should correspond to a reference element in [0, 1] x dim. Provide the correct bounds if this is not the case.

tp_fix : float, optional

The value that is used to multiply the tensor product element volume (= 1.0) to get the correct volume.

weight_fix : float, optional

The value that is used to multiply the weights to get the correct values.

symmetric : bool

If True, the integral is 1D and the given coordinates and weights are symmetric w.r.t. the centre of bounds; only the non-negative coordinates are given.

static from_table(geometry, order)[source]

Create a new QuadraturePoints instance, given reference element geometry name and polynomial order. For tensor product geometries, the polynomial order is the 1D (line) order.

sfepy.discrete.quadratures.get_actual_order(geometry, order)[source]

Return the actual integration order for given geometry.

Parameters:
geometry : str

The geometry key describing the integration domain, see the keys of quadrature_tables.

Returns:
order : int

If order is in quadrature tables it is this value. Otherwise it is the closest higher order. If no higher order is available, a warning is printed and the highest available order is used.