hp3D::Space Class Referenceabstract
A 3D hp FEM space with continuous, picewise polynomial basis functions. More...
#include <space.hh>
Public Types | |
typedef concepts::Scan< hp3D::Element< Real > > | Scan |
typedef Scan< ElementWithCell< Real > > | Scanner |
typedef hp3D::TraceSpace | TraceSpaceT |
typedef ElementWithCell< Real > | type |
Public Member Functions | |
virtual void | adjust (const concepts::Element< Real > &elm, const concepts::AdaptiveAdjustP< 3 > &a) |
virtual void | adjust (const Element< Real > &elm, const concepts::AdaptiveAdjustP< 3 > &a)=0 |
Adjusts the space in the next rebuild step for this element. More... | |
const concepts::BoundaryConditions * | bc () |
virtual const concepts::BoundaryConditions * | bc () const |
uint | dim () |
Returns the dimension of the space. More... | |
virtual uint | dim () const |
virtual uint | getOutputDimension () const |
Returns the default output dimension, when we consider plotting a real-valued operator on this space. More... | |
void | getPmax (const concepts::Hexahedron &cntr, ushort Pmax[3]) const |
const std::unique_ptr< concepts::SMatrixBase< Real > > * | getSmatrices (const concepts::Hex3dSubdivision *subdiv) const |
Returns the list of matrices matching the subdivision type subdiv . More... | |
uint | nelm () |
Returns the number of elements in the space. More... | |
virtual uint | nelm () const |
Space | rebuild () |
Rebuilds the space after an adjustment with adjust . More... | |
void | recomputeShapefunctions () |
Recompute shape functions, e.g. More... | |
void | recomputeSmatrices (const Hexahedron *elm) const |
Checks if the S matrices need to be recomputed and does so if necessary. More... | |
Scan * | scan () |
Returns a scanner to iterate over the elements of the space. More... | |
virtual Scan * | scan () const |
void | setDim (uint dim) |
Space (concepts::Mesh3 &msh, uint l, uint p, concepts::BoundaryConditions *bc=0, concepts::CellConditions *cc=0) | |
Constructor. More... | |
Space (const Space &spc) | |
Copy constructor. More... | |
virtual | ~Space () |
Static Public Member Functions | |
Timing Interface | |
These functions are used to get timings from class internal computations. The values are stored in a user defined concepts::InOutParameters structure in different arrays (see concepts::ResultsTable table;
std::ofstream ofs("table.gnuplot");
ofs << std::setprecision(20);
table.print<concepts::ResultsTable::GNUPLOT>(ofs);
| |
static void | setTimings (concepts::InOutParameters *timings) |
Sets the class to store the timing values in. More... | |
static bool | timings () |
Returns true if the class is able to do timings. More... | |
Friends | |
class | BuildDofsBase |
Strategies to Build the Degrees of Freedom | |
concepts::Mesh3 & | msh_ |
Mesh. More... | |
concepts::BoundaryConditions * | bc_ |
Boundary conditions. More... | |
concepts::CellConditions * | cc_ |
Cell conditions. More... | |
bool | rebuild_ |
If true: the elements have to be rebuilt. More... | |
uint | dim_ |
Dimension of the FE space. More... | |
uint | nelm_ |
Number of elements currently active in the mesh. More... | |
concepts::Joiner< Element< Real > *, 1 > * | elm_ |
Linked list of the elements. More... | |
concepts::Joiner< Element< Real > *, 1 > * | endElm_ |
Mesh. More... | |
std::unordered_map< uint, concepts::AdaptiveControl<> > | ctrl0_ |
Hash table of the control information for the vertices. More... | |
std::unique_ptr< std::unordered_map< uint, concepts::AdaptiveControlP< 1 > > > | ctrl1_ |
Hash table of the control information for the edges. More... | |
std::unique_ptr< std::unordered_map< uint, concepts::AdaptiveControlP< 2 > > > | ctrl2_ |
Hash table of the control information for the faces. More... | |
std::unordered_map< uint, concepts::AdaptiveControlP< 3 > > | ctrl3_ |
Hash table of the control information for the 3D elements. More... | |
std::unordered_map< uint, concepts::AdaptiveControlP< 3 > > | pMax_ |
Hash table of the dimension of the shape function space for the 3D elements. More... | |
std::map< uint, concepts::CellData > | cellList_ |
List of cells. More... | |
std::map< uint, concepts::FaceData > | faceList_ |
List of faces. More... | |
std::map< uint, concepts::EdgeData > | edgeList_ |
List of edges. More... | |
std::map< uint, concepts::VertexData > | vertexList_ |
List of vertices. More... | |
std::unordered_map< uint, concepts::AdaptiveAdjustP< 3 > > | adj_ |
Hash table of the adjustment information of the elements. More... | |
std::unique_ptr< concepts::SMatrix1D > | S1left_ |
The S matrices in 1D. More... | |
std::unique_ptr< concepts::SMatrix1D > | S1right_ |
Mesh. More... | |
std::unique_ptr< concepts::SMatrixBase< Real > > | Smatrices8_ [8] |
The S matrices in 3D for the subdivision into 8 children. More... | |
std::unique_ptr< concepts::SMatrixBase< Real > > | Smatrices4x_ [4] |
The S matrices in 3D for the subdivision into 4 children. More... | |
std::unique_ptr< concepts::SMatrixBase< Real > > | Smatrices4y_ [4] |
Mesh. More... | |
std::unique_ptr< concepts::SMatrixBase< Real > > | Smatrices4z_ [4] |
Mesh. More... | |
std::unique_ptr< concepts::SMatrixBase< Real > > | Smatrices2x_ [2] |
The S matrices in 3D for the subdivision into 2 children. More... | |
std::unique_ptr< concepts::SMatrixBase< Real > > | Smatrices2y_ [2] |
Mesh. More... | |
std::unique_ptr< concepts::SMatrixBase< Real > > | Smatrices2z_ [2] |
Mesh. More... | |
std::unique_ptr< BuildDofsBase > | buildVertexDofs_ |
Strategy to build the vertex degrees of freedom. More... | |
std::unique_ptr< BuildDofsBase > | buildEdgeDofs_ |
Strategy to build the edge degrees of freedom. More... | |
std::unique_ptr< BuildDofsBase > | buildFaceDofs_ |
Strategy to build the face degrees of freedom. More... | |
std::unique_ptr< BuildDofsBase > | buildInnerDofs_ |
Strategy to build the inner degrees of freedom. More... | |
static concepts::InOutParameters * | timings_ |
Place to store timing values. More... | |
static uint | timeCntr_ |
Counter for timing table. More... | |
void | buildVertexDofs (const BuildDofsBase &b) |
Change the strategy how the degrees of freedom for the vertices are built. More... | |
void | buildEdgeDofs (const BuildDofsBase &b) |
Change the strategy how the degrees of freedom for the edge are built. More... | |
void | buildFaceDofs (const BuildDofsBase &b) |
Change the strategy how the degrees of freedom for the face are built. More... | |
void | buildInnerDofs (const BuildDofsBase &b) |
Change the strategy how the degrees of freedom for the interior are built. More... | |
virtual std::ostream & | info (std::ostream &os) const |
Mesh. More... | |
void | meshAndPoly_ (concepts::Hexahedron3d &cell, concepts::Level< 3 > &L, int *P) |
Adjusts the mesh and the polynomial degree. More... | |
void | enforceBC_ (concepts::Hexahedron3d &cell) |
Enforce Dirichlet boundary conditions. More... | |
void | enforceCC_ (concepts::Hexahedron3d &cell) |
Enforce cell condition INACTIVEPLUS. More... | |
void | createCellList_ (const concepts::Connector &cntr, const concepts::CellData *father) |
Initially fills cellList_ (recusively). More... | |
void | createEdgFaceList_ (const concepts::Cell3 &cell) |
Initially fills faceList_ and edgeList_ . More... | |
void | createVtxEdgList_ (const concepts::Cell3 &cell) |
Initially fills vertexList_ . More... | |
void | enrichElm_ (const concepts::Hexahedron3d &cell, int *p) |
Enriches the polynomial degrees. More... | |
void | buildElements_ (concepts::Hexahedron3d &cell, ushort *Pmax, concepts::TColumn< Real > *T0=0) |
Creates the elements and their T matrices. More... | |
void | getPmax_ (const concepts::Hexahedron &cntr, ushort Pmax[3]) const |
Returns the polynomial degrees of the shape function space on connector. More... | |
void | computePmax_ (const concepts::Hexahedron &cntr, ushort Pmax[3]) const |
Computes the maximal polynomial degree of cell in each direction. More... | |
void | deactivate_ (const concepts::Connector2 &face) |
Deactivates the children of face (including the new edges and vertices). More... | |
void | deactivate_ (const concepts::Connector1 &edg) |
Deactivates the children of the edge edg (including the middle vertex). More... | |
void | deactivate_ (const concepts::Connector0 &vtx) |
Deactivates the children of the vertex vtx . More... | |
void | doCorrectRefinement_ (concepts::Hexahedron3d &cell, concepts::Level< 3 > &L) const |
Does apply the correct refinement based on the requested level change. More... | |
Detailed Description
A 3D hp FEM space with continuous, picewise polynomial basis functions.
Currently, only hexahedra are possible.
- Adaptivity
- You can adaptively refine this space by using adjust() on some elements to refine them or increase their polynomial degree. Please note that these extensions can be locally varying (ie. non-uniform) and anisotropic.
- Trunk Space
- Classically, the whole tensor product polynomial space is used on each element to build the degress of freedom. However, a considerable amount of degrees of freedom can be saved, when only a trunk of this space is used. In fact, only very little is lost in accuracy while the gain in terms of degrees of freedom is noteworthy. The methods buildVertexDofs(), buildEdgeDofs(), buildFaceDofs() and buildInnerDofs() can be used to change the way the degrees of freedom on each part of each element are built.
- Timing
- This class is equiped with an interface to get timings of internal computations if compiled accordingly (see
space.cc
file), seesetTimings()
andtimings()
.
- Examples
- hpFEM3d-EV.cc.
Member Typedef Documentation
◆ Scan
typedef concepts::Scan<hp3D::Element<Real> > hp3D::Space::Scan |
◆ Scanner
|
inherited |
◆ TraceSpaceT
◆ type
|
inherited |
Constructor & Destructor Documentation
◆ Space() [1/2]
hp3D::Space::Space | ( | concepts::Mesh3 & | msh, |
uint | l, | ||
uint | p, | ||
concepts::BoundaryConditions * | bc = 0 , |
||
concepts::CellConditions * | cc = 0 |
||
) |
Constructor.
Scans the mesh and sets the cells in the mesh active and the level of refinement and the polynomial degree in all cells to the given values. rebuild_
is set to true, ie. if the mesh is used it will firstly be rebuilt.
- Parameters
-
msh The domain of interest partitioned into a mesh. l Level of refinement p Degree of the polynomials to be used. bc Boundary conditions cc Cell conditions
◆ Space() [2/2]
hp3D::Space::Space | ( | const Space & | spc | ) |
Copy constructor.
◆ ~Space()
|
virtual |
Member Function Documentation
◆ adjust() [1/2]
|
virtual |
◆ adjust() [2/2]
|
pure virtualinherited |
Adjusts the space in the next rebuild step for this element.
◆ bc() [1/2]
|
inline |
◆ bc() [2/2]
|
inlinevirtual |
◆ buildEdgeDofs()
void hp3D::Space::buildEdgeDofs | ( | const BuildDofsBase & | b | ) |
Change the strategy how the degrees of freedom for the edge are built.
The default strategy is BuildEdgeDofs. You can change this strategy any time you chose.
buildEdgeDofs_ is reset with a clone of b
.
- Parameters
-
b New strategy
◆ buildElements_()
|
private |
Creates the elements and their T matrices.
This routine is called recursively until the finest level is reached. If there have degrees of freedom to be assigned on higher level, the according T columns are adapted using the S matrices. On the finest level, the elements are created and the T matrix is created from the T columns.
- Parameters
-
cell The current cell Pmax Request for maximal polynomial degree, array with 3 entries T0 Precomputed T columns which belong to a parent of the cell
- Precondition
Pmax
is an array of 3 integers
◆ buildFaceDofs()
void hp3D::Space::buildFaceDofs | ( | const BuildDofsBase & | b | ) |
Change the strategy how the degrees of freedom for the face are built.
The default strategy is BuildFaceDofs. You can change this strategy any time you chose.
buildFaceDofs_ is reset with a clone of b
.
- Parameters
-
b New strategy
◆ buildInnerDofs()
void hp3D::Space::buildInnerDofs | ( | const BuildDofsBase & | b | ) |
Change the strategy how the degrees of freedom for the interior are built.
The default strategy is BuildInnerDofs. You can change this strategy any time you chose.
buildInnerDofs_ is reset with a clone of b
.
- Parameters
-
b New strategy
◆ buildVertexDofs()
void hp3D::Space::buildVertexDofs | ( | const BuildDofsBase & | b | ) |
Change the strategy how the degrees of freedom for the vertices are built.
The default strategy is BuildVertexDofs. You can change this strategy any time you chose.
buildVertexDofs_ is reset with a clone of b
.
- Parameters
-
b New strategy
◆ computePmax_()
|
private |
Computes the maximal polynomial degree of cell
in each direction.
Takes into account higher polynomial degrees on faces and edges compared to the interior.
◆ createCellList_()
|
private |
Initially fills cellList_
(recusively).
Up to the active level, all cells are stored in the list with their key as list key. The father pointer in the CellData
structure is set during this process.
- Parameters
-
cntr Cell to work on father Parent cell of cell
◆ createEdgFaceList_()
|
private |
Initially fills faceList_
and edgeList_
.
In faceList_
every face of all active cells is stored together with a pointer to its cells. In edgeList_
every edge of all active cells is stored together with pointers to its cells and faces.
- Parameters
-
cell Cell to work on
◆ createVtxEdgList_()
|
private |
Initially fills vertexList_
.
In vertexLIst_
every vertex of all active cells is stored together with pointers to its cells and edges.
- Parameters
-
cell Cell to work on
◆ deactivate_() [1/3]
|
private |
Deactivates the children of the vertex vtx
.
◆ deactivate_() [2/3]
|
private |
Deactivates the children of the edge edg
(including the middle vertex).
◆ deactivate_() [3/3]
|
private |
Deactivates the children of face
(including the new edges and vertices).
◆ dim() [1/2]
|
inlinevirtual |
Returns the dimension of the space.
Implements concepts::SpaceOnCells< Real >.
◆ dim() [2/2]
◆ doCorrectRefinement_()
|
private |
Does apply the correct refinement based on the requested level change.
If the refinement is not readly applicable, upgraded refinements are tested to circumvent a geometric deadlock problem.
This routine is called by meshAndPoly_
.
◆ enforceBC_()
|
private |
Enforce Dirichlet boundary conditions.
This routine is called by meshAndPoly_
on every cell.
- Parameters
-
cell Cell to work on
◆ enforceCC_()
|
private |
Enforce cell condition INACTIVEPLUS.
This routine is called by meshAndPoly_
on every cell in order to deactivate all faces, edges and vertices of cells that are INACTIVEPLUS.
- Parameters
-
cell Cell to work on
◆ enrichElm_()
|
private |
Enriches the polynomial degrees.
This enrichement is necessary for the edges and faces to be able to create conforming global basis functions.
In refined cells, some of the small edges and faces might have a polynomial degree which is too low to represent a part of a basis function which is defined on a larger element.
To achieve this, the cells are traversed top down until the level of the elements is reached. During the traversal, the needed polynomial degree for the small elements is computed and enforced on the element level (only on the edges and faces!).
- Precondition
p
is an array of 27 integers
◆ getOutputDimension()
|
inlinevirtualinherited |
◆ getPmax()
|
inline |
◆ getPmax_()
|
private |
Returns the polynomial degrees of the shape function space on connector.
◆ getSmatrices()
const std::unique_ptr<concepts::SMatrixBase<Real> >* hp3D::Space::getSmatrices | ( | const concepts::Hex3dSubdivision * | subdiv | ) | const |
Returns the list of matrices matching the subdivision type subdiv
.
◆ info()
|
protectedvirtual |
Mesh.
Reimplemented from concepts::SpaceOnCells< Real >.
◆ meshAndPoly_()
|
private |
Adjusts the mesh and the polynomial degree.
On the old mesh, the adjustment information is taken. This leads to mesh refinements or coarsening and a certain polynomial degree distribution.
The polynomial degree is fixed on the finest level (where the elements are eventually built) and given to the coarser levels (for faces, edges and the interior of the cell). This is done bottom up.
Finally, the boundary conditions are enforced on every cell using enforceBC_
@param cell Cell
- Parameters
-
L The desired level P The desired polynomial degree
- Precondition
P
is an array of 3 integers
- See also
- enforceBC_
◆ nelm() [1/2]
|
inlinevirtual |
Returns the number of elements in the space.
Implements concepts::SpaceOnCells< Real >.
◆ nelm() [2/2]
◆ rebuild()
Space hp3D::Space::rebuild | ( | ) |
Rebuilds the space after an adjustment with adjust
.
Steps:
-
Calls
meshAndPoly_
for every cell in the coarse mesh. -
Creates the lists
cellList_
,faceList_
,faceListEdg_
,edgeList_
,edgeListVtx_
, andvertexList_
by callingcreateCellList_
,createEdgFaceList_
andcreateVtxEdgList_
for every cell in the coarse mesh. -
The lists which were created in the previous step are update and changed using
FaceData::checkRelatedFaces
,EdgeData::checkRelatedEdges
,EdgeData::checkRelations
,VertexData::checkRelations
in order to find the supports of the basis functions (which have to be continuous). -
Calls
enrichElm_
for every cell in the coarse mesh. -
Calls
buildElements_
for every cell in the coarse mesh.
They all work recursively and rebuild the space.
At the end of this routine, the information in adj_
, ctrl0_
, ctrl1_
, ctrl2_
, cellList_
, faceList_
, faceListEdg_
, edgeList_
, edgeListVtx_
and vertexList_
is deleted. Just the list of elements in elm_
and ctrl3_
is kept.
- See also
- concepts::TMatrixBase for more information on T matrices and how they are stored
- createCellList_
- createVtxEdgList_
- concepts::FaceData::checkRelatedFaces
- concepts::EdgeData::checkRelatedEdges
- concepts::EdgeData::checkRelations
- concepts::VertexData::checkRelations
- enrichElm_
- buildElements_
- Returns
- Copy of the space on the old level
- Examples
- hpFEM3d-EV.cc.
◆ recomputeShapefunctions()
void hp3D::Space::recomputeShapefunctions | ( | ) |
Recompute shape functions, e.g.
for other abscissas redefined through setIntegrationRule
◆ recomputeSmatrices()
void hp3D::Space::recomputeSmatrices | ( | const Hexahedron * | elm | ) | const |
Checks if the S matrices need to be recomputed and does so if necessary.
- Parameters
-
elm Element which the S matrices are needed for (the polynomial degrees are important)
◆ scan() [1/2]
|
inlinevirtual |
Returns a scanner to iterate over the elements of the space.
Implements concepts::SpaceOnCells< Real >.
◆ scan() [2/2]
|
inlinevirtual |
◆ setDim()
void hp3D::Space::setDim | ( | uint | dim | ) |
◆ setTimings()
|
static |
Sets the class to store the timing values in.
Additionally, the timeCntr_
is reset to 0. This counter is used to fill in the values into the arrays listed below in subsequent calls. The following timings are taken and stored in timings:
-
subroutine
meshAndPoly_
inmeshAndPoly
-
creation of the lists of vertices, edges and faces in
createTopologyLists
-
finding the support of faces modes in
checkRelatedFaces
-
finding the support of edge modes in
checkRelatedEdges
-
finding the support of vertex modes in
checkRelatedVertices
-
enriching the degree of elements in
enrichElements
-
building the elements and the T matrices in
buildElements
◆ timings()
|
static |
Returns true if the class is able to do timings.
The ability to do timings depends on a compiler switch in space.cc
file.
Friends And Related Function Documentation
◆ BuildDofsBase
|
friend |
Member Data Documentation
◆ adj_
|
private |
◆ bc_
|
private |
◆ buildEdgeDofs_
|
private |
◆ buildFaceDofs_
|
private |
◆ buildInnerDofs_
|
private |
◆ buildVertexDofs_
|
private |
◆ cc_
|
private |
◆ cellList_
|
private |
◆ ctrl0_
|
private |
◆ ctrl1_
|
private |
◆ ctrl2_
|
private |
◆ ctrl3_
|
private |
◆ dim_
◆ edgeList_
|
private |
◆ elm_
|
private |
◆ endElm_
|
private |
◆ faceList_
|
private |
◆ msh_
|
private |
◆ nelm_
|
private |
◆ pMax_
|
private |
◆ rebuild_
|
private |
◆ S1left_
|
mutableprivate |
◆ S1right_
|
private |
◆ Smatrices2x_
|
mutableprivate |
The S matrices in 3D for the subdivision into 2 children.
◆ Smatrices2y_
|
private |
◆ Smatrices2z_
|
private |
◆ Smatrices4x_
|
mutableprivate |
The S matrices in 3D for the subdivision into 4 children.
◆ Smatrices4y_
|
private |
◆ Smatrices4z_
|
private |
◆ Smatrices8_
|
mutableprivate |
The S matrices in 3D for the subdivision into 8 children.
- See also
- concepts::Hex3dSubdiv8
◆ timeCntr_
|
staticprivate |
◆ timings_
|
staticprivate |
◆ vertexList_
|
private |
The documentation for this class was generated from the following file:
- hp3D/space.hh