IPPL (Independent Parallel Particle Layer)
IPPL
Loading...
Searching...
No Matches
ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType > Class Template Reference

A class representing a Nedelec space for finite element methods on a structured, rectilinear grid. More...

#include <NedelecSpace.h>

Inheritance diagram for ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >:
Collaboration diagram for ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >:

Public Types

typedef FiniteElementSpace< T, Dim, numElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::indices_t indices_t
typedef FiniteElementSpace< T, Dim, numElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::point_t point_t
typedef FiniteElementSpace< T, Dim, numElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::vertex_points_t vertex_points_t
typedef FieldLayout< DimLayout_t
typedef detail::ViewType< T, 1 >::view_type ViewType
typedef detail::ViewType< T, 1, Kokkos::MemoryTraits< Kokkos::Atomic > >::view_type AtomicViewType
typedef Vector< size_t, numElementVerticesvertex_indices_t
typedef Vector< indices_t, numElementVerticesindices_list_t

Public Member Functions

 NedelecSpace (UniformCartesian< T, Dim > &mesh, ElementType &ref_element, const QuadratureType &quadrature, const Layout_t &layout)
 Construct a new NedelecSpace object.
 NedelecSpace (UniformCartesian< T, Dim > &mesh, ElementType &ref_element, const QuadratureType &quadrature)
 Construct a new NedelecSpace object (without layout) This constructor is made to work with the default constructor in FEMPoissonSolver.h such that it is compatible with alpine.
void initialize (UniformCartesian< T, Dim > &mesh, const Layout_t &layout)
 Initialize a NedelecSpace object created with the default constructor.
void initializeElementIndices (const Layout_t &layout)
 Initialize a Kokkos view containing the element indices.
KOKKOS_FUNCTION size_t numGlobalDOFs () const override
 Degree of Freedom operations //////////////////////////////////////.
KOKKOS_FUNCTION size_t getLocalDOFIndex (const size_t &elementIndex, const size_t &globalDOFIndex) const override
 Get the elements local DOF from the element index and global DOF index.
KOKKOS_FUNCTION size_t getGlobalDOFIndex (const size_t &elementIndex, const size_t &localDOFIndex) const override
 Get the global DOF index from the element index and local DOF.
KOKKOS_FUNCTION Vector< size_t, numElementDOFsgetLocalDOFIndices () const override
 Get the local DOF indices (vector of local DOF indices) They are independent of the specific element because it only depends on the reference element type.
KOKKOS_FUNCTION Vector< size_t, numElementDOFsgetGlobalDOFIndices (const size_t &elementIndex) const override
 Get the global DOF indices (vector of global DOF indices) of an element.
KOKKOS_FUNCTION Vector< size_t, numElementDOFsgetGlobalDOFIndices (const indices_t &elementIndex) const
 Get the global DOF indices (vector of global DOF indices) of an element.
KOKKOS_FUNCTION Vector< size_t, numElementDOFsgetFEMVectorDOFIndices (const size_t &elementIndex, NDIndex< Dim > ldom) const
 Get the DOF indices (vector of indices) corresponding to the position inside the FEMVector of an element.
KOKKOS_FUNCTION Vector< size_t, numElementDOFsgetFEMVectorDOFIndices (indices_t elementIndex, NDIndex< Dim > ldom) const
 Get the DOF indices (vector of indices) corresponding to the position inside the FEMVector of an element.
KOKKOS_FUNCTION point_t getLocalDOFPosition (size_t localDOFIndex) const
 Get the cartesion position of a local DOF in the reference element.
KOKKOS_FUNCTION point_t evaluateRefElementShapeFunction (const size_t &localDOF, const point_t &localPoint) const
 Basis functions and gradients /////////////////////////////////////.
KOKKOS_FUNCTION point_t evaluateRefElementShapeFunctionCurl (const size_t &localDOF, const point_t &localPoint) const
 Evaluate the curl of the shape function of a local degree of of freedom at ta given point in the reference element.
template<typename F>
FEMVector< TevaluateAx (FEMVector< T > &x, F &evalFunction) const
 Assembly operations ///////////////////////////////////////////////.
