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

Discrete Serendipity Hcurl space: local operators, L2 product and global interpolator. More...

#include <sxcurl.hpp>

Inheritance diagram for HArDCore2D::SXCurl:
Inheritance graph
[legend]
Collaboration diagram for HArDCore2D::SXCurl:
Collaboration graph
[legend]

Classes

struct  TransferOperators
 A structure to store the serendipity, extension and reduction operators. More...
 

Public Types

typedef std::function< Eigen::Vector2d(const Eigen::Vector2d &)> FunctionType
 

Public Member Functions

 SXCurl (const DDRCore &ddr_core, const SerendipityProblem &ser_pro, bool use_threads=true, std::ostream &output=std::cout)
 Constructor.
 
const Meshmesh () const
 Return the mesh.
 
const size_t & degree () const
 Return the polynomial degree.
 
Eigen::VectorXd interpolate (const FunctionType &v, const int deg_quad=-1) const
 Interpolator of a continuous function.
 
const Eigen::MatrixXd & ScurlCell (size_t iT) const
 Return the serendipity reconstruction for the cell of index iT.
 
const Eigen::MatrixXd & EcurlCell (size_t iT) const
 Return the extension for the cell of index iT.
 
const Eigen::MatrixXd & RcurlCell (size_t iT) const
 Return the reduction for the cell of index iT.
 
const Eigen::MatrixXd & ScurlCell (const Cell &T) const
 Return the serendipity reconstruction for cell T.
 
const Eigen::MatrixXd & EcurlCell (const Cell &T) const
 Return the extension for cell T.
 
const Eigen::MatrixXd & RcurlCell (const Cell &T) const
 Return the reduction for cell T.
 
const Eigen::MatrixXd cellCurl (size_t iT) const
 Return the full curl operator on the cell of index iT.
 
const Eigen::MatrixXd cellCurl (const Cell &T) const
 Return the full curl operator on cell T.
 
const Eigen::MatrixXd cellPotential (size_t iT) const
 Return the potential operator on the cell of index iT.
 
const Eigen::MatrixXd cellPotential (const Cell &T) const
 Return the potential operator on cell T.
 
Eigen::MatrixXd computeL2Product (const size_t iT, const double &penalty_factor=1., const Eigen::MatrixXd &mass_Pk2_T=Eigen::MatrixXd::Zero(1, 1), const IntegralWeight &weight=IntegralWeight(1.)) const
 Compute the matrix of the (weighted) L2-product.
 
Eigen::MatrixXd computeL2ProductGradient (const size_t iT, const SXGrad &sx_grad, const std::string &side, const double &penalty_factor=1., const Eigen::MatrixXd &mass_Pk2_T=Eigen::MatrixXd::Zero(1, 1), const IntegralWeight &weight=IntegralWeight(1.)) const
 Compute the matrix of the (weighted) L2-product as 'computeL2Product', with application of the discrete gradient on the left/right/both sides (depending on argument "side").
 
const DDRCore::CellBasescellBases (size_t iT) const
 Return cell bases for the face of index iT.
 
const DDRCore::CellBasescellBases (const Cell &T) const
 Return cell bases for cell T.
 
const DDRCore::EdgeBasesedgeBases (size_t iE) const
 Return edge bases for the edge of index iE.
 
const DDRCore::EdgeBasesedgeBases (const Edge &E) const
 Return edge bases for edge E.
 
- Public Member Functions inherited from HArDCore2D::VariableDOFSpace
 VariableDOFSpace (const Mesh &mesh, const Eigen::VectorXd n_local_vertex_dofs, const Eigen::VectorXd n_local_edge_dofs, const Eigen::VectorXd n_local_cell_dofs)
 Constructor.
 
 VariableDOFSpace (const Mesh &mesh, size_t n_local_vertex_dofs, const Eigen::VectorXd n_local_edge_dofs, const Eigen::VectorXd n_local_cell_dofs)
 Simpler constructor if all vertices have the same number of DOFs.
 
 VariableDOFSpace (const Mesh &mesh, size_t n_local_vertex_dofs, size_t n_local_edge_dofs, const Eigen::VectorXd n_local_cell_dofs)
 Simpler constructor if all vertices/edges have the same number of DOFs.
 
 VariableDOFSpace (const Mesh &mesh, size_t n_local_vertex_dofs, size_t n_local_edge_dofs, size_t n_local_cell_dofs)
 Simpler constructor if all vertices/edges/cells have the same number of DOFs.
 
const Meshmesh () const
 Returns the mesh.
 
size_t numLocalDofsVertex (const size_t iV) const
 Returns the number of local DOFs on vertex of index iV.
 
size_t numLocalDofsVertex (const Vertex &V) const
 Returns the number of local DOFs on vertex V.
 
