Basis based on Gauss-Lobatto quadrature. More...
#include <Gauss_lobatto.hpp>
Public Member Functions | |
Gauss_lobatto (int row_size_arg) | |
see hexed::Basis::Basis . | |
double | node (int i) const override |
the i th interpolation node (i.e. quadrature point). Should be in interval [0, 1]. | |
Eigen::VectorXd | node_weights () const override |
node_weights()(i) is the quadrature weight associated with node(i) Integrates in domain [0, 1]. | |
Eigen::MatrixXd | diff_mat () const override |
Differentiation matrix. | |
Eigen::MatrixXd | boundary () const override |
boundary()(i, j) is the j th basis polynomial evaluated at i for i in {0, 1}. | |
Eigen::VectorXd | orthogonal (int degree) const override |
orthogonal(degree)(i) is the degree th-degree Legendre polynomial evaluated at node(i) . | |
Eigen::MatrixXd | filter () const override |
Matrix that applies a low-pass filter to the Legendre modes of a polynomial. | |
Eigen::MatrixXd | prolong (int i_half) const override |
Matrix to prolong into polynomial space of one higher refinement level. | |
Eigen::MatrixXd | restrict (int i_half) const override |
Restrict from refined space above to polynomial space spanned by this basis. | |
double | min_eig_diffusion () const override |
nondimensional minimum eigenvalue of 1D diffusion operator | |
![]() | |
Basis (int row_size_arg) | |
Mat | nodes () const |
Get vector of nodes. | |
Mat< dyn, dyn > | interpolate (const Mat<> &sample) const |
Interpolate to arbitrary points. | |
double | max_cfl () const |
maximum stable CFL number for 1D convection | |
double | step_ratio () const |
ratio of time step for convection stage 2 to stage 1 | |
Protected Member Functions | |
double | min_eig_convection () const override |
nondimensional minimum real part of eigenvals of 1D convection operator | |
double | quadratic_safety () const override |
Additional Inherited Members | |
![]() | |
const int | row_size |
row size | |
Basis based on Gauss-Lobatto quadrature.
The node
s and node_weights
of this basis are those of the Gauss-Legendre quadrature, which is exact for polynomials of degree 2*row_size - 3
.
|
overridevirtual |
boundary()(i, j)
is the j
th basis polynomial evaluated at i
for i
in {0, 1}.
Implements hexed::Basis.
|
overridevirtual |
Differentiation matrix.
diff_mat()(i, j)
is the derivative of interpolating polynomial j
evaluated at node(i)
.
Implements hexed::Basis.
|
overridevirtual |
Matrix that applies a low-pass filter to the Legendre modes of a polynomial.
Implements hexed::Basis.
|
overrideprotectedvirtual |
nondimensional minimum real part of eigenvals of 1D convection operator
Implements hexed::Basis.
|
overridevirtual |
nondimensional minimum eigenvalue of 1D diffusion operator
Implements hexed::Basis.
|
overridevirtual |
the i
th interpolation node (i.e. quadrature point). Should be in interval [0, 1].
Implements hexed::Basis.
|
overridevirtual |
node_weights()(i)
is the quadrature weight associated with node(i)
Integrates in domain [0, 1].
Implements hexed::Basis.
|
overridevirtual |
orthogonal(degree)(i)
is the degree
th-degree Legendre polynomial evaluated at node(i)
.
Polynomial is scaled so that its norm is 1 with respect to the quadrature rule of the basis.
Implements hexed::Basis.
|
overridevirtual |
Matrix to prolong into polynomial space of one higher refinement level.
prolong(i_half)(i, j)
is the j
th basis polynomial evaluated at the i
th node in the refined space. If i_half
is 0, refined space is interval [0, 0.5]. If i_half
is 1, then [0.5, 1]. This is a basic, naive implementation. Derived classes may override for performance and/or precision.
Reimplemented from hexed::Basis.
|
inlineoverrideprotectedvirtual |
Implements hexed::Basis.
|
overridevirtual |
Restrict from refined space above to polynomial space spanned by this basis.
restrict(i_half)(i, j)
is the value at the i
th standard node of the orthogonal projection of the j
th polynomial in the refined space into the space of this basis.
Implements hexed::Basis.