FEMVector< TevaluateLoadVector (const FEMVector< point_t > &f) const
 Assemble the load vector b of the system Ax = b, given a field of the right hand side defined at the Nédélec DOF positions. If a functor instead of a field should be used, use the function NedelecSpace::evaluateLoadVectorFunctor.
template<typename F>
FEMVector< TevaluateLoadVectorFunctor (const F &f) const
 Assemble the load vector b of the system Ax = b, given a functional of the rhs. If a field instead of a functor should be used, use the function NedelecSpace::evaluateLoadVector.
FEMVector< TcreateFEMVector () const
 FEMVector conversion and creation//////////////////////////////////.
Kokkos::View< point_t * > reconstructToPoints (const Kokkos::View< point_t * > &positions, const FEMVector< T > &coef) const
 Reconstructs function values at arbitrary points in the mesh given the Nedelec DOF coefficients.
template<typename F>
T computeError (const FEMVector< T > &u_h, const F &u_sol) const
 Error norm computations ///////////////////////////////////////////.
KOKKOS_FUNCTION bool isDOFOnBoundary (const size_t &dofIdx) const
 Check if a DOF is on the boundary of the mesh.
KOKKOS_FUNCTION int getBoundarySide (const size_t &dofIdx) const
 Returns which side the boundary is on.
void setMesh (UniformCartesian< T, Dim > &mesh)
KOKKOS_FUNCTION size_t numElements () const
 Mesh and Element operations ///////////////////////////////////////.
KOKKOS_FUNCTION size_t numElementsInDim (const size_t &dim) const
 Get the number of elements in a given dimension.
KOKKOS_FUNCTION indices_t getMeshVertexNDIndex (const size_t &vertex_index) const
 Get the NDIndex of a mesh vertex.
KOKKOS_FUNCTION size_t getMeshVertexIndex (const indices_t &vertex_nd_index) const
 Get the global index of a mesh vertex given its NDIndex.
KOKKOS_FUNCTION indices_t getElementNDIndex (const size_t &elementIndex) const
 Get the NDIndex (vector of indices for each dimension) of a mesh element.
KOKKOS_FUNCTION size_t getElementIndex (const indices_t &ndindex) const
 Get the global index of a mesh element given the NDIndex.
KOKKOS_FUNCTION vertex_indices_t getElementMeshVertexIndices (const indices_t &elementNDIndex) const
 Get all the global vertex indices of an element (given by its NDIndex).
KOKKOS_FUNCTION indices_list_t getElementMeshVertexNDIndices (const indices_t &elementNDIndex) const
 Get all the NDIndices of the vertices of an element (given by its NDIndex).
KOKKOS_FUNCTION vertex_points_t getElementMeshVertexPoints (const indices_t &elementNDIndex) const
 Get all the global vertex points of an element (given by its NDIndex).

Public Attributes

UniformCartesian< T, Dim > & mesh_m
 Member variables //////////////////////////////////////////////////.
ElementType ref_element_m
const QuadratureType & quadrature_m
Vector< size_t, Dimnr_m
Vector< double, Dimhr_m
Vector< double, Dimorigin_m

Static Public Attributes

static constexpr unsigned numElementDOFs = getNedelecNumElementDOFs(Dim, Order)
static constexpr unsigned dim
static constexpr unsigned order = Order
static constexpr unsigned numElementVertices

Private Member Functions

FEMVector< TcreateFEMVector2d () const
 Implementation of the NedelecSpace::createFEMVector function for 2d.
FEMVector< TcreateFEMVector3d () const
 Implementation of the NedelecSpace::createFEMVector function for 3d.

Private Attributes

Kokkos::View< size_t * > elementIndices
 Stores which elements (squares or cubes) belong to the current MPI rank.
Vector< point_t, 12 > localDofPositions_m
 Stores the positions of the local Degrees of Freedoms on the reference elements.
Layout_t layout_m
 The layout of the MPI ranks over the mesh.

Detailed Description

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
class ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >

A class representing a Nedelec space for finite element methods on a structured, rectilinear grid.

