16 #ifndef Intrepid2_CellGeometry_h
17 #define Intrepid2_CellGeometry_h
31 #include "Intrepid2_ScalarView.hpp"
45 template<
class Po
intScalar,
int spaceDim,
typename DeviceType>
81 KOKKOS_INLINE_FUNCTION
110 Kokkos::Array<PointScalar,spaceDim> origin_;
111 Kokkos::Array<PointScalar,spaceDim> domainExtents_;
112 Kokkos::Array<int,spaceDim> gridCellCounts_;
119 ScalarView<int,DeviceType> cellToNodes_;
120 ScalarView<PointScalar,DeviceType> nodes_;
121 using BasisPtr = Teuchos::RCP< Basis<DeviceType,PointScalar,PointScalar> >;
123 unsigned numCells_ = 0;
124 unsigned numNodesPerCell_ = 0;
133 CellGeometry(
const Kokkos::Array<PointScalar,spaceDim> &origin,
134 const Kokkos::Array<PointScalar,spaceDim> &domainExtents,
135 const Kokkos::Array<int,spaceDim> &gridCellCounts,
147 ScalarView<int,DeviceType> cellToNodes,
148 ScalarView<PointScalar,DeviceType> nodes,
149 const bool claimAffine =
false,
157 ScalarView<PointScalar,DeviceType> cellNodes);
169 KOKKOS_INLINE_FUNCTION
193 KOKKOS_INLINE_FUNCTION
194 ordinal_type
cellToNode(ordinal_type cellIndex, ordinal_type cellNodeNumber)
const;
198 KOKKOS_INLINE_FUNCTION
202 KOKKOS_INLINE_FUNCTION
209 KOKKOS_INLINE_FUNCTION
210 size_t extent(
const int& r)
const;
213 template <
typename iType>
214 KOKKOS_INLINE_FUNCTION
215 typename std::enable_if<std::is_integral<iType>::value,
int>::type
219 KOKKOS_INLINE_FUNCTION
223 KOKKOS_INLINE_FUNCTION
227 KOKKOS_INLINE_FUNCTION
231 KOKKOS_INLINE_FUNCTION
235 KOKKOS_INLINE_FUNCTION
245 void orientations(ScalarView<Orientation,DeviceType> orientationsView,
const int &startCell = 0,
const int &endCell = -1);
248 KOKKOS_INLINE_FUNCTION
249 PointScalar
gridCellCoordinate(
const int &gridCellOrdinal,
const int &localNodeNumber,
const int &dim)
const;
252 KOKKOS_INLINE_FUNCTION
253 unsigned rank()
const;
256 KOKKOS_INLINE_FUNCTION
258 const int &subdivisionNodeNumber)
const;
261 KOKKOS_INLINE_FUNCTION
263 const int &subdivisionNodeNumber,
const int &d)
const;
266 KOKKOS_INLINE_FUNCTION
268 operator()(
const int& cell,
const int& node,
const int& dim)
const;
271 KOKKOS_INLINE_FUNCTION
318 const int startCell=0,
const int endCell=-1)
const;
328 const int startCell=0,
const int endCell=-1)
const;
340 #include <Intrepid2_CellGeometryDef.hpp>
void setJacobianDataPrivate(Data< PointScalar, DeviceType > &jacobianData, const int &pointsPerCell, const Data< PointScalar, DeviceType > &refData, const int startCell, const int endCell) const
Notionally-private method that provides a common interface for multiple public-facing setJacobianData...
geometry expressible in terms of vertices of the cell
BasisPtr basisForNodes() const
H^1 Basis used in the reference-to-physical transformation. Linear for straight-edged geometry; highe...
KOKKOS_INLINE_FUNCTION int numCellsPerGridCell(SubdivisionStrategy subdivisionStrategy) const
Helper method that returns the number of cells into which each grid cell will be subdivided based on ...
Data< PointScalar, DeviceType > allocateJacobianDataPrivate(const ScalarView< PointScalar, DeviceType > &pointComponentView, const int &pointsPerCell, const int startCell, const int endCell) const
Notionally-private method that provides a common interface for multiple public-facing allocateJacobia...
KOKKOS_INLINE_FUNCTION PointScalar subdivisionCoordinate(const int &gridCellOrdinal, const int &subdivisionOrdinal, const int &subdivisionNodeNumber, const int &d) const
returns coordinate in dimension d for the indicated subdivision of the indicated grid cell ...
CellGeometry provides the nodes for a set of cells; has options that support efficient definition of ...
View-like interface to tensor points; point components are stored separately; the appropriate coordin...
cube –> six tetrahedra
each grid division has the same dimensions
cube –> five tetrahedra
const shards::CellTopology & cellTopology() const
The shards CellTopology for each cell within the CellGeometry object. Note that this is always a lowe...
KOKKOS_INLINE_FUNCTION ordinal_type cellToNode(ordinal_type cellIndex, ordinal_type cellNodeNumber) const
Returns a global node number corresponding to the specified (cell, node) combination. If uniform grid (possibly subdivided), this number is computed dynamically; for more general meshes, this simply returns the result of a lookup in the cellToNodes container provided at construction.
KOKKOS_INLINE_FUNCTION PointScalar operator()(const int &cell, const int &node, const int &dim) const
Return the coordinate (weight) of the specified node. For straight-edged geometry, this is simply the physical coordinate of the vertex. For all geometries, this can be understood as a weight on the corresponding H^1 basis function used in the reference-to-physical map.
An abstract base class that defines interface for concrete basis implementations for Finite Element (...
void orientations(ScalarView< Orientation, DeviceType > orientationsView, const int &startCell=0, const int &endCell=-1)
Fills the provided container with the orientations for the specified cell range. Calls getOrientation...
View-like interface to tensor points; point components are stored separately; the appropriate coordin...
Defines the Data class, a wrapper around a Kokkos::View that allows data that is constant or repeatin...
CellGeometry(const Kokkos::Array< PointScalar, spaceDim > &origin, const Kokkos::Array< PointScalar, spaceDim > &domainExtents, const Kokkos::Array< int, spaceDim > &gridCellCounts, SubdivisionStrategy subdivisionStrategy=NO_SUBDIVISION, HypercubeNodeOrdering nodeOrdering=HYPERCUBE_NODE_ORDER_TENSOR)
Uniform grid constructor, with optional subdivision into simplices.
KOKKOS_INLINE_FUNCTION unsigned rank() const
Returns the logical rank of this container. This is always 3.
KOKKOS_INLINE_FUNCTION int hypercubeComponentNodeNumber(int hypercubeNodeNumber, int d) const
For hypercube vertex number hypercubeNodeNumber, returns the component node number in specified dimen...
Data< PointScalar, DeviceType > allocateJacobianData(const TensorPoints< PointScalar, DeviceType > &points, const int startCell=0, const int endCell=-1) const
Allocate a container into which Jacobians of the reference-to-physical mapping can be placed...
KOKKOS_INLINE_FUNCTION bool affine() const
Returns true if Jacobian is constant within each cell.
Wrapper around a Kokkos::View that allows data that is constant or repeating in various logical dimen...
KOKKOS_INLINE_FUNCTION HypercubeNodeOrdering nodeOrderingForHypercubes() const
Returns the node ordering used for hypercubes.
Header function for Intrepid2::Util class and other utility functions.
Header file for the Intrepid2::Orientation class.
KOKKOS_INLINE_FUNCTION int gridCellNodeForSubdivisionNode(const int &gridCellOrdinal, const int &subdivisionOrdinal, const int &subdivisionNodeNumber) const
returns coordinate in dimension d for the indicated subdivision of the indicated grid cell ...
TensorData< PointScalar, DeviceType > allocateCellMeasure(const Data< PointScalar, DeviceType > &jacobianDet, const TensorData< PointScalar, DeviceType > &cubatureWeights) const
Allocate a TensorData object appropriate for passing to computeCellMeasure().
KOKKOS_INLINE_FUNCTION PointScalar gridCellCoordinate(const int &gridCellOrdinal, const int &localNodeNumber, const int &dim) const
returns coordinate in dimension dim of the indicated node in the indicated grid cell ...
KOKKOS_INLINE_FUNCTION std::enable_if< std::is_integral< iType >::value, int >::type extent_int(const iType &r) const
Returns the logical extent of the container in the specified dimension as an int; the shape of CellGe...
lower-dimensional geometry that is orthogonally extruded in higher dimensions
KOKKOS_INLINE_FUNCTION int uniformJacobianModulus() const
Returns an integer indicating the number of distinct cell types vis-a-vis Jacobians.
Data< PointScalar, DeviceType > getJacobianRefData(const ScalarView< PointScalar, DeviceType > &points) const
Computes reference-space data for the specified points, to be used in setJacobian().
Orientation encoding and decoding.
geometry expressible in terms of a higher-order basis (must be specified)
Reference-space field values for a basis, designed to support typical vector-valued bases...
grid expressed as a Cartesian product of 1D grids (could be a Shishkin mesh, e.g.) ...
void initializeOrientations()
Initialize the internal orientations_ member with the orientations of each member cell...
square –> four triangles, with a new vertex at center
Data< Orientation, DeviceType > getOrientations()
Returns the orientations for all cells. Calls initializeOrientations() if it has not previously been ...
KOKKOS_INLINE_FUNCTION Orientation getOrientation(int &cellNumber) const
Returns the orientation for the specified cell. Requires that initializeOrientations() has been calle...
KOKKOS_INLINE_FUNCTION DataVariationType cellVariationType() const
KOKKOS_INLINE_FUNCTION int numNodesPerCell() const
Returns the number of nodes per cell; may be more than the number of vertices in the corresponding Ce...
KOKKOS_INLINE_FUNCTION int numCellsInDimension(const int &dim) const
For uniform grid and tensor grid CellGeometry, returns the number of cells in the specified component...
void setJacobian(Data< PointScalar, DeviceType > &jacobianData, const TensorPoints< PointScalar, DeviceType > &points, const Data< PointScalar, DeviceType > &refData, const int startCell=0, const int endCell=-1) const
Compute Jacobian values for the reference-to-physical transformation, and place them in the provided ...
a more natural tensor ordering
square –> two triangles, with a hypotenuse of slope -1
square –> two triangles, with a hypotenuse of slope 1
View-like interface to tensor data; tensor components are stored separately and multiplied together a...
KOKKOS_INLINE_FUNCTION int numCells() const
Returns the number of cells.
KOKKOS_INLINE_FUNCTION ~CellGeometry()
Destructor.
KOKKOS_INLINE_FUNCTION size_t extent(const int &r) const
Returns the logical extent of the container in the specified dimension; the shape of CellGeometry is ...
Stateless class that acts as a factory for a family of nodal bases (hypercube topologies only at this...
Header file for the abstract base class Intrepid2::Basis.
void computeCellMeasure(TensorData< PointScalar, DeviceType > &cellMeasure, const Data< PointScalar, DeviceType > &jacobianDet, const TensorData< PointScalar, DeviceType > &cubatureWeights) const
Compute cell measures that correspond to provided Jacobian determinants and.