Intrepid2
Intrepid2_CellTools.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ************************************************************************
3 //
4 // Intrepid2 Package
5 // Copyright (2007) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Kyungjoo Kim (kyukim@sandia.gov), or
38 // Mauro Perego (mperego@sandia.gov)
39 //
40 // ************************************************************************
41 // @HEADER
42 
49 #ifndef __INTREPID2_CELLTOOLS_HPP__
50 #define __INTREPID2_CELLTOOLS_HPP__
51 
52 #include "Intrepid2_ConfigDefs.hpp"
53 
54 #include "Shards_CellTopology.hpp"
55 #include "Shards_BasicTopologies.hpp"
56 
57 #include "Teuchos_RCP.hpp"
58 
59 #include "Intrepid2_Types.hpp"
60 #include "Intrepid2_Utils.hpp"
61 
63 
64 #include "Intrepid2_Basis.hpp"
65 
69 
74 
77 
81 
83 //#include "Intrepid2_HGRAD_WEDGE_I2_FEM.hpp"
84 
85 namespace Intrepid2 {
86 
87  //============================================================================================//
88  // //
89  // CellTools //
90  // //
91  //============================================================================================//
92 
103  template<typename ExecSpaceType>
104  class CellTools {
105  public:
106 
111  inline
112  static bool
113  hasReferenceCell( const shards::CellTopology cellTopo ) {
114  bool r_val = false;
115  switch ( cellTopo.getKey() ) {
116  case shards::Line<2>::key:
117  case shards::Line<3>::key:
118  case shards::ShellLine<2>::key:
119  case shards::ShellLine<3>::key:
120  case shards::Beam<2>::key:
121  case shards::Beam<3>::key:
122 
123  case shards::Triangle<3>::key:
124  // case shards::Triangle<4>::key:
125  case shards::Triangle<6>::key:
126  // case shards::ShellTriangle<3>::key:
127  // case shards::ShellTriangle<6>::key:
128 
129  case shards::Quadrilateral<4>::key:
130  case shards::Quadrilateral<8>::key:
131  case shards::Quadrilateral<9>::key:
132  //case shards::ShellQuadrilateral<4>::key:
133  //case shards::ShellQuadrilateral<8>::key:
134  //case shards::ShellQuadrilateral<9>::key:
135 
136  case shards::Tetrahedron<4>::key:
137  // case shards::Tetrahedron<8>::key:
138  case shards::Tetrahedron<10>::key:
139  // case shards::Tetrahedron<11>::key:
140 
141  case shards::Hexahedron<8>::key:
142  case shards::Hexahedron<20>::key:
143  case shards::Hexahedron<27>::key:
144 
145  case shards::Pyramid<5>::key:
146  // case shards::Pyramid<13>::key:
147  // case shards::Pyramid<14>::key:
148 
149  case shards::Wedge<6>::key:
150  // case shards::Wedge<15>::key:
151  case shards::Wedge<18>::key:
152  r_val = true;
153  }
154  return r_val;
155  }
156 
157  private:
158 
162  template<typename outputValueType,
163  typename pointValueType>
164  static Teuchos::RCP<Basis<ExecSpaceType,outputValueType,pointValueType> >
165  createHGradBasis( const shards::CellTopology cellTopo ) {
166  Teuchos::RCP<Basis<ExecSpaceType,outputValueType,pointValueType> > r_val;
167 
168  switch (cellTopo.getKey()) {
169  case shards::Line<2>::key: r_val = Teuchos::rcp(new Basis_HGRAD_LINE_C1_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
170  case shards::Triangle<3>::key: r_val = Teuchos::rcp(new Basis_HGRAD_TRI_C1_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
171  case shards::Quadrilateral<4>::key: r_val = Teuchos::rcp(new Basis_HGRAD_QUAD_C1_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
172  case shards::Tetrahedron<4>::key: r_val = Teuchos::rcp(new Basis_HGRAD_TET_C1_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
173  case shards::Hexahedron<8>::key: r_val = Teuchos::rcp(new Basis_HGRAD_HEX_C1_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
174  case shards::Wedge<6>::key: r_val = Teuchos::rcp(new Basis_HGRAD_WEDGE_C1_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
175  case shards::Pyramid<5>::key: r_val = Teuchos::rcp(new Basis_HGRAD_PYR_C1_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
176 
177  case shards::Triangle<6>::key: r_val = Teuchos::rcp(new Basis_HGRAD_TRI_C2_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
178  case shards::Quadrilateral<9>::key: r_val = Teuchos::rcp(new Basis_HGRAD_QUAD_C2_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
179  case shards::Tetrahedron<10>::key: r_val = Teuchos::rcp(new Basis_HGRAD_TET_C2_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
180  case shards::Tetrahedron<11>::key: r_val = Teuchos::rcp(new Basis_HGRAD_TET_COMP12_FEM<ExecSpaceType,outputValueType,pointValueType>()); break;
181  //case shards::Hexahedron<20>::key: r_val = Teuchos::rcp(new Basis_HGRAD_HEX_I2_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
182  case shards::Hexahedron<27>::key: r_val = Teuchos::rcp(new Basis_HGRAD_HEX_C2_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
183  //case shards::Wedge<15>::key: r_val = Teuchos::rcp(new Basis_HGRAD_WEDGE_I2_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
184  case shards::Wedge<18>::key: r_val = Teuchos::rcp(new Basis_HGRAD_WEDGE_C2_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
185  //case shards::Pyramid<13>::key: r_val = Teuchos::rcp(new Basis_HGRAD_PYR_I2_FEM <ExecSpaceType,outputValueType,pointValueType>()); break;
186 
187  case shards::Quadrilateral<8>::key:
188  case shards::Line<3>::key:
189  case shards::Beam<2>::key:
190  case shards::Beam<3>::key:
191  case shards::ShellLine<2>::key:
192  case shards::ShellLine<3>::key:
193  case shards::ShellTriangle<3>::key:
194  case shards::ShellTriangle<6>::key:
195  case shards::ShellQuadrilateral<4>::key:
196  case shards::ShellQuadrilateral<8>::key:
197  case shards::ShellQuadrilateral<9>::key:
198  default: {
199  INTREPID2_TEST_FOR_EXCEPTION( true, std::invalid_argument,
200  ">>> ERROR (Intrepid2::CellTools::createHGradBasis): Cell topology not supported.");
201  }
202  }
203  return r_val;
204  }
205 
206  // reference nodes initialized
210  typedef Kokkos::DynRankView<const double,Kokkos::LayoutRight,Kokkos::HostSpace> referenceNodeDataViewHostType;
212  double line[2][3], line_3[3][3];
213  double triangle[3][3], triangle_4[4][3], triangle_6[6][3];
214  double quadrilateral[4][3], quadrilateral_8[8][3], quadrilateral_9[9][3];
215  double tetrahedron[4][3], tetrahedron_8[8][3], tetrahedron_10[10][3], tetrahedron_11[10][3];
216  double hexahedron[8][3], hexahedron_20[20][3], hexahedron_27[27][3];
217  double pyramid[5][3], pyramid_13[13][3], pyramid_14[14][3];
218  double wedge[6][3], wedge_15[15][3], wedge_18[18][3];
219  };
220 
224  typedef Kokkos::DynRankView<double,Kokkos::LayoutRight,ExecSpaceType> referenceNodeDataViewType;
226  referenceNodeDataViewType line, line_3;
227  referenceNodeDataViewType triangle, triangle_4, triangle_6;
228  referenceNodeDataViewType quadrilateral, quadrilateral_8, quadrilateral_9;
229  referenceNodeDataViewType tetrahedron, tetrahedron_8, tetrahedron_10, tetrahedron_11;
230  referenceNodeDataViewType hexahedron, hexahedron_20, hexahedron_27;
231  referenceNodeDataViewType pyramid, pyramid_13, pyramid_14;
232  referenceNodeDataViewType wedge, wedge_15, wedge_18;
233  };
234 
235  static const ReferenceNodeDataStatic refNodeDataStatic_;
236  static ReferenceNodeData refNodeData_;
237 
238  static bool isReferenceNodeDataSet_;
241  static void setReferenceNodeData();
242 
243  //============================================================================================//
244  // //
245  // Parametrization coefficients of edges and faces of reference cells //
246  // //
247  //============================================================================================//
248 
249  // static variables
253  typedef Kokkos::DynRankView<double,ExecSpaceType> subcellParamViewType;
255  subcellParamViewType dummy;
256  subcellParamViewType lineEdges; // edge maps for 2d non-standard cells; shell line and beam
257  subcellParamViewType triEdges, quadEdges; // edge maps for 2d standard cells
258  subcellParamViewType shellTriEdges, shellQuadEdges; // edge maps for 3d non-standard cells; shell tri and quad
259  subcellParamViewType tetEdges, hexEdges, pyrEdges, wedgeEdges; // edge maps for 3d standard cells
260  subcellParamViewType shellTriFaces, shellQuadFaces; // face maps for 3d non-standard cells
261  subcellParamViewType tetFaces, hexFaces, pyrFaces, wedgeFaces; // face maps for 3d standard cells
262  };
263 
264  static SubcellParamData subcellParamData_;
265 
266  static bool isSubcellParametrizationSet_;
267 
302  static void setSubcellParametrization();
303 
316  static void
317  getSubcellParametrization( subcellParamViewType &subcellParam,
318  const ordinal_type subcellDim,
319  const shards::CellTopology parentCell );
320 
321 
332  static void
333  setSubcellParametrization( subcellParamViewType &subcellParam,
334  const ordinal_type subcellDim,
335  const shards::CellTopology parentCell );
336 
337  public:
338 
341  CellTools() = default;
342 
345  ~CellTools() = default;
346 
347  //============================================================================================//
348  // //
349  // Jacobian, inverse Jacobian and Jacobian determinant //
350  // //
351  //============================================================================================//
352 
386  template<typename jacobianValueType, class ...jacobianProperties,
387  typename pointValueType, class ...pointProperties,
388  typename worksetCellValueType, class ...worksetCellProperties,
389  typename HGradBasisPtrType>
390  static void
391  setJacobian( Kokkos::DynRankView<jacobianValueType,jacobianProperties...> jacobian,
392  const Kokkos::DynRankView<pointValueType,pointProperties...> points,
393  const Kokkos::DynRankView<worksetCellValueType,worksetCellProperties...> worksetCell,
394  const HGradBasisPtrType basis );
395 
430  template<typename jacobianValueType, class ...jacobianProperties,
431  typename pointValueType, class ...pointProperties,
432  typename worksetCellValueType, class ...worksetCellProperties>
433  static void
434  setJacobian( Kokkos::DynRankView<jacobianValueType,jacobianProperties...> jacobian,
435  const Kokkos::DynRankView<pointValueType,pointProperties...> points,
436  const Kokkos::DynRankView<worksetCellValueType,worksetCellProperties...> worksetCell,
437  const shards::CellTopology cellTopo ) {
438  auto basis = createHGradBasis<pointValueType,pointValueType>(cellTopo);
439  setJacobian(jacobian,
440  points,
441  worksetCell,
442  basis);
443  }
444 
455  template<typename jacobianInvValueType, class ...jacobianInvProperties,
456  typename jacobianValueType, class ...jacobianProperties>
457  static void
458  setJacobianInv( Kokkos::DynRankView<jacobianInvValueType,jacobianInvProperties...> jacobianInv,
459  const Kokkos::DynRankView<jacobianValueType,jacobianProperties...> jacobian );
460 
471  template<typename jacobianDetValueType, class ...jacobianDetProperties,
472  typename jacobianValueType, class ...jacobianProperties>
473  static void
474  setJacobianDet( Kokkos::DynRankView<jacobianDetValueType,jacobianDetProperties...> jacobianDet,
475  const Kokkos::DynRankView<jacobianValueType,jacobianProperties...> jacobian );
476 
477 
478  //============================================================================================//
479  // //
480  // Node information //
481  // //
482  //============================================================================================//
483 
484  // the node information can be used inside of kokkos functor and needs kokkos inline and
485  // exception should be an abort. for now, let's not decorate
486 
497  template<typename cellCenterValueType, class ...cellCenterProperties,
498  typename cellVertexValueType, class ...cellVertexProperties>
499  static void
500  getReferenceCellCenter( Kokkos::DynRankView<cellCenterValueType,cellCenterProperties...> cellCenter,
501  Kokkos::DynRankView<cellVertexValueType,cellVertexProperties...> cellVertex,
502  const shards::CellTopology cell );
503 
514  template<typename cellVertexValueType, class ...cellVertexProperties>
515  static void
516  getReferenceVertex( Kokkos::DynRankView<cellVertexValueType,cellVertexProperties...> cellVertex,
517  const shards::CellTopology cell,
518  const ordinal_type vertexOrd );
519 
520 
535  template<typename subcellVertexValueType, class ...subcellVertexProperties>
536  static void
537  getReferenceSubcellVertices( Kokkos::DynRankView<subcellVertexValueType,subcellVertexProperties...> subcellVertices,
538  const ordinal_type subcellDim,
539  const ordinal_type subcellOrd,
540  const shards::CellTopology parentCell );
541 
542 
543 
559  template<typename cellNodeValueType, class ...cellNodeProperties>
560  static void
561  getReferenceNode( Kokkos::DynRankView<cellNodeValueType,cellNodeProperties...> cellNode,
562  const shards::CellTopology cell,
563  const ordinal_type nodeOrd );
564 
565 
566 
580  template<typename subcellNodeValueType, class ...subcellNodeProperties>
581  static void
582  getReferenceSubcellNodes( Kokkos::DynRankView<subcellNodeValueType,subcellNodeProperties...> subcellNodes,
583  const ordinal_type subcellDim,
584  const ordinal_type subcellOrd,
585  const shards::CellTopology parentCell );
586 
612  template<typename refEdgeTangentValueType, class ...refEdgeTangentProperties>
613  static void
614  getReferenceEdgeTangent( Kokkos::DynRankView<refEdgeTangentValueType,refEdgeTangentProperties...> refEdgeTangent,
615  const ordinal_type edgeOrd,
616  const shards::CellTopology parentCell );
617 
654  template<typename refFaceTanUValueType, class ...refFaceTanUProperties,
655  typename refFaceTanVValueType, class ...refFaceTanVProperties>
656  static void
657  getReferenceFaceTangents( Kokkos::DynRankView<refFaceTanUValueType,refFaceTanUProperties...> refFaceTanU,
658  Kokkos::DynRankView<refFaceTanVValueType,refFaceTanVProperties...> refFaceTanV,
659  const ordinal_type faceOrd,
660  const shards::CellTopology parentCell );
661 
724  template<typename refSideNormalValueType, class ...refSideNormalProperties>
725  static void
726  getReferenceSideNormal( Kokkos::DynRankView<refSideNormalValueType,refSideNormalProperties...> refSideNormal,
727  const ordinal_type sideOrd,
728  const shards::CellTopology parentCell );
729 
768  template<typename refFaceNormalValueType, class ...refFaceNormalProperties>
769  static void
770  getReferenceFaceNormal( Kokkos::DynRankView<refFaceNormalValueType,refFaceNormalProperties...> refFaceNormal,
771  const ordinal_type faceOrd,
772  const shards::CellTopology parentCell );
773 
803  template<typename edgeTangentValueType, class ...edgeTangentProperties,
804  typename worksetJacobianValueType, class ...worksetJacobianProperties>
805  static void
806  getPhysicalEdgeTangents( Kokkos::DynRankView<edgeTangentValueType,edgeTangentProperties...> edgeTangents,
807  const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
808  const ordinal_type worksetEdgeOrd,
809  const shards::CellTopology parentCell );
810 
850  template<typename faceTanUValueType, class ...faceTanUProperties,
851  typename faceTanVValueType, class ...faceTanVProperties,
852  typename worksetJacobianValueType, class ...worksetJacobianProperties>
853  static void
854  getPhysicalFaceTangents( Kokkos::DynRankView<faceTanUValueType,faceTanUProperties...> faceTanU,
855  Kokkos::DynRankView<faceTanVValueType,faceTanVProperties...> faceTanV,
856  const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
857  const ordinal_type worksetFaceOrd,
858  const shards::CellTopology parentCell );
859 
920  template<typename sideNormalValueType, class ...sideNormalProperties,
921  typename worksetJacobianValueType, class ...worksetJacobianProperties>
922  static void
923  getPhysicalSideNormals( Kokkos::DynRankView<sideNormalValueType,sideNormalProperties...> sideNormals,
924  const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
925  const ordinal_type worksetSideOrd,
926  const shards::CellTopology parentCell );
927 
966  template<typename faceNormalValueType, class ...faceNormalProperties,
967  typename worksetJacobianValueType, class ...worksetJacobianProperties>
968  static void
969  getPhysicalFaceNormals( Kokkos::DynRankView<faceNormalValueType,faceNormalProperties...> faceNormals,
970  const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
971  const ordinal_type worksetFaceOrd,
972  const shards::CellTopology parentCell );
973 
974  //============================================================================================//
975  // //
976  // Reference-to-physical frame mapping and its inverse //
977  // //
978  //============================================================================================//
979 
1016  template<typename physPointValueType, class ...physPointProperties,
1017  typename refPointValueType, class ...refPointProperties,
1018  typename worksetCellValueType, class ...worksetCellProperties,
1019  typename HGradBasisPtrType>
1020  static void
1021  mapToPhysicalFrame( Kokkos::DynRankView<physPointValueType,physPointProperties...> physPoints,
1022  const Kokkos::DynRankView<refPointValueType,refPointProperties...> refPoints,
1023  const Kokkos::DynRankView<worksetCellValueType,worksetCellProperties...> worksetCell,
1024  const HGradBasisPtrType basis );
1025 
1066  template<typename physPointValueType, class ...physPointProperties,
1067  typename refPointValueType, class ...refPointProperties,
1068  typename worksetCellValueType, class ...worksetCellProperties>
1069  static void
1070  mapToPhysicalFrame( Kokkos::DynRankView<physPointValueType,physPointProperties...> physPoints,
1071  const Kokkos::DynRankView<refPointValueType,refPointProperties...> refPoints,
1072  const Kokkos::DynRankView<worksetCellValueType,worksetCellProperties...> worksetCell,
1073  const shards::CellTopology cellTopo ) {
1074  auto basis = createHGradBasis<refPointValueType,refPointValueType>(cellTopo);
1075  mapToPhysicalFrame(physPoints,
1076  refPoints,
1077  worksetCell,
1078  basis);
1079  }
1080 
1131  template<typename refSubcellPointValueType, class ...refSubcellPointProperties,
1132  typename paramPointValueType, class ...paramPointProperties>
1133  static void
1134  mapToReferenceSubcell( Kokkos::DynRankView<refSubcellPointValueType,refSubcellPointProperties...> refSubcellPoints,
1135  const Kokkos::DynRankView<paramPointValueType,paramPointProperties...> paramPoints,
1136  const ordinal_type subcellDim,
1137  const ordinal_type subcellOrd,
1138  const shards::CellTopology parentCell );
1139 
1140 
1141  //============================================================================================//
1142  // //
1143  // Physical-to-reference frame mapping and its inverse //
1144  // //
1145  //============================================================================================//
1146 
1147 
1191  template<typename refPointValueType, class ...refPointProperties,
1192  typename physPointValueType, class ...physPointProperties,
1193  typename worksetCellValueType, class ...worksetCellProperties>
1194  static void
1195  mapToReferenceFrame( Kokkos::DynRankView<refPointValueType,refPointProperties...> refPoints,
1196  const Kokkos::DynRankView<physPointValueType,physPointProperties...> physPoints,
1197  const Kokkos::DynRankView<worksetCellValueType,worksetCellProperties...> worksetCell,
1198  const shards::CellTopology cellTopo );
1199 
1226  template<typename refPointValueType, class ...refPointProperties,
1227  typename initGuessValueType, class ...initGuessProperties,
1228  typename physPointValueType, class ...physPointProperties,
1229  typename worksetCellValueType, class ...worksetCellProperties,
1230  typename HGradBasisPtrType>
1231  static void
1232  mapToReferenceFrameInitGuess( Kokkos::DynRankView<refPointValueType,refPointProperties...> refPoints,
1233  const Kokkos::DynRankView<initGuessValueType,initGuessProperties...> initGuess,
1234  const Kokkos::DynRankView<physPointValueType,physPointProperties...> physPoints,
1235  const Kokkos::DynRankView<worksetCellValueType,worksetCellProperties...> worksetCell,
1236  const HGradBasisPtrType basis );
1237 
1238 
1269  template<typename refPointValueType, class ...refPointProperties,
1270  typename initGuessValueType, class ...initGuessProperties,
1271  typename physPointValueType, class ...physPointProperties,
1272  typename worksetCellValueType, class ...worksetCellProperties>
1273  static void
1274  mapToReferenceFrameInitGuess( Kokkos::DynRankView<refPointValueType,refPointProperties...> refPoints,
1275  const Kokkos::DynRankView<initGuessValueType,initGuessProperties...> initGuess,
1276  const Kokkos::DynRankView<physPointValueType,physPointProperties...> physPoints,
1277  const Kokkos::DynRankView<worksetCellValueType,worksetCellProperties...> worksetCell,
1278  const shards::CellTopology cellTopo ) {
1279  auto basis = createHGradBasis<refPointValueType,refPointValueType>(cellTopo);
1280  mapToReferenceFrameInitGuess(refPoints,
1281  initGuess,
1282  physPoints,
1283  worksetCell,
1284  basis);
1285  }
1286 
1287 
1288  //============================================================================================//
1289  // //
1290  // Control Volume Coordinates //
1291  // //
1292  //============================================================================================//
1293 
1367  template<typename subcvCoordValueType, class ...subcvCoordProperties,
1368  typename cellCoordValueType, class ...cellCoordProperties>
1369  static void
1370  getSubcvCoords( Kokkos::DynRankView<subcvCoordValueType,subcvCoordProperties...> subcvCoords,
1371  const Kokkos::DynRankView<cellCoordValueType,cellCoordProperties...> cellCoords,
1372  const shards::CellTopology primaryCell );
1373 
1374  //============================================================================================//
1375  // //
1376  // Inclusion tests //
1377  // //
1378  //============================================================================================//
1379 
1389  template<typename pointValueType, class ...pointProperties>
1390  static bool
1391  checkPointInclusion( const Kokkos::DynRankView<pointValueType,pointProperties...> point,
1392  const shards::CellTopology cellTopo,
1393  const double thres = threshold() );
1394 
1395  // template<class ArrayPoint>
1396  // static int checkPointsetInclusion(const ArrayPoint & points,
1397  // const shards::CellTopology & cellTopo,
1398  // const double & threshold = INTREPID2_THRESHOLD);
1399 
1400 
1401 
1428  template<typename inCellValueType, class ...inCellProperties,
1429  typename pointValueType, class ...pointProperties>
1430  static void checkPointwiseInclusion( Kokkos::DynRankView<inCellValueType,inCellProperties...> inCell,
1431  const Kokkos::DynRankView<pointValueType,pointProperties...> points,
1432  const shards::CellTopology cellTopo,
1433  const double thres = threshold() );
1434 
1456  template<typename inCellValueType, class ...inCellProperties,
1457  typename pointValueType, class ...pointProperties,
1458  typename cellWorksetValueType, class ...cellWorksetProperties>
1459  static void checkPointwiseInclusion( Kokkos::DynRankView<inCellValueType,inCellProperties...> inCell,
1460  const Kokkos::DynRankView<pointValueType,pointProperties...> points,
1461  const Kokkos::DynRankView<cellWorksetValueType,cellWorksetProperties...> cellWorkset,
1462  const shards::CellTopology cellTopo,
1463  const double thres = threshold() );
1464 
1465 
1466  // //============================================================================================//
1467  // // //
1468  // // Debug //
1469  // // //
1470  // //============================================================================================//
1471 
1472 
1473  // /** \brief Prints the reference space coordinates of the vertices of the specified subcell
1474  // \param subcellDim [in] - dimension of the subcell where points are mapped to
1475  // \param subcellOrd [in] - subcell ordinal
1476  // \param parentCell [in] - cell topology of the parent cell.
1477  // */
1478  // static void printSubcellVertices(const int subcellDim,
1479  // const int subcellOrd,
1480  // const shards::CellTopology & parentCell);
1481 
1482 
1483 
1484  // /** \brief Prints the nodes of a subcell from a cell workset
1485 
1486  // */
1487  // template<class ArrayCell>
1488  // static void printWorksetSubcell(const ArrayCell & cellWorkset,
1489  // const shards::CellTopology & parentCell,
1490  // const int& pCellOrd,
1491  // const int& subcellDim,
1492  // const int& subcellOrd,
1493  // const int& fieldWidth = 3);
1494  };
1495 
1496  //============================================================================================//
1497  // //
1498  // Validation of input/output arguments for CellTools methods //
1499  // //
1500  //============================================================================================//
1501 
1508  template<typename jacobianViewType,
1509  typename pointViewType,
1510  typename worksetCellViewType>
1511  static void
1512  CellTools_setJacobianArgs( const jacobianViewType jacobian,
1513  const pointViewType points,
1514  const worksetCellViewType worksetCell,
1515  const shards::CellTopology cellTopo );
1516 
1521  template<typename jacobianInvViewType,
1522  typename jacobianViewType>
1523  static void
1524  CellTools_setJacobianInvArgs( const jacobianInvViewType jacobianInv,
1525  const jacobianViewType jacobian );
1526 
1527 
1532  template<typename jacobianDetViewType,
1533  typename jacobianViewType>
1534  static void
1535  CellTools_setJacobianDetArgs( const jacobianDetViewType jacobianDet,
1536  const jacobianViewType jacobian );
1537 
1538 
1545  template<typename physPointViewType,
1546  typename refPointViewType,
1547  typename worksetCellViewType>
1548  static void
1549  CellTools_mapToPhysicalFrameArgs( const physPointViewType physPoints,
1550  const refPointViewType refPoints,
1551  const worksetCellViewType worksetCell,
1552  const shards::CellTopology cellTopo );
1553 
1554 
1561  template<typename refPointViewType,
1562  typename physPointViewType,
1563  typename worksetCellViewType>
1564  static void
1565  CellTools_mapToReferenceFrameArgs( const refPointViewType refPoints,
1566  const physPointViewType physPoints,
1567  const worksetCellViewType worksetCell,
1568  const shards::CellTopology cellTopo );
1569 
1570 
1571 
1579  template<typename refPointViewType,
1580  typename initGuessViewType,
1581  typename physPointViewType,
1582  typename worksetCellViewType>
1583  static void
1584  CellTools_mapToReferenceFrameInitGuess( const refPointViewType refPoints,
1585  const initGuessViewType initGuess,
1586  const physPointViewType physPoints,
1587  const worksetCellViewType worksetCell,
1588  const shards::CellTopology cellTopo );
1589 
1590  // /** \brief Validates arguments to Intrepid2::CellTools::checkPointwiseInclusion
1591  // \param inCell [out] - rank-1 (P) array required
1592  // \param physPoints [in] - rank-2 (P,D) array required
1593  // \param cellWorkset [in] - rank-3 (C,N,D) array required
1594  // \param whichCell [in] - 0 <= whichCell < C required
1595  // \param cellTopo [in] - cell topology with a reference cell required
1596  // */
1597  // template<class ArrayIncl, class ArrayPoint, class ArrayCell>
1598  // static void
1599  // validateArguments_checkPointwiseInclusion(ArrayIncl & inCell,
1600  // const ArrayPoint & physPoints,
1601  // const ArrayCell & cellWorkset,
1602  // const int & whichCell,
1603  // const shards::CellTopology & cell);
1604 
1605 
1606 
1607 }
1608 
1610 
1613 
1616 
1619 
1621 
1622 // not yet converted ...
1624 
1625 // #include "Intrepid2_CellToolsDefDebug.hpp"
1626 
1627 
1628 #endif
1629 
static void getReferenceEdgeTangent(Kokkos::DynRankView< refEdgeTangentValueType, refEdgeTangentProperties...> refEdgeTangent, const ordinal_type edgeOrd, const shards::CellTopology parentCell)
Computes constant tangent vectors to edges of 2D or 3D reference cells.
Implementation of the default H(grad)-compatible FEM basis of degree 2 on Wedge cell.
Implementation of the default H(grad)-compatible FEM basis of degree 1 on Quadrilateral cell...
static void getReferenceSubcellNodes(Kokkos::DynRankView< subcellNodeValueType, subcellNodeProperties...> subcellNodes, const ordinal_type subcellDim, const ordinal_type subcellOrd, const shards::CellTopology parentCell)
Retrieves the Cartesian coordinates of all nodes of a reference subcell.
Definition file for the reference to physical mappings in the Intrepid2::CellTools class...
static bool checkPointInclusion(const Kokkos::DynRankView< pointValueType, pointProperties...> point, const shards::CellTopology cellTopo, const double thres=threshold())
Checks if a point belongs to a reference cell.
Header file for the Intrepid2::Basis_HGRAD_TRI_C1_FEM class.
Implementation of the default H(grad)-compatible FEM basis of degree 1 on Wedge cell.
Header file for the Intrepid2::Basis_HGRAD_QUAD_C1_FEM class.
static void mapToPhysicalFrame(Kokkos::DynRankView< physPointValueType, physPointProperties...> physPoints, const Kokkos::DynRankView< refPointValueType, refPointProperties...> refPoints, const Kokkos::DynRankView< worksetCellValueType, worksetCellProperties...> worksetCell, const HGradBasisPtrType basis)
Computes F, the reference-to-physical frame map.
Header file for the Intrepid2::Basis_HGRAD_QUAD_C2_FEM class.
Definition file for the control volume functions of the Intrepid2::CellTools class.
Implementation of the default H(grad)-compatible FEM basis of degree 2 on Hexahedron cell...
static void getSubcvCoords(Kokkos::DynRankView< subcvCoordValueType, subcvCoordProperties...> subcvCoords, const Kokkos::DynRankView< cellCoordValueType, cellCoordProperties...> cellCoords, const shards::CellTopology primaryCell)
Computes coordinates of sub-control volumes in each primary cell.
Implementation of the default H(grad)-compatible FEM basis of degree 1 on Triangle cell...
static void getReferenceVertex(Kokkos::DynRankView< cellVertexValueType, cellVertexProperties...> cellVertex, const shards::CellTopology cell, const ordinal_type vertexOrd)
Retrieves the Cartesian coordinates of a reference cell vertex.
Implementation of the default H(grad)-compatible FEM basis of degree 1 on Pyramid cell...
CellTools()=default
Default constructor.
static void getReferenceSideNormal(Kokkos::DynRankView< refSideNormalValueType, refSideNormalProperties...> refSideNormal, const ordinal_type sideOrd, const shards::CellTopology parentCell)
Computes constant normal vectors to sides of 2D or 3D reference cells.
static void setJacobian(Kokkos::DynRankView< jacobianValueType, jacobianProperties...> jacobian, const Kokkos::DynRankView< pointValueType, pointProperties...> points, const Kokkos::DynRankView< worksetCellValueType, worksetCellProperties...> worksetCell, const shards::CellTopology cellTopo)
Computes the Jacobian matrix DF of the reference-to-physical frame map F.
Header file for the Intrepid2::Basis_HGRAD_HEX_C2_FEM class.
Header file for the Intrepid2::Basis_HGRAD_PYR_C1_FEM class.
static bool hasReferenceCell(const shards::CellTopology cellTopo)
Checks if a cell topology has reference cell.
static void setJacobianDet(Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties...> jacobianDet, const Kokkos::DynRankView< jacobianValueType, jacobianProperties...> jacobian)
Computes the determinant of the Jacobian matrix DF of the reference-to-physical frame map F...
Reference node data for each supported topology.
A stateless class for operations on cell data. Provides methods for:
Implementation of the default H(grad)-compatible FEM basis of degree 2 on Triangle cell...
Header file for the Intrepid2::Basis_HGRAD_TRI_C2_FEM class.
static void setJacobianInv(Kokkos::DynRankView< jacobianInvValueType, jacobianInvProperties...> jacobianInv, const Kokkos::DynRankView< jacobianValueType, jacobianProperties...> jacobian)
Computes the inverse of the Jacobian matrix DF of the reference-to-physical frame map F...
static void mapToPhysicalFrame(Kokkos::DynRankView< physPointValueType, physPointProperties...> physPoints, const Kokkos::DynRankView< refPointValueType, refPointProperties...> refPoints, const Kokkos::DynRankView< worksetCellValueType, worksetCellProperties...> worksetCell, const shards::CellTopology cellTopo)
Computes F, the reference-to-physical frame map.
Header file for the Intrepid2::Basis_HGRAD_WEDGE_C2_FEM class.
Definition file for point inclusion functions of the Intrepid2::CellTools class.
Header function for Intrepid2::Util class and other utility functions.
Header file for the Intrepid2::Basis_HGRAD_WEDGE_C1_FEM class.
Header file for the Intrepid2::Basis_HGRAD_TET_C2_FEM class.
static void mapToReferenceSubcell(Kokkos::DynRankView< refSubcellPointValueType, refSubcellPointProperties...> refSubcellPoints, const Kokkos::DynRankView< paramPointValueType, paramPointProperties...> paramPoints, const ordinal_type subcellDim, const ordinal_type subcellOrd, const shards::CellTopology parentCell)
Computes parameterization maps of 1- and 2-subcells of reference cells.
static void mapToReferenceFrameInitGuess(Kokkos::DynRankView< refPointValueType, refPointProperties...> refPoints, const Kokkos::DynRankView< initGuessValueType, initGuessProperties...> initGuess, const Kokkos::DynRankView< physPointValueType, physPointProperties...> physPoints, const Kokkos::DynRankView< worksetCellValueType, worksetCellProperties...> worksetCell, const HGradBasisPtrType basis)
Computation of , the inverse of the reference-to-physical frame map using user-supplied initial guess...
Header file for the Intrepid2::Basis_HGRAD_HEX_C1_FEM class.
static void getPhysicalFaceTangents(Kokkos::DynRankView< faceTanUValueType, faceTanUProperties...> faceTanU, Kokkos::DynRankView< faceTanVValueType, faceTanVProperties...> faceTanV, const Kokkos::DynRankView< worksetJacobianValueType, worksetJacobianProperties...> worksetJacobians, const ordinal_type worksetFaceOrd, const shards::CellTopology parentCell)
Computes non-normalized tangent vector pairs to physical faces in a face workset ; (see Subcell works...
Definition file for the validate arguments functions of the Intrepid2::CellTools class.
Definition file for the physical to reference mappings in the Intrepid2::CellTools class...
static void getPhysicalEdgeTangents(Kokkos::DynRankView< edgeTangentValueType, edgeTangentProperties...> edgeTangents, const Kokkos::DynRankView< worksetJacobianValueType, worksetJacobianProperties...> worksetJacobians, const ordinal_type worksetEdgeOrd, const shards::CellTopology parentCell)
Computes non-normalized tangent vectors to physical edges in an edge workset ; (see Subcell worksets ...
static void setJacobian(Kokkos::DynRankView< jacobianValueType, jacobianProperties...> jacobian, const Kokkos::DynRankView< pointValueType, pointProperties...> points, const Kokkos::DynRankView< worksetCellValueType, worksetCellProperties...> worksetCell, const HGradBasisPtrType basis)
Computes the Jacobian matrix DF of the reference-to-physical frame map F.
Definition file for the parameterization functions of the Intrepid2::CellTools class.
Definition file for the Jacobian functions in the Intrepid2::CellTools class.
Definition file for node data and subcell functions of the Intrepid2::CellTools class.
Header file for the Intrepid2::Basis_HGRAD_TET_C1_FEM class.
Contains definitions of custom data types in Intrepid2.
static void mapToReferenceFrame(Kokkos::DynRankView< refPointValueType, refPointProperties...> refPoints, const Kokkos::DynRankView< physPointValueType, physPointProperties...> physPoints, const Kokkos::DynRankView< worksetCellValueType, worksetCellProperties...> worksetCell, const shards::CellTopology cellTopo)
Computes , the inverse of the reference-to-physical frame map using a default initial guess...
static void setReferenceNodeData()
Set reference node coordinates for supported topologies.
Header file with additional documentation for the Intrepid2::CellTools class.
Implementation of the default H(grad)-compatible FEM basis of degree 2 on Quadrilateral cell...
Header file for the Intrepid2::Basis_HGRAD_LINE_C1_FEM class.
static void checkPointwiseInclusion(Kokkos::DynRankView< inCellValueType, inCellProperties...> inCell, const Kokkos::DynRankView< pointValueType, pointProperties...> points, const shards::CellTopology cellTopo, const double thres=threshold())
Checks every point in a set for inclusion in a reference cell.
Implementation of the default H(grad)-compatible FEM basis of degree 2 on Tetrahedron cell...
Implementation of the default H(grad)-compatible FEM basis of degree 2 on Tetrahedron cell...
static Teuchos::RCP< Basis< ExecSpaceType, outputValueType, pointValueType > > createHGradBasis(const shards::CellTopology cellTopo)
Generates default HGrad basis based on cell topology.
Header file for the Intrepid2::Basis_HGRAD_TET_COMP12_FEM class.
static void getReferenceFaceTangents(Kokkos::DynRankView< refFaceTanUValueType, refFaceTanUProperties...> refFaceTanU, Kokkos::DynRankView< refFaceTanVValueType, refFaceTanVProperties...> refFaceTanV, const ordinal_type faceOrd, const shards::CellTopology parentCell)
Computes pairs of constant tangent vectors to faces of a 3D reference cells.
Implementation of the default H(grad)-compatible FEM basis of degree 1 on Line cell.
Parametrization coefficients of edges and faces of reference cells.
static void getPhysicalSideNormals(Kokkos::DynRankView< sideNormalValueType, sideNormalProperties...> sideNormals, const Kokkos::DynRankView< worksetJacobianValueType, worksetJacobianProperties...> worksetJacobians, const ordinal_type worksetSideOrd, const shards::CellTopology parentCell)
Computes non-normalized normal vectors to physical sides in a side workset .
Implementation of the default H(grad)-compatible FEM basis of degree 1 on Tetrahedron cell...
~CellTools()=default
Destructor.
Header file for Intrepid2::RealSpaceTools class providing basic linear algebra functionality in 1D...
static void setSubcellParametrization()
Defines orientation-preserving parametrizations of reference edges and faces of cell topologies with ...
static void mapToReferenceFrameInitGuess(Kokkos::DynRankView< refPointValueType, refPointProperties...> refPoints, const Kokkos::DynRankView< initGuessValueType, initGuessProperties...> initGuess, const Kokkos::DynRankView< physPointValueType, physPointProperties...> physPoints, const Kokkos::DynRankView< worksetCellValueType, worksetCellProperties...> worksetCell, const shards::CellTopology cellTopo)
Computation of , the inverse of the reference-to-physical frame map using user-supplied initial guess...
static void getSubcellParametrization(subcellParamViewType &subcellParam, const ordinal_type subcellDim, const shards::CellTopology parentCell)
Returns array with the coefficients of the parametrization maps for the edges or faces of a reference...
static void getReferenceCellCenter(Kokkos::DynRankView< cellCenterValueType, cellCenterProperties...> cellCenter, Kokkos::DynRankView< cellVertexValueType, cellVertexProperties...> cellVertex, const shards::CellTopology cell)
Computes the Cartesian coordinates of reference cell center.
Implementation of the default H(grad)-compatible FEM basis of degree 1 on Hexahedron cell...
static void getPhysicalFaceNormals(Kokkos::DynRankView< faceNormalValueType, faceNormalProperties...> faceNormals, const Kokkos::DynRankView< worksetJacobianValueType, worksetJacobianProperties...> worksetJacobians, const ordinal_type worksetFaceOrd, const shards::CellTopology parentCell)
Computes non-normalized normal vectors to physical faces in a face workset ; (see Subcell worksets fo...
Reference node containers for each supported topology.
static void getReferenceFaceNormal(Kokkos::DynRankView< refFaceNormalValueType, refFaceNormalProperties...> refFaceNormal, const ordinal_type faceOrd, const shards::CellTopology parentCell)
Computes constant normal vectors to faces of 3D reference cell.
Header file for the abstract base class Intrepid2::Basis.
static void getReferenceNode(Kokkos::DynRankView< cellNodeValueType, cellNodeProperties...> cellNode, const shards::CellTopology cell, const ordinal_type nodeOrd)
Retrieves the Cartesian coordinates of a reference cell node.
static void getReferenceSubcellVertices(Kokkos::DynRankView< subcellVertexValueType, subcellVertexProperties...> subcellVertices, const ordinal_type subcellDim, const ordinal_type subcellOrd, const shards::CellTopology parentCell)
Retrieves the Cartesian coordinates of all vertices of a reference subcell.