Template Parameters
TThe floating point number type of the field values
DimThe dimension of the mesh
OrderThe order of the Nedelec space
QuadratureTypeThe type of the quadrature rule
FieldTypeThe type of field to use.

Definition at line 34 of file NedelecSpace.h.

Member Typedef Documentation

◆ AtomicViewType

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
typedef detail::ViewType<T,1,Kokkos::MemoryTraits<Kokkos::Atomic>>::view_type ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::AtomicViewType

Definition at line 68 of file NedelecSpace.h.

◆ indices_list_t

typedef Vector<indices_t, numElementVertices> ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::indices_list_t
inherited

Definition at line 59 of file FiniteElementSpace.h.

◆ indices_t

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
typedef FiniteElementSpace<T,Dim,numElementDOFs,ElementType,QuadratureType,FEMVector<T>,FEMVector<T>>::indices_t ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::indices_t

Definition at line 53 of file NedelecSpace.h.

◆ Layout_t

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
typedef FieldLayout<Dim> ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::Layout_t

Definition at line 63 of file NedelecSpace.h.

◆ point_t

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
typedef FiniteElementSpace<T,Dim,numElementDOFs,ElementType,QuadratureType,FEMVector<T>,FEMVector<T>>::point_t ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::point_t

Definition at line 57 of file NedelecSpace.h.

◆ vertex_indices_t

typedef Vector<size_t, numElementVertices> ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::vertex_indices_t
inherited

Definition at line 57 of file FiniteElementSpace.h.

◆ vertex_points_t

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
typedef FiniteElementSpace<T,Dim,numElementDOFs,ElementType,QuadratureType,FEMVector<T>,FEMVector<T>>::vertex_points_t ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::vertex_points_t

Definition at line 60 of file NedelecSpace.h.

◆ ViewType

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
typedef detail::ViewType<T,1>::view_type ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::ViewType

Definition at line 66 of file NedelecSpace.h.

Constructor & Destructor Documentation

◆ NedelecSpace() [1/2]

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::NedelecSpace ( UniformCartesian< T, Dim > & mesh,
ElementType & ref_element,
const QuadratureType & quadrature,
const Layout_t & layout )

Construct a new NedelecSpace object.

Parameters
meshReference to the mesh
ref_elementReference to the reference element
quadratureReference to the quadrature rule
layoutReference to the field layout

Definition at line 9 of file NedelecSpace.hpp.

References Dim, ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::FiniteElementSpace(), getNedelecNumElementDOFs(), and initializeElementIndices().

Referenced by evaluateLoadVector(), evaluateRefElementShapeFunction(), evaluateRefElementShapeFunctionCurl(), and reconstructToPoints().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ NedelecSpace() [2/2]

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::NedelecSpace ( UniformCartesian< T, Dim > & mesh,
ElementType & ref_element,
const QuadratureType & quadrature )

Construct a new NedelecSpace object (without layout) This constructor is made to work with the default constructor in FEMPoissonSolver.h such that it is compatible with alpine.

Parameters
meshReference to the mesh
ref_elementReference to the reference element
quadratureReference to the quadrature rule

Definition at line 26 of file NedelecSpace.hpp.

References Dim, ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::FiniteElementSpace(), and getNedelecNumElementDOFs().

Here is the call graph for this function:

Member Function Documentation

◆ computeError()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
template<typename F>
T ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::computeError ( const FEMVector< T > & u_h,
const F & u_sol ) const

Error norm computations ///////////////////////////////////////////.

Given the Nedelec space DoF coefficients and an analytical solution computes the L2 norm error.

Parameters
u_hThe basis function coefficients obtained via FEM.  
u_solThe analytical solution (functor)
Returns
error - The error ||u_h - u_sol||_L2

Definition at line 901 of file NedelecSpace.hpp.

References ippl::Comm, ippl::Vector< T, Dim >::dot(), elementIndices, evaluateRefElementShapeFunction(), ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::getElementMeshVertexPoints(), getFEMVectorDOFIndices(), getGlobalDOFIndices(), ippl::FEMVector< T >::getView(), layout_m, numElementDOFs, ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::quadrature_m, and ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::ref_element_m.

