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

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

#include <LagrangeSpace.h>

Inheritance diagram for ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >:
Collaboration diagram for ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >:

Classes

struct  DeviceStruct
 Device struct for copies //////////////////////////////////////////. More...

Public Types

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

Public Member Functions

 LagrangeSpace (UniformCartesian< T, Dim > &mesh, ElementType &ref_element, const QuadratureType &quadrature, const Layout_t &layout)
 Construct a new LagrangeSpace object.
 LagrangeSpace (UniformCartesian< T, Dim > &mesh, ElementType &ref_element, const QuadratureType &quadrature)
 Construct a new LagrangeSpace 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 LagrangeSpace object created with the default constructor.
void initializeElementIndices (const Layout_t &layout)
 Initialize a Kokkos view containing the element indices. This distributes the elements among MPI ranks.
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 &element_index) const override
 Get the global DOF indices (vector of global DOF indices) of an element.
KOKKOS_FUNCTION T evaluateRefElementShapeFunction (const size_t &localDOF, const point_t &localPoint) const
 Basis functions and gradients /////////////////////////////////////.
KOKKOS_FUNCTION point_t evaluateRefElementShapeFunctionGradient (const size_t &localDOF, const point_t &localPoint) const
 Evaluate the gradient of the shape function of a local degree of freedom at a given point in the reference element.
template<typename F>
FieldLHS evaluateAx (FieldLHS &field, F &evalFunction) const
 Assembly operations ///////////////////////////////////////////////.
template<typename F>
FieldLHS evaluateAx_lower (FieldLHS &field, F &evalFunction) const
template<typename F>
FieldLHS evaluateAx_upper (FieldLHS &field, F &evalFunction) const
template<typename F>
FieldLHS evaluateAx_upperlower (FieldLHS &field, F &evalFunction) const
template<typename F>
FieldLHS evaluateAx_inversediag (FieldLHS &field, F &evalFunction) const
template<typename F>
FieldLHS evaluateAx_diag (FieldLHS &field, F &evalFunction) const
template<typename F>
FieldLHS evaluateAx_lift (FieldLHS &field, F &evalFunction) const
 Assemble the left stiffness matrix A of the system but only for the boundary values, so that they can be subtracted from the RHS for treatment of Dirichlet BCs.
void evaluateLoadVector (FieldRHS &field) const
 Assemble the load vector b of the system Ax = b.
template<typename F>
T computeErrorL2 (const FieldLHS &u_h, const F &u_sol) const
 Error norm computations ///////////////////////////////////////////.
T computeAvg (const FieldLHS &u_h) const
 Given a field, compute the average.
DeviceStruct getDeviceMirror () const
 Device struct definitions /////////////////////////////////////////.
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 = getLagrangeNumElementDOFs(Dim, Order)
static constexpr unsigned dim
static constexpr unsigned order = Order
static constexpr unsigned numElementVertices

Private Member Functions

KOKKOS_FUNCTION bool isDOFOnBoundary (const indices_t &ndindex) const
 Check if a DOF is on the boundary of the mesh.

Private Attributes

Kokkos::View< size_t * > elementIndices

Detailed Description

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
class ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >

A class representing a Lagrange 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 Lagrange space
QuadratureTypeThe type of the quadrature rule
FieldLHSThe type of the left hand side field
FieldRHSThe type of the right hand side field

Definition at line 33 of file LagrangeSpace.h.

Member Typedef Documentation

◆ AtomicViewType

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
typedef detail::ViewType<T,Dim,Kokkos::MemoryTraits<Kokkos::Atomic>>::view_type ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::AtomicViewType

Definition at line 69 of file LagrangeSpace.h.

◆ indices_list_t

typedef Vector<indices_t, numElementVertices> ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FieldLHS, FieldRHS >::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 FieldLHS, typename FieldRHS>
typedef FiniteElementSpace<T,Dim,numElementDOFs,ElementType,QuadratureType,FieldLHS,FieldRHS>::indices_t ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::indices_t

Definition at line 54 of file LagrangeSpace.h.

◆ Layout_t

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
typedef FieldLayout<Dim> ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::Layout_t

Definition at line 64 of file LagrangeSpace.h.

◆ point_t

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
typedef FiniteElementSpace<T,Dim,numElementDOFs,ElementType,QuadratureType,FieldLHS,FieldRHS>::point_t ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::point_t

