HArD::Core2D
Hybrid Arbitrary Degree::Core 2D - Library to implement 2D schemes with edge and cell polynomials as unknowns
Classes | Public Types | Public Member Functions | List of all members
HArDCore2D::HybridCore Class Reference

#include <hybridcore.hpp>

Classes

struct  qrule
 description of one node and one weight from a quadrature rule More...
 

Public Types

using cell_basis_type = std::function< double(double, double)>
 type for cell basis
 
using cell_gradient_type = std::function< Eigen::Vector2d(double, double)>
 type for gradients of cell basis
 
using edge_basis_type = std::function< double(double, double)>
 type for edge basis
 
using tensor_function_type = std::function< Eigen::Matrix2d(double, double)>
 type for 2D tensors basis
 

Public Member Functions

 HybridCore (const Mesh *mesh_ptr, const size_t K, const size_t L, const std::string orth="noON")
 Class constructor: initialises the data structure with the given mesh, and desired polynomial degrees of the basis functions. More...
 
size_t dim_Pcell (const size_t m) const
 Compute the size of the basis of 2-variate polynomials up to degree m. More...
 
size_t dim_Pedge (const size_t m) const
 Compute the size of the basis of 1-variate polynomials up to degree m. More...
 
const cell_basis_typecell_basis (size_t iT, size_t i) const
 Return a reference to the i'th basis function of the cell iT. More...
 
const edge_basis_typeedge_basis (size_t iF, size_t i) const
 Return a reference to the i'th basis function of the edge iF. More...
 
const cell_gradient_typecell_gradient (size_t iT, size_t i) const
 Return a reference to the gradient of the i'th basis function of the celliT. More...
 
const Meshget_mesh_ptr () const
 returns a pointer to the mesh
 
std::vector< HybridCore::qrulecell_qrule (const size_t iT, const size_t doe) const
 Compute a quadrature rule of a given degree of exactness on a cell. More...
 
std::vector< HybridCore::qruleedge_qrule (const size_t iE, const size_t doe) const
 Compute a quadrature rule of a given degree of exactness on an edge. More...
 
Eigen::VectorXd restr (const Eigen::VectorXd &Xh, size_t iT) const
 Extract from a global vector Xh of unknowns the unknowns corresponding to cell iT.
 
double L2norm (const Eigen::VectorXd &Xh) const
 Compute L2 norm of a discrete function (using cell values)
 
double H1norm (const Eigen::VectorXd &Xh) const
 Compute discrete H1 norm of a discrete function (using cell values)
 
double Linf_edge (const Eigen::VectorXd &Xh) const
 Compute maximum of the coefficients on the edge basis functions.
 
template<typename Function >
Eigen::VectorXd interpolate (const Function &f, size_t doe) const
 Compute the interpolant in the discrete space of a continuous function. More...
 
Eigen::MatrixXd gram_matrix (const std::vector< Eigen::ArrayXd > &f_quad, const std::vector< Eigen::ArrayXd > &g_quad, const size_t &nrows, const size_t &ncols, const std::vector< HybridCore::qrule > &quad, const bool &sym, std::vector< double > L2weight={}) const
 
Eigen::MatrixXd gram_matrix (const std::vector< Eigen::ArrayXXd > &F_quad, const std::vector< Eigen::ArrayXXd > &G_quad, const size_t &nrows, const size_t &ncols, const std::vector< HybridCore::qrule > &quad, const bool &sym, std::vector< Eigen::Matrix2d > L2Weight={}) const
 Overloaded version of the previous one for vector-valued functions: the functions (F_i) and (G_j) are vector-valued functions. More...
 
Eigen::VectorXd compute_weights (size_t iT) const
 Weights to compute cell unknowns from edge unknowns when l=-1.
 
const std::vector< Eigen::ArrayXd > basis_quad (const std::string celledge, const size_t iTF, const std::vector< qrule > quad, const size_t nbasis) const
 Computes (cell or edge) basis functions at the given quadrature nodes. More...
 
const std::vector< Eigen::ArrayXXd > grad_basis_quad (const size_t iT, const std::vector< qrule > quad, const size_t nbasis) const
 Compute $(\nabla \phi_i)_{i\in I}$ at the given quadrature nodes, where $(\phi_i)_{i\in I}$ are the cell basis functions. More...
 
double integral (const Eigen::VectorXd &XTF) const
 Computes the integral of a discrete function using cell unknowns.
 
double evaluate_in_cell (const Eigen::VectorXd XTF, size_t iT, double x, double y) const
 Evaluates a discrete function in the cell iT at point (x,y)
 