Here is the call graph for this function:

◆ createFEMVector()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
FEMVector< T > ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::createFEMVector ( ) const

FEMVector conversion and creation//////////////////////////////////.

FEMVector conversion //////////////////////////////////////////////.

Creates and empty FEMVector.

Creates and empty FEMVector which corresponds to the domain this MPI rank owns (according to the ippl layout created for this mesh). To this extend it will also setup all the information needed to exchange halo cells.

Returns
An empty FEMVector for this domain.

Definition at line 804 of file NedelecSpace.hpp.

References createFEMVector2d(), createFEMVector3d(), and Dim.

Referenced by evaluateLoadVector(), and evaluateLoadVectorFunctor().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ createFEMVector2d()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
FEMVector< T > ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::createFEMVector2d ( ) const
private

Implementation of the NedelecSpace::createFEMVector function for 2d.

Definition at line 1169 of file NedelecSpace.hpp.

References ippl::Comm, ippl::NDIndex< Dim >::first(), getFEMVectorDOFIndices(), ippl::NDIndex< Dim >::last(), and layout_m.

Referenced by createFEMVector().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ createFEMVector3d()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
FEMVector< T > ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::createFEMVector3d ( ) const
private

Implementation of the NedelecSpace::createFEMVector function for 3d.

Definition at line 1417 of file NedelecSpace.hpp.

References ippl::Comm, ippl::NDIndex< Dim >::first(), getFEMVectorDOFIndices(), ippl::NDIndex< Dim >::last(), and layout_m.

Referenced by createFEMVector().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ evaluateAx()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
template<typename F>
FEMVector< T > ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateAx ( FEMVector< T > & x,
F & evalFunction ) const

Assembly operations ///////////////////////////////////////////////.

Assemble the left stiffness matrix A of the system Ax = b

Parameters
xThe vector which we want to multiply
Returns
The vector containing A*x

Definition at line 364 of file NedelecSpace.hpp.

References elementIndices, evaluateRefElementShapeFunction(), evaluateRefElementShapeFunctionCurl(), getFEMVectorDOFIndices(), getGlobalDOFIndices(), IpplTimings::getTimer(), ippl::FEMVector< T >::getView(), isDOFOnBoundary(), layout_m, numElementDOFs, ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::quadrature_m, IpplTimings::startTimer(), and IpplTimings::stopTimer().

Here is the call graph for this function:

◆ evaluateLoadVector()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
FEMVector< T > ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateLoadVector ( const FEMVector< point_t > & f) const

Assemble the load vector b of the system Ax = b, given a field of the right hand side defined at the Nédélec DOF positions. If a functor instead of a field should be used, use the function NedelecSpace::evaluateLoadVectorFunctor.

Parameters
fThe source field defined at the Nédélec degrees fo freedom.
Returns
The resulting rhs b of the Galerkin discretization.

Definition at line 483 of file NedelecSpace.hpp.

References createFEMVector(), Dim, ippl::Vector< T, Dim >::dot(), elementIndices, evaluateRefElementShapeFunction(), getFEMVectorDOFIndices(), getGlobalDOFIndices(), getLocalDOFPosition(), ippl::FEMVector< T >::getView(), isDOFOnBoundary(), layout_m, NedelecSpace(), numElementDOFs, ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::quadrature_m, and ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::ref_element_m.

Here is the call graph for this function:

◆ evaluateLoadVectorFunctor()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
template<typename F>
FEMVector< T > ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateLoadVectorFunctor ( const F & f) const

◆ evaluateRefElementShapeFunction()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
KOKKOS_FUNCTION NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::point_t ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateRefElementShapeFunction ( const size_t & localDOF,
const point_t & localPoint ) const

Basis functions and gradients /////////////////////////////////////.

Evaluate the shape function of a local degree of freedom at a given point in the reference element

Parameters
localDOFsize_t - The local degree of freedom index
localPointpoint_t (Vector<T, Dim>) - The point in the reference element
Returns
T - The value of the shape function at the given point

Definition at line 695 of file NedelecSpace.hpp.

