|
IPPL (Independent Parallel Particle Layer)
IPPL
|
A class representing a Nedelec space for finite element methods on a structured, rectilinear grid. More...
#include <NedelecSpace.h>
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< Dim > | Layout_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, numElementVertices > | vertex_indices_t |
| typedef Vector< indices_t, numElementVertices > | indices_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, numElementDOFs > | getLocalDOFIndices () 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, numElementDOFs > | getGlobalDOFIndices (const size_t &elementIndex) const override |
| Get the global DOF indices (vector of global DOF indices) of an element. | |
| KOKKOS_FUNCTION Vector< size_t, numElementDOFs > | getGlobalDOFIndices (const indices_t &elementIndex) const |
| Get the global DOF indices (vector of global DOF indices) of an element. | |
| KOKKOS_FUNCTION Vector< size_t, numElementDOFs > | 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. | |
| KOKKOS_FUNCTION Vector< size_t, numElementDOFs > | 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. | |
| 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< T > | evaluateAx (FEMVector< T > &x, F &evalFunction) const |
| Assembly operations ///////////////////////////////////////////////. | |
| FEMVector< T > | 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. | |
| template<typename F> | |
| FEMVector< T > | evaluateLoadVectorFunctor (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< T > | createFEMVector () 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, Dim > | nr_m |
| Vector< double, Dim > | hr_m |
| Vector< double, Dim > | origin_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< T > | createFEMVector2d () const |
Implementation of the NedelecSpace::createFEMVector function for 2d. | |
| FEMVector< T > | createFEMVector3d () 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. | |
A class representing a Nedelec space for finite element methods on a structured, rectilinear grid.
| T | The floating point number type of the field values |
| Dim | The dimension of the mesh |
| Order | The order of the Nedelec space |
| QuadratureType | The type of the quadrature rule |
| FieldType | The type of field to use. |
Definition at line 34 of file NedelecSpace.h.
| 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.
|
inherited |
Definition at line 59 of file FiniteElementSpace.h.
| 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.
| typedef FieldLayout<Dim> ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::Layout_t |
Definition at line 63 of file NedelecSpace.h.
| 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.
|
inherited |
Definition at line 57 of file FiniteElementSpace.h.
| 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.
| typedef detail::ViewType<T,1>::view_type ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::ViewType |
Definition at line 66 of file NedelecSpace.h.
| 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.
| mesh | Reference to the mesh |
| ref_element | Reference to the reference element |
| quadrature | Reference to the quadrature rule |
| layout | Reference 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().
| 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.
| mesh | Reference to the mesh |
| ref_element | Reference to the reference element |
| quadrature | Reference 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().
| 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.
| u_h | The basis function coefficients obtained via FEM. |
| u_sol | The analytical solution (functor) |
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.
| 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.
Definition at line 804 of file NedelecSpace.hpp.
References createFEMVector2d(), createFEMVector3d(), and Dim.
Referenced by evaluateLoadVector(), and evaluateLoadVectorFunctor().
|
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().
|
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().
| 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
| x | The vector which we want to multiply |
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().
| 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.
| f | The source field defined at the Nédélec degrees fo freedom. |
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.
| FEMVector< T > ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateLoadVectorFunctor | ( | 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.
| f | The source function, which can be evaluated at arbitrary points. |
| F | The functor type. |
Definition at line 599 of file NedelecSpace.hpp.
References createFEMVector(), ippl::Vector< T, Dim >::dot(), elementIndices, evaluateRefElementShapeFunction(), ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::getElementMeshVertexPoints(), ippl::FiniteElementSpace< T, Dim, getNedelecNumElementDOFs(Dim, Order), ElementType, QuadratureType, FEMVector< T >, FEMVector< T > >::getElementNDIndex(), getFEMVectorDOFIndices(), getGlobalDOFIndices(), ippl::FEMVector< T >::getView(), isDOFOnBoundary(), 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.
| 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
| localDOF | size_t - The local degree of freedom index |
| localPoint | point_t (Vector<T, Dim>) - The point in the reference element |
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().
| 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.
| localDOF | size_t - The local degree of freedom index |
| localPoint | point_t (Vector<T, Dim>) - The point in the reference element |
Definition at line 750 of file NedelecSpace.hpp.
References Dim, and NedelecSpace().
Referenced by evaluateAx().
| 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.
| dofIdx | the global DoF index for which the boundary side should be retrieved. |
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.
|
inherited |
Get the global index of a mesh element given the NDIndex.
| ndindex | indices_t (Vector<size_t, Dim>) - vector of indices for each direction |
Definition at line 135 of file FiniteElementSpace.hpp.
Referenced by ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::initializeElementIndices().
|
inherited |
Get all the global vertex indices of an element (given by its NDIndex).
| elementNDIndex | The NDIndex of the element |
Definition at line 146 of file FiniteElementSpace.hpp.
|
inherited |
Get all the NDIndices of the vertices of an element (given by its NDIndex).
| elementNDIndex | The NDIndex of the element |
Definition at line 157 of file FiniteElementSpace.hpp.
|
inherited |
Get all the global vertex points of an element (given by its NDIndex).
| elementNDIndex | The NDIndex of the element |
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().
|
inherited |
Get the NDIndex (vector of indices for each dimension) of a mesh element.
| elementIndex | indices_t (Vector<size_t, Dim>) - The index of the element |
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().
| 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.
| elementIndex | size_t - The index of the element |
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().
| 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.
| elementIndex | indices_t - The index of the element |
|
overridevirtual |
Get the global DOF index from the element index and local DOF.
| elementIndex | size_t - The index of the element |
| localDOFIndex | size_t - The local DOF index |
Definition at line 180 of file NedelecSpace.hpp.
References getGlobalDOFIndices().
| 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.
| elementIndex | indices_t - The index of the element |
|
overridevirtual |
Get the global DOF indices (vector of global DOF indices) of an element.
| elementIndex | size_t - The index of the element |
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().
|
overridevirtual |
Get the elements local DOF from the element index and global DOF index.
| elementIndex | size_t - The index of the element |
| globalDOFIndex | size_t - The global DOF index |
Definition at line 149 of file NedelecSpace.hpp.
References Dim, ippl::Vector< T, Dim >::dim, and getGlobalDOFIndices().
|
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.
Definition at line 193 of file NedelecSpace.hpp.
References numElementDOFs.
| 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().
|
inherited |
Get the global index of a mesh vertex given its NDIndex.
| vertex_nd_index | indices_t (Vector<size_t, Dim>) - The NDIndex of the vertex (vector of indices for each dimension). |
Definition at line 117 of file FiniteElementSpace.hpp.
|
inherited |
Get the NDIndex of a mesh vertex.
| vertex_index | size_t - The index of the vertex |
Definition at line 107 of file FiniteElementSpace.hpp.
| 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.
| mesh | Reference to the mesh |
| layout | Reference 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().
| 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().
| 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.
| dofIdx | The global DoF index for which should be checked if it is on the 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().
|
inherited |
Mesh and Element operations ///////////////////////////////////////.
Get the number of elements in the mesh of the space
Definition at line 88 of file FiniteElementSpace.hpp.
|
inherited |
Get the number of elements in a given dimension.
| dim | size_t - representing the dimension |
Definition at line 97 of file FiniteElementSpace.hpp.
|
overridevirtual |
Degree of Freedom operations //////////////////////////////////////.
Get the number of global degrees of freedom in the space
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.
| 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.
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.| positions | The points at which the function should be evaluated. A Kokkos::View which stores in each element a 2D/3D point. |
| coef | The basis function coefficients obtained via FEM. |
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.
|
inherited |
Definition at line 77 of file FiniteElementSpace.hpp.
Referenced by ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::initialize().
|
staticconstexpr |
Definition at line 41 of file NedelecSpace.h.
|
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().
|
inherited |
Definition at line 231 of file FiniteElementSpace.h.
Referenced by ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::reconstructToPoints().
|
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().
|
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().
|
inherited |
Member variables //////////////////////////////////////////////////.
Definition at line 227 of file FiniteElementSpace.h.
|
inherited |
Definition at line 230 of file FiniteElementSpace.h.
Referenced by ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::getBoundarySide(), ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::initializeElementIndices(), ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::isDOFOnBoundary(), ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::numGlobalDOFs(), and ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::reconstructToPoints().
|
staticconstexpr |
Definition at line 38 of file NedelecSpace.h.
Referenced by computeError(), evaluateAx(), evaluateLoadVector(), evaluateLoadVectorFunctor(), evaluateRefElementShapeFunction(), getLocalDOFIndices(), and reconstructToPoints().
|
staticconstexpr |
Definition at line 48 of file NedelecSpace.h.
|
staticconstexpr |
Definition at line 45 of file NedelecSpace.h.
|
inherited |
Definition at line 232 of file FiniteElementSpace.h.
Referenced by ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::reconstructToPoints().
|
inherited |
Definition at line 229 of file FiniteElementSpace.h.
Referenced by ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::computeError(), ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateAx(), ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateLoadVector(), and ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateLoadVectorFunctor().
|
inherited |
Definition at line 228 of file FiniteElementSpace.h.
Referenced by ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::computeError(), ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateLoadVector(), ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateLoadVectorFunctor(), and ippl::NedelecSpace< T, Dim, Order, ElementType, QuadratureType, FieldType >::evaluateRefElementShapeFunction().