Definition at line 58 of file LagrangeSpace.h.

◆ vertex_indices_t

typedef Vector<size_t, numElementVertices> ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FieldLHS, FieldRHS >::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 FieldLHS, typename FieldRHS>
typedef FiniteElementSpace<T,Dim,numElementDOFs,ElementType,QuadratureType,FieldLHS,FieldRHS>::vertex_points_t ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::vertex_points_t

Definition at line 61 of file LagrangeSpace.h.

◆ ViewType

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
typedef detail::ViewType<T,Dim>::view_type ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::ViewType

Definition at line 67 of file LagrangeSpace.h.

Constructor & Destructor Documentation

◆ LagrangeSpace() [1/2]

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

Construct a new LagrangeSpace 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 8 of file LagrangeSpace.hpp.

References Dim, ippl::FiniteElementSpace< T, Dim, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >::FiniteElementSpace(), getLagrangeNumElementDOFs(), and initializeElementIndices().

Referenced by ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::DeviceStruct::evaluateRefElementShapeFunction().

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

◆ LagrangeSpace() [2/2]

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

Construct a new LagrangeSpace 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 24 of file LagrangeSpace.hpp.

References Dim, ippl::FiniteElementSpace< T, Dim, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >::FiniteElementSpace(), and getLagrangeNumElementDOFs().

Here is the call graph for this function:

Member Function Documentation

◆ computeAvg()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
T ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::computeAvg ( const FieldLHS & u_h) const

◆ computeErrorL2()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
template<typename F>
T ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::computeErrorL2 ( const FieldLHS & u_h,
const F & u_sol ) const

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

Functions for error computations, etc. ////////////////////////////.

Given two fields, compute the L2 norm error

Parameters
u_hThe numerical solution found using FEM  
u_solThe analytical solution (functor)
Returns
error - The error ||u_h - u_sol||_L2

Definition at line 1428 of file LagrangeSpace.hpp.

References ippl::apply(), ippl::Comm, Dim, elementIndices, evaluateRefElementShapeFunction(), getGlobalDOFIndices(), ippl::FiniteElementSpace< T, Dim, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >::getMeshVertexNDIndex(), numElementDOFs, ippl::FiniteElementSpace< T, Dim, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >::quadrature_m, and ippl::FiniteElementSpace< T, Dim, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >::ref_element_m.

Here is the call graph for this function:

◆ evaluateAx()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
template<typename F>
FieldLHS ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::evaluateAx ( FieldLHS & field,
F & evalFunction ) const

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

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

Parameters
fieldThe field to assemble the matrix for
evalFunctionThe lambda telling us the form which A takes
Returns
FieldLHS - The LHS field containing A*x

Definition at line 350 of file LagrangeSpace.hpp.

References ippl::apply(), ippl::BConds< Field, Dim >::apply(), ippl::BConds< Field, Dim >::assignGhostToPhysical(), ippl::CONSTANT_FACE, Dim, elementIndices, evaluateAx(), evaluateRefElementShapeFunctionGradient(), getGlobalDOFIndices(), ippl::FiniteElementSpace< T, Dim, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >::getMeshVertexNDIndex(), IpplTimings::getTimer(), isDOFOnBoundary(), numElementDOFs, ippl::PERIODIC_FACE, ippl::FiniteElementSpace< T, Dim, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >::quadrature_m, IpplTimings::startTimer(), IpplTimings::stopTimer(), and ippl::ZERO_FACE.

Referenced by evaluateAx().

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

◆ evaluateAx_diag()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
template<typename F>
FieldLHS ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::evaluateAx_diag ( FieldLHS & field,
F & evalFunction ) const

◆ evaluateAx_inversediag()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
template<typename F>
FieldLHS ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::evaluateAx_inversediag ( FieldLHS & field,
F & evalFunction ) const

◆ evaluateAx_lift()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
template<typename F>
FieldLHS ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::evaluateAx_lift ( FieldLHS & field,
F & evalFunction ) const

Assemble the left stiffness matrix A of the system but only for the boundary values, so that they can be subtracted from the RHS for treatment of Dirichlet BCs.

Parameters
fieldThe field to assemble the matrix for
evalFunctionThe lambda telling us the form which A takes
Returns
FieldLHS - The LHS field containing A*x