References Dim, NedelecSpace(), numElementDOFs, and ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::ref_element_m.

Referenced by computeError(), evaluateAx(), evaluateLoadVector(), evaluateLoadVectorFunctor(), and reconstructToPoints().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ evaluateRefElementShapeFunctionCurl()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
KOKKOS_FUNCTION NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::point_t ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateRefElementShapeFunctionCurl ( const size_t & localDOF,
const point_t & localPoint ) const

Evaluate the curl of the shape function of a local degree of of freedom at ta given point in the reference element.

Parameters
localDOFsize_t - The local degree of freedom index
localPointpoint_t (Vector<T, Dim>) - The point in the reference element
Returns
point_t (Vector<T, Dim>) - The curl of the shape function at the given point

Definition at line 750 of file NedelecSpace.hpp.

References Dim, and NedelecSpace().

Referenced by evaluateAx().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ getBoundarySide()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
KOKKOS_FUNCTION int ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::getBoundarySide ( const size_t & dofIdx) const

Returns which side the boundary is on.

This function takes as input the global index of a DoF and then returns on which side of the domain boundary it is on, in 2d that would be either south, north, west, east and in 3d space and ground is added. The mapping is as follows: 0 = south 1 = west 2 = north 3 = east 4 = ground 5 = space -1 = not on a boundary.

Parameters
dofIdxthe global DoF index for which the boundary side should be retrieved.
Returns
Which boundary side the DoF is on or -1 if on no boundary.

Definition at line 1080 of file NedelecSpace.hpp.

References Dim, and ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::nr_m.

◆ getElementIndex()

KOKKOS_FUNCTION size_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::getElementIndex ( const indices_t & ndindex) const
inherited

Get the global index of a mesh element given the NDIndex.

Parameters
ndindexindices_t (Vector<size_t, Dim>) - vector of indices for each direction
Returns
size_t - the index of the element

Definition at line 135 of file FiniteElementSpace.hpp.

Referenced by ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::initializeElementIndices().

Here is the caller graph for this function:

◆ getElementMeshVertexIndices()

KOKKOS_FUNCTION FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::vertex_indices_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::getElementMeshVertexIndices ( const indices_t & elementNDIndex) const
inherited

Get all the global vertex indices of an element (given by its NDIndex).

Parameters
elementNDIndexThe NDIndex of the element
Returns
vertex_indices_t (Vector<size_t, numElementVertices>) - vector of vertex indices

Definition at line 146 of file FiniteElementSpace.hpp.

◆ getElementMeshVertexNDIndices()

KOKKOS_FUNCTION FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::indices_list_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::getElementMeshVertexNDIndices ( const indices_t & elementNDIndex) const
inherited

Get all the NDIndices of the vertices of an element (given by its NDIndex).

Parameters
elementNDIndexThe NDIndex of the element
Returns
indices_list_t (Vector<Vector<size_t, Dim>, numElementVertices>) - vector of vertex NDIndices

Definition at line 157 of file FiniteElementSpace.hpp.

◆ getElementMeshVertexPoints()

KOKKOS_FUNCTION FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::vertex_points_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::getElementMeshVertexPoints ( const indices_t & elementNDIndex) const
inherited

Get all the global vertex points of an element (given by its NDIndex).

Parameters
elementNDIndexThe NDIndex of the element
Returns
vertex_points_t (Vector<Vector<T, Dim>, numElementVertices>) -

Definition at line 167 of file FiniteElementSpace.hpp.

Referenced by ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::computeError(), and ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateLoadVectorFunctor().

Here is the caller graph for this function:

◆ getElementNDIndex()

KOKKOS_FUNCTION FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::indices_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::getElementNDIndex ( const size_t & elementIndex) const
inherited

Get the NDIndex (vector of indices for each dimension) of a mesh element.

Parameters
elementIndexindices_t (Vector<size_t, Dim>) - The index of the element
Returns
indices_t (Vector<size_t, Dim>) - vector of indices for each dimension

Definition at line 126 of file FiniteElementSpace.hpp.

Referenced by ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateLoadVectorFunctor(), ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::getFEMVectorDOFIndices(), and ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::getGlobalDOFIndices().

