modepy.
Quadrature
(nodes, weights)¶The basic interface for a quadrature rule.
nodes
¶An array of shape (d, nnodes), where d is the dimension of the qudrature rule. In 1D, the shape is just (nnodes,).
weights
¶A vector of length nnodes.
__call__
(f)¶Evaluate the callable f at the quadrature nodes and return its integral.
f is assumed to accept arrays of shape (dims, npts), or of shape (npts,) for 1D quadrature.
modepy.
Quadrature
(nodes, weights)The basic interface for a quadrature rule.
nodes
An array of shape (d, nnodes), where d is the dimension of the qudrature rule. In 1D, the shape is just (nnodes,).
weights
A vector of length nnodes.
modepy.
JacobiGaussQuadrature
(alpha, beta, N)¶An Gauss quadrature of order N associated with the Jacobi weight \((1-x)^\alpha(1+x)^\beta\).
Assumes \(\alpha,\beta > -1\) with \((\alpha,\beta)\not\in\{(-1/2,-1/2)\}\).
Integrates on the interval (-1,1). The quadrature rule is exact up to degree \(2N+1\).
Inherits from modepy.Quadrature
. See there for the interface
to obtain nodes and weights.
modepy.
LegendreGaussQuadrature
(N)¶An Gauss quadrature associated with weight 1.
Integrates on the interval (-1,1). The quadrature rule is exact up to degree \(2N+1\).
Inherits from modepy.Quadrature
. See there for the interface
to obtain nodes and weights.
modepy.quadrature.jacobi_gauss.
jacobi_gauss_lobatto_nodes
(alpha, beta, N)¶Compute the Gauss-Lobatto quadrature
nodes corresponding to the JacobiGaussQuadrature
with the same parameters.
Exact to degree \(2N-3\).
modepy.quadrature.jacobi_gauss.
legendre_gauss_lobatto_nodes
(N)¶Compute the Legendre-Gauss-Lobatto quadrature nodes.
Exact to degree \(2N-1\).
New in version 2013.3.
modepy.
GrundmannMoellerSimplexQuadrature
(order, dimension)¶Cubature on an n-simplex.
This cubature rule has both negative and positive weights. It is exact for polynomials up to order \(2s+1\), where \(s\) is given as order.
The integration domain is the unit simplex. (see Coordinates on the triangle and Coordinates on the tetrahedron)
exact_to
¶The total degree up to which the quadrature is exact.
See
Parameters: |
|
---|
__call__
(f)¶Evaluate the callable f at the quadrature nodes and return its integral.
f is assumed to accept arrays of shape (dims, npts), or of shape (npts,) for 1D quadrature.
modepy.
XiaoGimbutasSimplexQuadrature
(order, dims)¶A (nearly) Gaussian simplicial quadrature with very few quadrature nodes, available for low-to-moderate orders.
Raises modepy.QuadratureRuleUnavailable
if no quadrature rule for the
requested parameters is available.
The integration domain is the unit simplex. (see Coordinates on the triangle and Coordinates on the tetrahedron)
Inherits from modepy.Quadrature
. See there for the interface
to obtain nodes and weights.
exact_to
¶The total degree up to which the quadrature is exact.
See
H. Xiao and Z. Gimbutas, “A numerical algorithm for the construction of efficient quadrature rules in two and higher dimensions,” Computers & Mathematics with Applications, vol. 59, no. 2, pp. 663-676, 2010. http://dx.doi.org/10.1016/j.camwa.2009.10.027
Parameters: |
|
---|
__call__
(f)¶Evaluate the callable f at the quadrature nodes and return its integral.
f is assumed to accept arrays of shape (dims, npts), or of shape (npts,) for 1D quadrature.
modepy.
VioreanuRokhlinSimplexQuadrature
(order, dims)¶Simplicial quadratures with symmetric node sets and positive weights suitable for well-conditioned interpolation. The integration domain is the unit simplex. (see Coordinates on the triangle and Coordinates on the tetrahedron)
Raises modepy.QuadratureRuleUnavailable
if no quadrature rule for the
requested parameters is available.
Inherits from modepy.Quadrature
. See there for the interface
to obtain nodes and weights.
exact_to
¶The total degree up to which the quadrature is exact.
When using these nodes, please acknowledge Zydrunas Gimbutas, who generated them as follows:
The 2D nodes are based on the interpolation node set derived in the article
B. Vioreanu and V. Rokhlin, “Spectra of Multiplication Operators as a Numerical Tool,” Yale CS Tech Report 1443
Note that in Vioreanu’s tables, only orders 5,6,9, and 12 are rotationally symmetric, which gives one extra order for integration and better interpolation conditioning. Also note that since the tables have been re-generated independently, the nodes and weights may be different.
The 3D nodes were derived from the modepy.warp_and_blend_nodes()
.
A tightening algorithm was then applied, as described in
B. Vioreanu, “Spectra of Multiplication Operators as a Numerical Tool”, Yale University, 2012. Dissertation
New in version 2013.3.
Parameters: |
|
---|
__call__
(f)¶Evaluate the callable f at the quadrature nodes and return its integral.
f is assumed to accept arrays of shape (dims, npts), or of shape (npts,) for 1D quadrature.