Definition at line 1191 of file LagrangeSpace.hpp.

References ippl::apply(), Dim, elementIndices, evaluateAx_lift(), evaluateRefElementShapeFunctionGradient(), getGlobalDOFIndices(), ippl::FiniteElementSpace< T, Dim, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >::getMeshVertexNDIndex(), isDOFOnBoundary(), numElementDOFs, and ippl::FiniteElementSpace< T, Dim, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >::quadrature_m.

Referenced by evaluateAx_lift().

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

◆ evaluateAx_lower()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
template<typename F>
FieldLHS ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::evaluateAx_lower ( FieldLHS & field,
F & evalFunction ) const

◆ evaluateAx_upper()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
template<typename F>
FieldLHS ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::evaluateAx_upper ( FieldLHS & field,
F & evalFunction ) const

◆ evaluateAx_upperlower()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
template<typename F>
FieldLHS ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::evaluateAx_upperlower ( FieldLHS & field,
F & evalFunction ) const

◆ evaluateLoadVector()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
void ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::evaluateLoadVector ( FieldRHS & field) const

◆ evaluateRefElementShapeFunction()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
KOKKOS_FUNCTION T ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::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 256 of file LagrangeSpace.hpp.

Referenced by computeAvg(), computeErrorL2(), and evaluateLoadVector().

Here is the caller graph for this function:

◆ evaluateRefElementShapeFunctionGradient()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
KOKKOS_FUNCTION LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::point_t ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::evaluateRefElementShapeFunctionGradient ( const size_t & localDOF,
const point_t & localPoint ) const

Evaluate the gradient of 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
point_t (Vector<T, Dim>) - The gradient of the shape function at the given point

Definition at line 291 of file LagrangeSpace.hpp.

Referenced by evaluateAx(), evaluateAx_diag(), evaluateAx_inversediag(), evaluateAx_lift(), evaluateAx_lower(), evaluateAx_upper(), and evaluateAx_upperlower().

Here is the caller graph for this function:

◆ getDeviceMirror()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::DeviceStruct ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::getDeviceMirror ( ) const

◆ getElementIndex()

KOKKOS_FUNCTION size_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FieldLHS, FieldRHS >::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::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::initializeElementIndices().

Here is the caller graph for this function:

◆ getElementMeshVertexIndices()

KOKKOS_FUNCTION FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FieldLHS, FieldRHS >::vertex_indices_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FieldLHS, FieldRHS >::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, FieldLHS, FieldRHS >::indices_list_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FieldLHS, FieldRHS >::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.

References Dim, and nr_m.

◆ getElementMeshVertexPoints()

KOKKOS_FUNCTION FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FieldLHS, FieldRHS >::vertex_points_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FieldLHS, FieldRHS >::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.

◆ getElementNDIndex()

KOKKOS_FUNCTION FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FieldLHS, FieldRHS >::indices_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FieldLHS, FieldRHS >::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.

References nr_m.

◆ getGlobalDOFIndex()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
KOKKOS_FUNCTION size_t ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::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, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >.

Definition at line 151 of file LagrangeSpace.hpp.

References getGlobalDOFIndex().

Referenced by getGlobalDOFIndex().

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

◆ getGlobalDOFIndices()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
KOKKOS_FUNCTION Vector< size_t, LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::numElementDOFs > ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::getGlobalDOFIndices ( const size_t & element_index) 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, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >.

Definition at line 178 of file LagrangeSpace.hpp.

Referenced by computeAvg(), computeErrorL2(), evaluateAx(), evaluateAx_diag(), evaluateAx_inversediag(), evaluateAx_lift(), evaluateAx_lower(), evaluateAx_upper(), evaluateAx_upperlower(), evaluateLoadVector(), and getLocalDOFIndex().

Here is the caller graph for this function:

◆ getLocalDOFIndex()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
KOKKOS_FUNCTION size_t ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::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, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >.

Definition at line 122 of file LagrangeSpace.hpp.

References Dim, and getGlobalDOFIndices().

Here is the call graph for this function:

◆ getLocalDOFIndices()

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
KOKKOS_FUNCTION Vector< size_t, LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::numElementDOFs > ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::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, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >.

Definition at line 163 of file LagrangeSpace.hpp.

References getLocalDOFIndices(), and numElementDOFs.

Referenced by getLocalDOFIndices().

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