Here is the caller graph for this function:

◆ getFEMVectorDOFIndices() [1/2]

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
KOKKOS_FUNCTION Vector< size_t, NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::numElementDOFs > ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::getFEMVectorDOFIndices ( const size_t & elementIndex,
NDIndex< Dim > ldom ) const

Get the DOF indices (vector of indices) corresponding to the position inside the FEMVector of an element.

Parameters
elementIndexsize_t - The index of the element
Returns
Vector<size_t, NumElementDOFs> - The DOF indices

Definition at line 335 of file NedelecSpace.hpp.

References ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::getElementNDIndex(), and getFEMVectorDOFIndices().

Referenced by computeError(), createFEMVector2d(), createFEMVector3d(), evaluateAx(), evaluateLoadVector(), evaluateLoadVectorFunctor(), getFEMVectorDOFIndices(), and reconstructToPoints().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ getFEMVectorDOFIndices() [2/2]

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
KOKKOS_FUNCTION Vector< size_t, numElementDOFs > ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::getFEMVectorDOFIndices ( indices_t elementIndex,
NDIndex< Dim > ldom ) const

Get the DOF indices (vector of indices) corresponding to the position inside the FEMVector of an element.

Parameters
elementIndexindices_t - The index of the element
Returns
Vector<size_t, NumElementDOFs> - The DOF indices

◆ getGlobalDOFIndex()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
KOKKOS_FUNCTION size_t ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::getGlobalDOFIndex ( const size_t & elementIndex,
const size_t & localDOFIndex ) const
overridevirtual

Get the global DOF index from the element index and local DOF.

Parameters
elementIndexsize_t - The index of the element
localDOFIndexsize_t - The local DOF index
Returns
size_t - The global DOF index

Implements ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >.

Definition at line 180 of file NedelecSpace.hpp.

References getGlobalDOFIndices().

Here is the call graph for this function:

◆ getGlobalDOFIndices() [1/2]

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
KOKKOS_FUNCTION Vector< size_t, numElementDOFs > ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::getGlobalDOFIndices ( const indices_t & elementIndex) const

Get the global DOF indices (vector of global DOF indices) of an element.

Parameters
elementIndexindices_t - The index of the element
Returns
Vector<size_t, NumElementDOFs> - The global DOF indices

◆ getGlobalDOFIndices() [2/2]

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
KOKKOS_FUNCTION Vector< size_t, NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::numElementDOFs > ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::getGlobalDOFIndices ( const size_t & elementIndex) const
overridevirtual

Get the global DOF indices (vector of global DOF indices) of an element.

Parameters
elementIndexsize_t - The index of the element
Returns
Vector<size_t, NumElementDOFs> - The global DOF indices

Implements ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >.

Definition at line 259 of file NedelecSpace.hpp.

References ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::getElementNDIndex(), and getGlobalDOFIndices().

Referenced by computeError(), evaluateAx(), evaluateLoadVector(), evaluateLoadVectorFunctor(), getGlobalDOFIndex(), getGlobalDOFIndices(), and getLocalDOFIndex().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ getLocalDOFIndex()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
KOKKOS_FUNCTION size_t ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::getLocalDOFIndex ( const size_t & elementIndex,
const size_t & globalDOFIndex ) const
overridevirtual

Get the elements local DOF from the element index and global DOF index.

Parameters
elementIndexsize_t - The index of the element
globalDOFIndexsize_t - The global DOF index
Returns
size_t - The local DOF index

Implements ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >.

Definition at line 149 of file NedelecSpace.hpp.

References Dim, ippl::Vector< T, Dim >::dim, and getGlobalDOFIndices().

Here is the call graph for this function:

◆ getLocalDOFIndices()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
KOKKOS_FUNCTION Vector< size_t, NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::numElementDOFs > ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::getLocalDOFIndices ( ) const
overridevirtual

Get the local DOF indices (vector of local DOF indices) They are independent of the specific element because it only depends on the reference element type.

Returns
Vector<size_t, NumElementDOFs> - The local DOF indices

Implements ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >.