double evaluate_in_edge (const Eigen::VectorXd XTF, size_t iF, double x, double y) const
 Evaluates a discrete function on the edge iF at point (x,y)
 
const size_t K ()
 polynomial degree of edge unknowns
 
const int L ()
 polynomial degree of cell unknowns
 
const int Ldeg ()
 usually equal to L, but put at 0 if L=-1
 
size_t ntotal_dofs () const
 Total number of degrees of freedom.
 
size_t nlocal_cell_dofs () const
 number of degrees of freedom in each cell (dimension of polynomial space)
 
size_t ntotal_cell_dofs () const
 total number of cell degrees of freedom
 
size_t nlocal_edge_dofs () const
 number of degrees of freedom on each edge (dimension of polynomial space)
 
size_t ntotal_edge_dofs () const
 total number of edge degrees of freedom
 
size_t ninternal_edge_dofs () const
 total number of edge degrees of freedom for internal edges
 
size_t nboundary_edge_dofs () const
 total number of edge degrees of freedom for boundary edges
 
size_t nhighorder_dofs () const
 total number of cell degrees of freedom with polynomials up to order k+1
 
size_t ngradient_dofs () const
 total number of degrees of freedom for gradients
 
template<typename Function >
void quadrature_over_cell (const size_t iT, const Function &f) const
 To integrate a function over a cell.
 
template<typename Function >
void quadrature_over_edge (const size_t iF, const Function &f) const
 To integrate a function over an edge.
 
template<typename Function >
double integrate_over_cell (const size_t iT, const Function &f) const
 Integrates a function over a cell. Use with parcimony, expensive (re-compute quadratures)
 
template<typename Function >
double integrate_over_edge (const size_t iF, const Function &f) const
 Integrates a function over an edge. Use with parcimony, expensive (re-compute quadratures)
 
template<typename Function >
double integrate_over_domain (const Function &f) const
 Integrates a function over the domaine. Use with parcimony, expensive (re-compute quadratures)
 
Eigen::VectorXd VertexValues (const Eigen::VectorXd Xh, const std::string from_dofs)
 From a hybrid function, computes a vector of values at the vertices of the mesh. More...
 

Detailed Description

The HybridCore class provides convenient interfaces for performing integration over mesh cells and edges and handling polynomial basis functions The class also provides convenient interfaces for dealing with solutions to Hybrid High-Order schemes, such as the computation of integrals, norms and interpolants in the HHO space.

Constructor & Destructor Documentation

◆ HybridCore()

HybridCore::HybridCore ( const Mesh mesh_ptr,
const size_t  K,
const size_t  L,
const std::string  orth = "noON" 
)

Class constructor: initialises the data structure with the given mesh, and desired polynomial degrees of the basis functions.

The orthonormalisation comes at a cost in terms of manipulation of the basis functions. This should only be use when the polynomial degree is large and/or the cell is very distorted. However, in these cases, it can make a huge difference on the observed convergence rate.

Parameters
mesh_ptrA pointer to the loaded mesh
KThe degree of the edge polynomials
LThe degree of the cell polynomials
orth=ON to orthonormalise the basis functions.

Member Function Documentation

◆ basis_quad()

const std::vector< Eigen::ArrayXd > HybridCore::basis_quad ( const std::string  celledge,
const size_t  iTF,
const std::vector< qrule quad,
const size_t  nbasis 
) const

Computes (cell or edge) basis functions at the given quadrature nodes.

Returns
phi_quad[i] = array listing the nbq (=nb of quadrature nodes) values of phi_i at the quadrature nodes
Parameters
celledgedetermines the type of basis function (cell or edge) we want the values of
iTFglobal index of the cell/edge
quadquadrature nodes and weights on the cell/edge
nbasisnumber of basis functions to consider

◆ cell_basis()

const HybridCore::cell_basis_type & HybridCore::cell_basis ( size_t  iT,
size_t  i 
) const

Return a reference to the i'th basis function of the cell iT.

Parameters
iTThe global cell number of the cell
iThe index of the desired basis function

◆ cell_gradient()

const HybridCore::cell_gradient_type & HybridCore::cell_gradient ( size_t  iT,
size_t  i 
) const

Return a reference to the gradient of the i'th basis function of the celliT.

Note that the gradient functions are indexed the same as the basis functions. In particular, this means that the first gradient function will always be identically zero, as it is the gradient of the constant basis function.

Parameters
iTThe global cell number of the cell
iThe index of the desired basis function

◆ cell_qrule()

std::vector< HybridCore::qrule > HybridCore::cell_qrule ( const size_t  iT,
const size_t  doe 
) const

Compute a quadrature rule of a given degree of exactness on a cell.

Returns
list of quadrature nodes and weights
Parameters
iTGlobal index of the cell
doeRequired degree of exactness

