Basix
Functions
basix::quadrature Namespace Reference

basix More...

Functions

Eigen::ArrayXXd compute_jacobi_deriv (double a, int n, int nderiv, const Eigen::ArrayXd &x)
 
Eigen::ArrayXd compute_gauss_jacobi_points (double a, int m)
 
std::pair< Eigen::ArrayXd, Eigen::ArrayXd > compute_gauss_jacobi_rule (double a, int m)
 Gauss-Jacobi quadrature rule (points and weights) More...
 
std::pair< Eigen::ArrayXd, Eigen::ArrayXd > make_quadrature_line (int m)
 
std::pair< Eigen::ArrayX2d, Eigen::ArrayXd > make_quadrature_triangle_collapsed (int m)
 
std::pair< Eigen::ArrayX3d, Eigen::ArrayXd > make_quadrature_tetrahedron_collapsed (int m)
 
std::pair< Eigen::ArrayXXd, Eigen::ArrayXd > make_quadrature (const std::string &rule, cell::type celltype, int m)
 
std::pair< Eigen::ArrayXd, Eigen::ArrayXd > make_gll_line (int m)
 
std::pair< Eigen::ArrayXd, Eigen::ArrayXd > compute_gll_rule (int m)
 GLL quadrature rule (points and weights)
 

Detailed Description

basix

Integration using Gauss-Jacobi quadrature on simplices. Other shapes can be obtained by using a product.

Todo:
  • pyramid

Function Documentation

◆ compute_gauss_jacobi_points()

Eigen::ArrayXd basix::quadrature::compute_gauss_jacobi_points ( double  a,
int  m 
)

Finds the m roots of \(P_{m}^{a,0}\) on [-1,1] by Newton's method.

Parameters
[in]aweight in Jacobi (b=0)
[in]morder
Returns
list of points in 1D

Computes the m roots of \(P_{m}^{a,0}\) on [-1,1] by Newton's method. The initial guesses are the Chebyshev points. Algorithm implemented from the pseudocode given by Karniadakis and Sherwin

◆ compute_gauss_jacobi_rule()

std::pair< Eigen::ArrayXd, Eigen::ArrayXd > basix::quadrature::compute_gauss_jacobi_rule ( double  a,
int  m 
)

Gauss-Jacobi quadrature rule (points and weights)

Note
Computes on [-1, 1]

◆ compute_jacobi_deriv()

Eigen::ArrayXXd basix::quadrature::compute_jacobi_deriv ( double  a,
int  n,
int  nderiv,
const Eigen::ArrayXd &  x 
)

Evaluate the nth Jacobi polynomial and derivatives with weight parameters (a, 0) at points x

Parameters
[in]aJacobi weight a
[in]nOrder of polynomial
[in]nderivNumber of derivatives (if zero, just compute polynomial itself)
[in]xPoints at which to evaluate
Returns
Array of polynomial derivative values (rows) at points (columns)

◆ make_gll_line()

std::pair< Eigen::ArrayXd, Eigen::ArrayXd > basix::quadrature::make_gll_line ( int  m)

Compute GLL line quadrature rule on [0, 1]

Parameters
morder
Returns
list of 1D points, list of weights

◆ make_quadrature()

std::pair< Eigen::ArrayXXd, Eigen::ArrayXd > basix::quadrature::make_quadrature ( const std::string &  rule,
cell::type  celltype,
int  m 
)

Utility for quadrature rule on reference cell

Parameters
[in]ruleName of quadrature rule (or use "default")
[in]celltype
[in]mMaximum degree of polynomial that this quadrature rule will integrate exactly
Returns
List of points and list of weights. The number of points arrays has shape (num points, gdim)

◆ make_quadrature_line()

std::pair< Eigen::ArrayXd, Eigen::ArrayXd > basix::quadrature::make_quadrature_line ( int  m)

Compute line quadrature rule on [0, 1]

Parameters
morder
Returns
list of points, list of weights

◆ make_quadrature_tetrahedron_collapsed()

std::pair< Eigen::ArrayX3d, Eigen::ArrayXd > basix::quadrature::make_quadrature_tetrahedron_collapsed ( int  m)

Compute tetrahedron quadrature rule on [0, 1]x[0, 1]x[0, 1]

Parameters
[in]morder
Returns
List of points, list of weights. The number of points arrays has shape (num points, gdim)

◆ make_quadrature_triangle_collapsed()

std::pair< Eigen::ArrayX2d, Eigen::ArrayXd > basix::quadrature::make_quadrature_triangle_collapsed ( int  m)

Compute triangle quadrature rule on [0, 1]x[0, 1]

Parameters
[in]morder
Returns
list of points, list of weights