Definition at line 193 of file NedelecSpace.hpp.

References numElementDOFs.

◆ getLocalDOFPosition()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
KOKKOS_FUNCTION NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::point_t ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::getLocalDOFPosition ( size_t localDOFIndex) const

Get the cartesion position of a local DOF in the reference element.

Given the local DOF index this function will return the cartesian position of this DOF with respect to the reference element.

Definition at line 347 of file NedelecSpace.hpp.

References localDofPositions_m.

Referenced by evaluateLoadVector().

Here is the caller graph for this function:

◆ getMeshVertexIndex()

KOKKOS_FUNCTION size_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::getMeshVertexIndex ( const indices_t & vertex_nd_index) const
inherited

Get the global index of a mesh vertex given its NDIndex.

Parameters
vertex_nd_indexindices_t (Vector<size_t, Dim>) - The NDIndex of the vertex (vector of indices for each dimension).
Returns
size_t - unsigned integer index of the mesh vertex

Definition at line 117 of file FiniteElementSpace.hpp.

◆ getMeshVertexNDIndex()

KOKKOS_FUNCTION FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::indices_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::getMeshVertexNDIndex ( const size_t & vertex_index) const
inherited

Get the NDIndex of a mesh vertex.

Parameters
vertex_indexsize_t - The index of the vertex
Returns
indices_t (Vector<size_t, Dim>) - Returns the NDIndex (vector of indices for each dimension)

Definition at line 107 of file FiniteElementSpace.hpp.

◆ initialize()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
void ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::initialize ( UniformCartesian< T, Dim > & mesh,
const Layout_t & layout )

Initialize a NedelecSpace object created with the default constructor.

Parameters
meshReference to the mesh
layoutReference to the field layout

Definition at line 43 of file NedelecSpace.hpp.

References Dim, ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::FiniteElementSpace(), getNedelecNumElementDOFs(), initializeElementIndices(), localDofPositions_m, and ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::setMesh().

Here is the call graph for this function:

◆ initializeElementIndices()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
void ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::initializeElementIndices ( const Layout_t & layout)

Initialize a Kokkos view containing the element indices.

Definition at line 73 of file NedelecSpace.hpp.

References Dim, elementIndices, first(), ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::getElementIndex(), ippl::FieldLayout< Dim >::getLocalNDIndex(), layout_m, ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::nr_m, and ippl::NDIndex< Dim >::size().

Referenced by initialize(), and NedelecSpace().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ isDOFOnBoundary()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
KOKKOS_FUNCTION bool ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::isDOFOnBoundary ( const size_t & dofIdx) const

Check if a DOF is on the boundary of the mesh.

This function takes as input the global index of a DoF and returns if this DoF is on. If one would like to know which boundary this is the function NedelecSpace::getBoundarySide can be used.

Parameters
dofIdxThe global DoF index for which should be checked if it is on the boundary.
Returns
true - If the DOF is on the domain boundary
false - If the DOF is not on the domain boundary

Definition at line 992 of file NedelecSpace.hpp.

References Dim, and ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::nr_m.

Referenced by evaluateAx(), evaluateLoadVector(), and evaluateLoadVectorFunctor().

Here is the caller graph for this function:

◆ numElements()

KOKKOS_FUNCTION size_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::numElements ( ) const
inherited

Mesh and Element operations ///////////////////////////////////////.

Get the number of elements in the mesh of the space

Returns
size_t - unsigned integer number of elements

Definition at line 88 of file FiniteElementSpace.hpp.

◆ numElementsInDim()

KOKKOS_FUNCTION size_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::numElementsInDim ( const size_t & dim) const
inherited

Get the number of elements in a given dimension.

Parameters
dimsize_t - representing the dimension
Returns
size_t - unsigned integer number of elements in the given dimension

Definition at line 97 of file FiniteElementSpace.hpp.

◆ numGlobalDOFs()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
KOKKOS_FUNCTION size_t ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::numGlobalDOFs ( ) const
overridevirtual

Degree of Freedom operations //////////////////////////////////////.

Get the number of global degrees of freedom in the space

Returns
size_t - unsigned integer number of global degrees of freedom