◆ dim_Pcell()

size_t HybridCore::dim_Pcell ( const size_t  m) const

Compute the size of the basis of 2-variate polynomials up to degree m.

Parameters
mThe maximum degree of the polynomial

◆ dim_Pedge()

size_t HybridCore::dim_Pedge ( const size_t  m) const

Compute the size of the basis of 1-variate polynomials up to degree m.

Parameters
mThe maximum degree of the polynomial

◆ edge_basis()

const HybridCore::edge_basis_type & HybridCore::edge_basis ( size_t  iF,
size_t  i 
) const

Return a reference to the i'th basis function of the edge iF.

Parameters
iFThe global edge number of the edge
iThe index of the desired basis function

◆ edge_qrule()

std::vector< HybridCore::qrule > HybridCore::edge_qrule ( const size_t  iE,
const size_t  doe 
) const

Compute a quadrature rule of a given degree of exactness on an edge.

Returns
list of quadrature nodes and weights
Parameters
iEGlobal index of the edge
doeRequired degree of exactness

◆ grad_basis_quad()

const std::vector< Eigen::ArrayXXd > HybridCore::grad_basis_quad ( const size_t  iT,
const std::vector< qrule quad,
const size_t  nbasis 
) const

Compute $(\nabla \phi_i)_{i\in I}$ at the given quadrature nodes, where $(\phi_i)_{i\in I}$ are the cell basis functions.

Returns
dphi_quad[i]: array of size 2*nbq (where nbq=nb of quadrature nodes), with each column being $\nabla \phi_i$ at the corresponding quadrature node
Parameters
iTglobal index of the cell
quadquadrature rules in the cell
nbasisnumber of basis functions $phi_i$ to consider

◆ gram_matrix() [1/2]

Eigen::MatrixXd HybridCore::gram_matrix ( const std::vector< Eigen::ArrayXd > &  f_quad,
const std::vector< Eigen::ArrayXd > &  g_quad,
const size_t &  nrows,
const size_t &  ncols,
const std::vector< HybridCore::qrule > &  quad,
const bool &  sym,
std::vector< double >  L2weight = {} 
) const

Create the matrix of L2 products of two families (f_i) and (g_j) of functions (this is not really a Gram matrix, unless the two families are the same)

Returns
The matrix $(\int f_i g_j)_{i=1\ldots nrows; j=1\ldots ncols}$
Parameters
f_quadValues of functions (f1,f2,...) at the quadrature nodes
g_quadValues of functions (g1,g2,...) at the quadrature nodes
nrowsNumber of rows of the matrix - typically number of functions f_i (but could be less)
ncolsNumber of columns of the matrix - typically number of functions g_j (but could be less)
quadQuadrature nodes for integration
symTrue if the matrix is pseudo-symmetric (that is, #f<=#g and f_i=g_i if i<=#f)
L2weightOptional weight for the L2 product. If provided, should be a std::vector<double> of the weight at the quadrature nodes

◆ gram_matrix() [2/2]

Eigen::MatrixXd HybridCore::gram_matrix ( const std::vector< Eigen::ArrayXXd > &  F_quad,
const std::vector< Eigen::ArrayXXd > &  G_quad,
const size_t &  nrows,
const size_t &  ncols,
const std::vector< HybridCore::qrule > &  quad,
const bool &  sym,
std::vector< Eigen::Matrix2d >  L2Weight = {} 
) const

Overloaded version of the previous one for vector-valued functions: the functions (F_i) and (G_j) are vector-valued functions.

Returns
The matrix $(\int F_i \cdot G_j)_{i=1\ldots nrows; j=1\ldots ncols}$
Parameters
F_quadValues of functions (F1,F2,...) at the quadrature nodes
G_quadValues of functions (G1,G2,...) at the quadrature nodes
nrowsNumber of rows of the matrix - typically number of functions F_i (but could be less)
ncolsNumber of rows of the matrix - typically number of functions G_j (but could be less)
quadQuadrature nodes for integration
symTrue if the matrix is pseudo-symmetric (that is, #F<=#G and F_i=G_i if i<=#F)
L2WeightOptional weight for the L2 product. If provided, should be a std::vector<Eigen::Matrix2d> of the weight at the quadrature nodes

◆ VertexValues()

Eigen::VectorXd HybridCore::VertexValues ( const Eigen::VectorXd  Xh,
const std::string  from_dofs 
)

From a hybrid function, computes a vector of values at the vertices of the mesh.

Parameters
Xhhybrid function (cell and edge polynomials)
from_dofsType of unknowns to use: "cell" or "edge"

The documentation for this class was generated from the following files: