Intrepid2
Intrepid2_CellGeometry.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Intrepid2 Package
4 //
5 // Copyright 2007 NTESS and the Intrepid2 contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
16 #ifndef Intrepid2_CellGeometry_h
17 #define Intrepid2_CellGeometry_h
18 
19 #include "Intrepid2_Basis.hpp"
20 #include "Intrepid2_CellTools.hpp"
21 #include "Intrepid2_Data.hpp"
25 #include "Intrepid2_TensorData.hpp"
28 #include "Intrepid2_Utils.hpp"
29 #include "Intrepid2_VectorData.hpp"
30 
31 #include "Intrepid2_ScalarView.hpp"
32 
33 namespace Intrepid2
34 {
45  template<class PointScalar, int spaceDim, typename DeviceType>
47  {
48  public:
56  };
57 
60  {
68  };
69 
72  {
75  };
76 
81  KOKKOS_INLINE_FUNCTION
82  int numCellsPerGridCell(SubdivisionStrategy subdivisionStrategy) const;
83 
84  public:
92  Data<PointScalar,DeviceType> allocateJacobianDataPrivate(const ScalarView<PointScalar,DeviceType> &pointComponentView, const int &pointsPerCell, const int startCell, const int endCell) const;
93 
101  void setJacobianDataPrivate(Data<PointScalar,DeviceType> &jacobianData, const int &pointsPerCell, const Data<PointScalar,DeviceType> &refData, const int startCell, const int endCell) const;
102  protected:
103  HypercubeNodeOrdering nodeOrdering_;
104  CellGeometryType cellGeometryType_;
105  SubdivisionStrategy subdivisionStrategy_ = NO_SUBDIVISION;
106  bool affine_; // if true, each cell has constant Jacobian across the cell
107  Data<Orientation, DeviceType> orientations_; // for grid types, this could have either a single entry or one matching numCellsPerGridCell(). For other types, it has as many entries as there are cells.
108 
109  // uniform grid data -- used for UNIFORM_GRID type
110  Kokkos::Array<PointScalar,spaceDim> origin_; // point specifying a corner of the mesh
111  Kokkos::Array<PointScalar,spaceDim> domainExtents_; // how far the domain extends in each dimension
112  Kokkos::Array<int,spaceDim> gridCellCounts_; // how many grid cells wide the mesh is in each dimension
113 
114  // tensor grid data -- only used for TENSOR_GRID type
116 
117  // arbitrary cell node data, used for both higher-order and first-order
118  // (here, nodes are understood as geometry degrees of freedom)
119  ScalarView<int,DeviceType> cellToNodes_; // (C,N) -- N is the number of nodes per cell; values are global node ordinals
120  ScalarView<PointScalar,DeviceType> nodes_; // (GN,D) or (C,N,D) -- GN is the number of global nodes; (C,N,D) used only if cellToNodes_ is empty.
121  using BasisPtr = Teuchos::RCP< Basis<DeviceType,PointScalar,PointScalar> >;
122 
123  unsigned numCells_ = 0;
124  unsigned numNodesPerCell_ = 0;
125  public:
133  CellGeometry(const Kokkos::Array<PointScalar,spaceDim> &origin,
134  const Kokkos::Array<PointScalar,spaceDim> &domainExtents,
135  const Kokkos::Array<int,spaceDim> &gridCellCounts,
136  SubdivisionStrategy subdivisionStrategy = NO_SUBDIVISION,
138 
146  CellGeometry(const shards::CellTopology &cellTopo,
147  ScalarView<int,DeviceType> cellToNodes,
148  ScalarView<PointScalar,DeviceType> nodes,
149  const bool claimAffine = false,
151 
157  ScalarView<PointScalar,DeviceType> cellNodes);
158 
162  KOKKOS_INLINE_FUNCTION CellGeometry(const CellGeometry &cellGeometry);
163 
166  KOKKOS_INLINE_FUNCTION ~CellGeometry();
167 
169  KOKKOS_INLINE_FUNCTION
170  bool affine() const;
171 
178 
184  void computeCellMeasure( TensorData<PointScalar,DeviceType> &cellMeasure, const Data<PointScalar,DeviceType> & jacobianDet, const TensorData<PointScalar,DeviceType> & cubatureWeights ) const;
185 
187  BasisPtr basisForNodes() const;
188 
190  const shards::CellTopology & cellTopology() const;
191 
193  KOKKOS_INLINE_FUNCTION
194  ordinal_type cellToNode(ordinal_type cellIndex, ordinal_type cellNodeNumber) const;
195 
198  KOKKOS_INLINE_FUNCTION
199  DataVariationType cellVariationType() const;
200 
202  KOKKOS_INLINE_FUNCTION
203  int hypercubeComponentNodeNumber(int hypercubeNodeNumber, int d) const;
204 
206  void initializeOrientations();
207 
209  KOKKOS_INLINE_FUNCTION
210  size_t extent(const int& r) const;
211 
213  template <typename iType>
214  KOKKOS_INLINE_FUNCTION
215  typename std::enable_if<std::is_integral<iType>::value, int>::type
216  extent_int(const iType& r) const;
217 
219  KOKKOS_INLINE_FUNCTION
221 
223  KOKKOS_INLINE_FUNCTION
224  int numCells() const;
225 
227  KOKKOS_INLINE_FUNCTION
228  int numCellsInDimension(const int &dim) const;
229 
231  KOKKOS_INLINE_FUNCTION
232  int numNodesPerCell() const;
233 
235  KOKKOS_INLINE_FUNCTION
236  Orientation getOrientation(int &cellNumber) const;
237 
240 
245  void orientations(ScalarView<Orientation,DeviceType> orientationsView, const int &startCell = 0, const int &endCell = -1);
246 
248  KOKKOS_INLINE_FUNCTION
249  PointScalar gridCellCoordinate(const int &gridCellOrdinal, const int &localNodeNumber, const int &dim) const;
250 
252  KOKKOS_INLINE_FUNCTION
253  unsigned rank() const;
254 
256  KOKKOS_INLINE_FUNCTION
257  int gridCellNodeForSubdivisionNode(const int &gridCellOrdinal, const int &subdivisionOrdinal,
258  const int &subdivisionNodeNumber) const;
259 
261  KOKKOS_INLINE_FUNCTION
262  PointScalar subdivisionCoordinate(const int &gridCellOrdinal, const int &subdivisionOrdinal,
263  const int &subdivisionNodeNumber, const int &d) const;
264 
266  KOKKOS_INLINE_FUNCTION
267  PointScalar
268  operator()(const int& cell, const int& node, const int& dim) const;
269 
271  KOKKOS_INLINE_FUNCTION
272  int uniformJacobianModulus() const;
273 
280  Data<PointScalar,DeviceType> allocateJacobianData(const TensorPoints<PointScalar,DeviceType> &points, const int startCell=0, const int endCell=-1) const;
281 
288  Data<PointScalar,DeviceType> allocateJacobianData(const ScalarView<PointScalar,DeviceType> &points, const int startCell=0, const int endCell=-1) const;
289 
296  Data<PointScalar,DeviceType> allocateJacobianData(const int &numPoints, const int startCell=0, const int endCell=-1) const;
297 
302  Data<PointScalar,DeviceType> getJacobianRefData(const ScalarView<PointScalar,DeviceType> &points) const;
303 
309 
318  const int startCell=0, const int endCell=-1) const;
319 
327  void setJacobian(Data<PointScalar,DeviceType> &jacobianData, const ScalarView<PointScalar,DeviceType> &points, const Data<PointScalar,DeviceType> &refData,
328  const int startCell=0, const int endCell=-1) const;
329 
336  void setJacobian(Data<PointScalar,DeviceType> &jacobianData, const int &numPoints, const int startCell=0, const int endCell=-1) const;
337  };
338 } // namespace Intrepid2
339 
340 #include <Intrepid2_CellGeometryDef.hpp>
341 
342 #endif /* Intrepid2_CellGeometry_h */
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...
each grid division has the same dimensions
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.
Header file for the Intrepid2::OrientationTools and Intrepid2::Impl::OrientationTools classes...
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 –&gt; 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 ...
Structure-preserving representation of transformed basis values; reference space values and transform...
square –&gt; two triangles, with a hypotenuse of slope -1
square –&gt; 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 ...
Header file for the Intrepid2::CellTools class.
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.