size_t numLocalDofsEdge (const size_t iE) const
 Returns the number of local DOFs on edge of index iE.
 
size_t numLocalDofsEdge (const Edge &E) const
 Returns the number of local DOFs on edge E.
 
size_t numLocalDofsCell (const size_t iT) const
 Returns the number of local DOFs on cell of index iT.
 
size_t numLocalDofsCell (const Cell &T) const
 Returns the number of local DOFs on cell T.
 
size_t nDOFs_vertices () const
 Total number of vertices DOFs.
 
size_t nDOFs_edges () const
 Total number of edges DOFs.
 
size_t nDOFs_cells () const
 Total number of cells DOFs.
 
size_t dimension () const
 Returns the dimension of the global space (all DOFs for all geometric entities)
 
size_t dimensionVertex (const Vertex &V) const
 Returns the dimension of the local space on the vertex V.
 
size_t dimensionVertex (size_t iV) const
 Returns the dimension of the local space on the vertex of index iV.
 
size_t dimensionEdge (const Edge &E) const
 Returns the dimension of the local space on the edge E (including vertices)
 
size_t dimensionEdge (size_t iE) const
 Returns the dimension of the local space on the edge of index iE (including vertices)
 
size_t dimensionCell (const Cell &T) const
 Returns the dimension of the local space on the cell T (including faces, edges and vertices)
 
size_t dimensionCell (size_t iT) const
 Returns the dimension of the local space on the cell of index iT (including faces, edges and vertices)
 
size_t localOffset (const Edge &E, const Vertex &V) const
 Returns the local offset of the vertex V with respect to the edge E.
 
size_t localOffset (const Edge &E) const
 Returns the local offset of the unknowns attached to the edge E.
 
size_t localOffset (const Cell &T, const Vertex &V) const
 Returns the local offset of the vertex V with respect to the cell T.
 
size_t localOffset (const Cell &T, const Edge &E) const
 Returns the local offset of the edge E with respect to the cell T.
 
size_t localOffset (const Cell &T) const
 Returns the local offset of the unknowns attached to the element T.
 
size_t globalOffset (const Vertex &V) const
 Return the global offset for the unknowns on the vertex V.
 
size_t globalOffset (const Edge &E) const
 Return the global offset for the unknowns on the edge E.
 
size_t globalOffset (const Cell &T) const
 Return the global offset for the unknowns on the cell T.
 
Eigen::VectorXd restrictEdge (size_t iE, const Eigen::VectorXd &vh) const
 Restrict to the edge (including its vertices) of index iE.
 
Eigen::VectorXd restrictCell (size_t iT, const Eigen::VectorXd &vh) const
 Restrict to the cell (including vertices, edges and faces) of index iT.
 
Eigen::VectorXd restrict (const Edge &E, const Eigen::VectorXd vh) const
 Restrict to an edge.
 
Eigen::VectorXd restrict (const Cell &T, const Eigen::VectorXd vh) const
 Restrict to a cell.
 
Eigen::MatrixXd extendOperator (const Cell &T, const Edge &E, const Eigen::MatrixXd &opE) const
 Extend an edge operator to a cell.
 
std::vector< size_t > globalDOFIndices (const Cell &T) const
 Returns a vector listing the global DOFs attached to the element T: vertex DOFs, edge DOFs, face DOFs and element DOFs.
 
std::vector< size_t > globalDOFIndices (const Edge &E) const
 Returns a vector listing the global DOFs attached to the edge E: vertex DOFs, edge DOFs,.
 

Detailed Description

Discrete Serendipity Hcurl space: local operators, L2 product and global interpolator.

On each edge, the DOFs correspond to the polynomial bases on the edge provided by m_ddr_core. On each face/element, the DOFs are first those of the \(\mathcal{R}^{k-1}\) component and then of the \(\mathcal{R}^{c,\ell+1}\) component, each one of them following the bases of these spaces provided by m_ddr_core or restriction thereof.

Member Typedef Documentation

◆ FunctionType

typedef std::function<Eigen::Vector2d(const Eigen::Vector2d &)> HArDCore2D::SXCurl::FunctionType

Constructor & Destructor Documentation

◆ SXCurl()

SXCurl::SXCurl ( const DDRCore ddr_core,
const SerendipityProblem ser_pro,
bool  use_threads = true,
std::ostream &  output = std::cout 
)

Constructor.

Member Function Documentation

◆ cellBases() [1/2]

const DDRCore::CellBases & HArDCore2D::SXCurl::cellBases ( const Cell &  T) const
inline

Return cell bases for cell T.

◆ cellBases() [2/2]

const DDRCore::CellBases & HArDCore2D::SXCurl::cellBases ( size_t  iT) const
inline