◆ getMeshVertexIndex()

KOKKOS_FUNCTION size_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FieldLHS, FieldRHS >::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.

References Dim, and FiniteElementSpace().

Here is the call graph for this function:

◆ getMeshVertexNDIndex()

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

◆ initialize()

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

Initialize a LagrangeSpace object created with the default constructor.

Parameters
meshReference to the mesh
layoutReference to the field layout

Definition at line 38 of file LagrangeSpace.hpp.

References Dim, ippl::FiniteElementSpace< T, Dim, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >::FiniteElementSpace(), getLagrangeNumElementDOFs(), initializeElementIndices(), and ippl::FiniteElementSpace< T, Dim, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >::setMesh().

Here is the call graph for this function:

◆ initializeElementIndices()

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

Initialize a Kokkos view containing the element indices. This distributes the elements among MPI ranks.

Definition at line 52 of file LagrangeSpace.hpp.

References Dim, first(), ippl::FiniteElementSpace< T, Dim, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >::getElementIndex(), ippl::FieldLayout< Dim >::getLocalNDIndex(), initializeElementIndices(), ippl::FiniteElementSpace< T, Dim, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >::nr_m, and ippl::NDIndex< Dim >::size().

Referenced by initialize(), initializeElementIndices(), and LagrangeSpace().

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 FieldLHS, typename FieldRHS>
KOKKOS_FUNCTION bool ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::isDOFOnBoundary ( const indices_t & ndindex) const
inlineprivate

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

Parameters
ndindexThe NDIndex of the global DOF
Returns
true - If the DOF is on the boundary
false - If the DOF is not on the boundary

Definition at line 303 of file LagrangeSpace.h.

Referenced by evaluateAx(), evaluateAx_diag(), evaluateAx_inversediag(), evaluateAx_lift(), evaluateAx_lower(), evaluateAx_upper(), evaluateAx_upperlower(), and evaluateLoadVector().

Here is the caller graph for this function:

◆ numElements()

KOKKOS_FUNCTION size_t ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FieldLHS, FieldRHS >::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, FieldLHS, FieldRHS >::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 FieldLHS, typename FieldRHS>
KOKKOS_FUNCTION size_t ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::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, getLagrangeNumElementDOFs(Dim, Order), ElementType, QuadratureType, FieldLHS, FieldRHS >.

Definition at line 110 of file LagrangeSpace.hpp.

References numGlobalDOFs().

Referenced by numGlobalDOFs().

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

◆ setMesh()

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

Definition at line 77 of file FiniteElementSpace.hpp.

References nr_m.

Referenced by ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::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 FieldLHS, typename FieldRHS>
unsigned ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::dim
staticconstexpr

◆ elementIndices

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
Kokkos::View<size_t*> ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::elementIndices
private

Private member containing the element indices owned by //////////// my MPI rank. //////////////////////////////////////////////////////

Definition at line 316 of file LagrangeSpace.h.

Referenced by computeAvg(), computeErrorL2(), evaluateAx(), evaluateAx_diag(), evaluateAx_inversediag(), evaluateAx_lift(), evaluateAx_lower(), evaluateAx_upper(), evaluateAx_upperlower(), and evaluateLoadVector().

◆ hr_m

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

Definition at line 231 of file FiniteElementSpace.h.

◆ mesh_m

UniformCartesian<T, Dim>& ippl::FiniteElementSpace< T, Dim, NumElementDOFs, ElementType, QuadratureType, FieldLHS, FieldRHS >::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 FieldLHS, typename FieldRHS>
unsigned ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::numElementDOFs = getLagrangeNumElementDOFs(Dim, Order)
staticconstexpr

◆ numElementVertices

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
unsigned ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::numElementVertices
staticconstexpr
Initial value:

Definition at line 48 of file LagrangeSpace.h.

◆ order

template<typename T, unsigned Dim, unsigned Order, typename ElementType, typename QuadratureType, typename FieldLHS, typename FieldRHS>
unsigned ippl::LagrangeSpace< T, Dim, Order, ElementType, QuadratureType, FieldLHS, FieldRHS >::order = Order
staticconstexpr

Definition at line 45 of file LagrangeSpace.h.

◆ origin_m

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

Definition at line 232 of file FiniteElementSpace.h.

◆ quadrature_m

◆ ref_element_m


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