Implements ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >.

Definition at line 130 of file NedelecSpace.hpp.

References Dim, and ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::nr_m.

◆ reconstructToPoints()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
Kokkos::View< typename NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::point_t * > ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::reconstructToPoints ( const Kokkos::View< point_t * > & positions,
const FEMVector< T > & coef ) const

Reconstructs function values at arbitrary points in the mesh given the Nedelec DOF coefficients.

This function can be used to retrieve the values of a solution function at arbitrary points inside of the mesh given the Nedelec DOF coefficients which solved the problem using FEM.

Note
Currently the function is able to handle both cases, where we have that positions only contains positions which are inside of local domain of this MPI rank (i.e. each rank gets its own unique position ) and where positions contains the positions of all ranks (i.e. positions is the same for all ranks). If in the future it can be guaranteed, that each rank will get its own positions then certain parts of the function implementation can be removed. Instructions for this are given in the implementation itself.
Parameters
positionsThe points at which the function should be evaluated. A Kokkos::View which stores in each element a 2D/3D point.
coefThe basis function coefficients obtained via FEM.
Returns
The function evaluated at the given points, stored inside of Kokkos::View where each element corresponts to the function value at the point described by the same element inside of positions.

Definition at line 819 of file NedelecSpace.hpp.

References Dim, evaluateRefElementShapeFunction(), getFEMVectorDOFIndices(), ippl::FEMVector< T >::getView(), ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::hr_m, ippl::NDIndex< Dim >::last(), layout_m, NedelecSpace(), ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::nr_m, numElementDOFs, and ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::origin_m.

Here is the call graph for this function:

◆ setMesh()

void ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::setMesh ( UniformCartesian< T, Dim > & mesh)
inherited

Definition at line 77 of file FiniteElementSpace.hpp.

Referenced by ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::initialize().

Here is the caller graph for this function:

Member Data Documentation

◆ dim

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
unsigned ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::dim
staticconstexpr
Initial value:

Definition at line 41 of file NedelecSpace.h.

◆ elementIndices

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
Kokkos::View<size_t*> ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::elementIndices
private

Stores which elements (squares or cubes) belong to the current MPI rank.

Definition at line 407 of file NedelecSpace.h.

Referenced by computeError(), evaluateAx(), evaluateLoadVector(), evaluateLoadVectorFunctor(), and initializeElementIndices().

◆ hr_m

Vector<double, Dim> ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::hr_m
inherited

◆ layout_m

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
Layout_t ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::layout_m
private

The layout of the MPI ranks over the mesh.

Standart ippl layout which dictates how the MPI ranks are layed out over the mesh. It is used in order to be able to create FEMVectors, retreive correct DOF indices and intitalize the elementIndices.

Definition at line 425 of file NedelecSpace.h.

Referenced by computeError(), createFEMVector2d(), createFEMVector3d(), evaluateAx(), evaluateLoadVector(), evaluateLoadVectorFunctor(), initializeElementIndices(), and reconstructToPoints().

◆ localDofPositions_m

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
Vector<point_t, 12> ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::localDofPositions_m
private

Stores the positions of the local Degrees of Freedoms on the reference elements.

We are saying that the local degree of freedom positions are simply the centers of the edges.

Definition at line 416 of file NedelecSpace.h.

Referenced by getLocalDOFPosition(), and initialize().

◆ mesh_m

UniformCartesian<T, Dim>& ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::mesh_m
inherited

Member variables //////////////////////////////////////////////////.

Definition at line 227 of file FiniteElementSpace.h.

◆ nr_m

◆ numElementDOFs

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
unsigned ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::numElementDOFs = getNedelecNumElementDOFs(Dim, Order)
staticconstexpr

◆ numElementVertices

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
unsigned ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::numElementVertices
staticconstexpr
Initial value:

Definition at line 48 of file NedelecSpace.h.

◆ order

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldType>
unsigned ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::order = Order
staticconstexpr

Definition at line 45 of file NedelecSpace.h.

◆ origin_m

Vector<double, Dim> ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::origin_m
inherited

◆ quadrature_m

◆ ref_element_m


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