Return cell bases for the face of index iT.

◆ cellCurl() [1/2]

const Eigen::MatrixXd HArDCore2D::SXCurl::cellCurl ( const Cell &  T) const
inline

Return the full curl operator on cell T.

◆ cellCurl() [2/2]

const Eigen::MatrixXd HArDCore2D::SXCurl::cellCurl ( size_t  iT) const
inline

Return the full curl operator on the cell of index iT.

◆ cellPotential() [1/2]

const Eigen::MatrixXd HArDCore2D::SXCurl::cellPotential ( const Cell &  T) const
inline

Return the potential operator on cell T.

◆ cellPotential() [2/2]

const Eigen::MatrixXd HArDCore2D::SXCurl::cellPotential ( size_t  iT) const
inline

Return the potential operator on the cell of index iT.

◆ computeL2Product()

Eigen::MatrixXd HArDCore2D::SXCurl::computeL2Product ( const size_t  iT,
const double &  penalty_factor = 1.,
const Eigen::MatrixXd &  mass_Pk2_T = Eigen::MatrixXd::Zero(1,1),
const IntegralWeight weight = IntegralWeight(1.) 
) const
inline

Compute the matrix of the (weighted) L2-product.

Parameters
iTindex of the cell
penalty_factorpre-factor for stabilisation term
mass_Pk2_Tif pre-computed, the mass matrix of (P^k(T))^2; if none is pre-computed, passing Eigen::MatrixXd::Zero(1,1) will force the calculation
weightweight function in the L2 product, defaults to 1

◆ computeL2ProductGradient()

Eigen::MatrixXd SXCurl::computeL2ProductGradient ( const size_t  iT,
const SXGrad sx_grad,
const std::string &  side,
const double &  penalty_factor = 1.,
const Eigen::MatrixXd &  mass_Pk2_T = Eigen::MatrixXd::Zero(1,1),
const IntegralWeight weight = IntegralWeight(1.) 
) const

Compute the matrix of the (weighted) L2-product as 'computeL2Product', with application of the discrete gradient on the left/right/both sides (depending on argument "side").

Parameters
iTindex of the cell
sx_gradinstance of SXGrad to access the full gradients
sidewhich side (left, right, both) we apply the gradient to
penalty_factorpre-factor for stabilisation term
mass_Pk2_Tif pre-computed, the mass matrix of (P^k(T))^2; if none is pre-computed, passing Eigen::MatrixXd::Zero(1,1) will force the calculation
weightweight function in the L2 product, defaults to constant 1.

◆ degree()

const size_t & HArDCore2D::SXCurl::degree ( ) const
inline

Return the polynomial degree.

◆ EcurlCell() [1/2]

const Eigen::MatrixXd & HArDCore2D::SXCurl::EcurlCell ( const Cell &  T) const
inline

Return the extension for cell T.

◆ EcurlCell() [2/2]

const Eigen::MatrixXd & HArDCore2D::SXCurl::EcurlCell ( size_t  iT) const
inline

Return the extension for the cell of index iT.

◆ edgeBases() [1/2]

const DDRCore::EdgeBases & HArDCore2D::SXCurl::edgeBases ( const Edge &  E) const
inline

Return edge bases for edge E.

◆ edgeBases() [2/2]

const DDRCore::EdgeBases & HArDCore2D::SXCurl::edgeBases ( size_t  iE) const
inline

Return edge bases for the edge of index iE.

◆ interpolate()

Eigen::VectorXd SXCurl::interpolate ( const FunctionType v,
const int  deg_quad = -1 
) const

Interpolator of a continuous function.

Parameters
vThe function to interpolate
deg_quadThe optional degre of quadrature rules to compute the interpolate. If negative, then 2*degree()+3 will be used.

◆ mesh()

const Mesh & HArDCore2D::SXCurl::mesh ( ) const
inline

Return the mesh.

◆ RcurlCell() [1/2]

const Eigen::MatrixXd & HArDCore2D::SXCurl::RcurlCell ( const Cell &  T) const
inline

Return the reduction for cell T.

◆ RcurlCell() [2/2]

const Eigen::MatrixXd & HArDCore2D::SXCurl::RcurlCell ( size_t  iT) const
inline

Return the reduction for the cell of index iT.

◆ ScurlCell() [1/2]

const Eigen::MatrixXd & HArDCore2D::SXCurl::ScurlCell ( const Cell &  T) const
inline

Return the serendipity reconstruction for cell T.

◆ ScurlCell() [2/2]

const Eigen::MatrixXd & HArDCore2D::SXCurl::ScurlCell ( size_t  iT) const
inline

Return the serendipity reconstruction for the cell of index iT.


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