Intrepid
Intrepid_CellTools.hpp
Go to the documentation of this file.
1 #ifndef INTREPID_CELTOOLS_HPP
2 #define INTREPID_CELTOOLS_HPP
3 
4 // @HEADER
5 // ************************************************************************
6 //
7 // Intrepid Package
8 // Copyright (2007) Sandia Corporation
9 //
10 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
11 // license for use of this work by or on behalf of the U.S. Government.
12 //
13 // Redistribution and use in source and binary forms, with or without
14 // modification, are permitted provided that the following conditions are
15 // met:
16 //
17 // 1. Redistributions of source code must retain the above copyright
18 // notice, this list of conditions and the following disclaimer.
19 //
20 // 2. Redistributions in binary form must reproduce the above copyright
21 // notice, this list of conditions and the following disclaimer in the
22 // documentation and/or other materials provided with the distribution.
23 //
24 // 3. Neither the name of the Corporation nor the names of the
25 // contributors may be used to endorse or promote products derived from
26 // this software without specific prior written permission.
27 //
28 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
29 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
30 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
31 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
32 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
33 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
34 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
35 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
36 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
37 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
38 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
39 //
40 // Questions? Contact Pavel Bochev (pbboche@sandia.gov)
41 // Denis Ridzal (dridzal@sandia.gov), or
42 // Kara Peterson (kjpeter@sandia.gov)
43 //
44 // ************************************************************************
45 // @HEADER
46 
52 #ifndef INTREPID_CELLTOOLS_HPP
53 #define INTREPID_CELLTOOLS_HPP
54 
55 
58 #include "Intrepid_ConfigDefs.hpp"
59 #include "Intrepid_Types.hpp"
60 #include "Intrepid_Utils.hpp"
61 #include "Intrepid_Basis.hpp"
63 #include "Intrepid_HGRAD_QUAD_C1_FEM.hpp"
64 #include "Intrepid_HGRAD_TET_C1_FEM.hpp"
65 #include "Intrepid_HGRAD_WEDGE_C1_FEM.hpp"
66 #include "Intrepid_HGRAD_PYR_C1_FEM.hpp"
67 #include "Intrepid_HGRAD_HEX_C1_FEM.hpp"
68 
70 
72 #include "Intrepid_HGRAD_QUAD_C2_FEM.hpp"
80 
81 #include "Shards_CellTopology.hpp"
82 #include "Shards_BasicTopologies.hpp"
83 
84 #include "Teuchos_Assert.hpp"
85 #include "Teuchos_RCP.hpp"
86 
87 #include <Intrepid_Rank.hpp>
88 
89 namespace Intrepid {
90 
91 //nn
92  //============================================================================================//
93  // //
94  // CellTools //
95  // //
96  //============================================================================================//
97 
108 template<class Scalar>
109 class CellTools {
110 private:
111 
112  //============================================================================================//
113  // //
114  // Parametrization coefficients of edges and faces of reference cells //
115  // //
116  //============================================================================================//
117 
118 
131  static const FieldContainer<double>& getSubcellParametrization(const int subcellDim,
132  const shards::CellTopology& parentCell);
133 
134 
135 
175  static void setSubcellParametrization(FieldContainer<double>& subcellParametrization,
176  const int subcellDim,
177  const shards::CellTopology& parentCell);
178 
179  //============================================================================================//
180  // //
181  // Validation of input/output arguments for CellTools methods //
182  // //
183  //============================================================================================//
184 
192  template<class ArrayJac, class ArrayPoint, class ArrayCell>
193  static void validateArguments_setJacobian(const ArrayJac & jacobian,
194  const ArrayPoint & points,
195  const ArrayCell & cellWorkset,
196  const int & whichCell,
197  const shards::CellTopology & cellTopo);
198 
199 
200 
205  template<class ArrayJacInv, class ArrayJac>
206  static void validateArguments_setJacobianInv(const ArrayJacInv & jacobianInv,
207  const ArrayJac & jacobian);
208 
209 
210 
215  template<class ArrayJacDet, class ArrayJac>
216  static void validateArguments_setJacobianDetArgs(const ArrayJacDet & jacobianDet,
217  const ArrayJac & jacobian);
218 
219 
220 
228  template<class ArrayPhysPoint, class ArrayRefPoint, class ArrayCell>
229  static void validateArguments_mapToPhysicalFrame(const ArrayPhysPoint & physPoints,
230  const ArrayRefPoint & refPoints,
231  const ArrayCell & cellWorkset,
232  const shards::CellTopology & cellTopo,
233  const int& whichCell);
234 
235 
236 
244  template<class ArrayRefPoint, class ArrayPhysPoint, class ArrayCell>
245  static void validateArguments_mapToReferenceFrame(const ArrayRefPoint & refPoints,
246  const ArrayPhysPoint & physPoints,
247  const ArrayCell & cellWorkset,
248  const shards::CellTopology & cellTopo,
249  const int& whichCell);
250 
251 
252 
261  template<class ArrayRefPoint, class ArrayInitGuess, class ArrayPhysPoint, class ArrayCell>
262  static void validateArguments_mapToReferenceFrame(const ArrayRefPoint & refPoints,
263  const ArrayInitGuess & initGuess,
264  const ArrayPhysPoint & physPoints,
265  const ArrayCell & cellWorkset,
266  const shards::CellTopology & cellTopo,
267  const int& whichCell);
268 
269 
270 
278  template<class ArrayIncl, class ArrayPoint, class ArrayCell>
279  static void validateArguments_checkPointwiseInclusion(ArrayIncl & inCell,
280  const ArrayPoint & physPoints,
281  const ArrayCell & cellWorkset,
282  const int & whichCell,
283  const shards::CellTopology & cell);
284 public:
285 
288  CellTools(){ };
289 
290 
294 
295  //============================================================================================//
296  // //
297  // Jacobian, inverse Jacobian and Jacobian determinant //
298  // //
299  //============================================================================================//
300 
349  template<class ArrayJac, class ArrayPoint, class ArrayCell, bool typecheck>
351 
352 
353  /*template<class ArrayJac, class ArrayPoint, class ArrayCell>
354  static void setJacobianTemp(ArrayJac & jacobian,
355  const ArrayPoint & points,
356  const ArrayCell & cellWorkset,
357  const shards::CellTopology & cellTopo,
358  const int & whichCell = -1);*/
359  template<class ArrayJac, class ArrayPoint, class ArrayCell>
360  static void setJacobian(ArrayJac & jacobian,
361  const ArrayPoint & points,
362  const ArrayCell & cellWorkset,
363  const shards::CellTopology & cellTopo,
364  const int & whichCell = -1);
365 
366  template<class ArrayJac, class ArrayPoint, class ArrayCell>
367  static void setJacobian(ArrayJac & jacobian,
368  const ArrayPoint & points,
369  const ArrayCell & cellWorkset,
370  const Teuchos::RCP< Basis< Scalar, FieldContainer<Scalar> > > HGRAD_Basis,
371  const int & whichCell = -1);
372 
373 
386  template<class ArrayJacInv, class ArrayJac>
387  static void setJacobianInv(ArrayJacInv & jacobianInv,
388  const ArrayJac & jacobian);
389 
390  /* template<class ArrayJacInv, class ArrayJac>
391  static void setJacobianInvTemp(ArrayJacInv & jacobianInv,
392  const ArrayJac & jacobian);
393  */
394 
407  template<class ArrayJacDet, class ArrayJac>
408  static void setJacobianDet(ArrayJacDet & jacobianDet,
409  const ArrayJac & jacobian);
410 
411  /* template<class ArrayJacDet, class ArrayJac>
412  static void setJacobianDetTemp(ArrayJacDet & jacobianDet,
413  const ArrayJac & jacobian);*/
414 
415  //============================================================================================//
416  // //
417  // Reference-to-physical frame mapping and its inverse //
418  // //
419  //============================================================================================//
420 
476  template<class ArrayPhysPoint, class ArrayRefPoint, class ArrayCell>
477  static void mapToPhysicalFrame(ArrayPhysPoint & physPoints,
478  const ArrayRefPoint & refPoints,
479  const ArrayCell & cellWorkset,
480  const shards::CellTopology & cellTopo,
481  const int & whichCell = -1);
482 
483  template<class ArrayPhysPoint, class ArrayRefPoint, class ArrayCell>
484  static void mapToPhysicalFrame(ArrayPhysPoint & physPoints,
485  const ArrayRefPoint & refPoints,
486  const ArrayCell & cellWorkset,
487  const Teuchos::RCP< Basis< Scalar, FieldContainer<Scalar> > > HGRAD_Basis,
488  const int & whichCell = -1);
489 
545  /* template<class ArrayPhysPoint, class ArrayRefPoint, class ArrayCell>
546  static void mapToPhysicalFrameTemp(ArrayPhysPoint & physPoints,
547  const ArrayRefPoint & refPoints,
548  const ArrayCell & cellWorkset,
549  const shards::CellTopology & cellTopo,
550  const int & whichCell = -1);
551  */
552  template<class ArrayPhysPoint, class ArrayRefPoint, class ArrayCell, int refRank,int phyptsrank>
554 
613  template<class ArrayRefPoint, class ArrayPhysPoint, class ArrayCell>
614  static void mapToReferenceFrame(ArrayRefPoint & refPoints,
615  const ArrayPhysPoint & physPoints,
616  const ArrayCell & cellWorkset,
617  const shards::CellTopology & cellTopo,
618  const int & whichCell = -1);
619 
620 
621  template<class ArrayRefPoint, class ArrayPhysPoint, class ArrayCell>
622  static void mapToReferenceFrame(ArrayRefPoint & refPoints,
623  const ArrayPhysPoint & physPoints,
624  const ArrayCell & cellWorkset,
625  const Teuchos::RCP< Basis< Scalar, FieldContainer<Scalar> > > HGRAD_Basis,
626  const int & whichCell = -1);
627 
628 
674  template<class ArrayRefPoint, class ArrayInitGuess, class ArrayPhysPoint, class ArrayCell>
675  static void mapToReferenceFrameInitGuess(ArrayRefPoint & refPoints,
676  const ArrayInitGuess & initGuess,
677  const ArrayPhysPoint & physPoints,
678  const ArrayCell & cellWorkset,
679  const shards::CellTopology & cellTopo,
680  const int & whichCell = -1);
681 
682  template<class ArrayRefPoint, class ArrayInitGuess, class ArrayPhysPoint, class ArrayCell>
683  static void mapToReferenceFrameInitGuess(ArrayRefPoint & refPoints,
684  const ArrayInitGuess & initGuess,
685  const ArrayPhysPoint & physPoints,
686  const ArrayCell & cellWorkset,
687  const Teuchos::RCP< Basis< Scalar, FieldContainer<Scalar> > > HGRAD_Basis,
688  const int & whichCell = -1);
689 
690 
691 
742  template<class ArraySubcellPoint, class ArrayParamPoint>
743  static void mapToReferenceSubcell(ArraySubcellPoint & refSubcellPoints,
744  const ArrayParamPoint & paramPoints,
745  const int subcellDim,
746  const int subcellOrd,
747  const shards::CellTopology & parentCell);
748 
749 
750 
776  template<class ArrayEdgeTangent>
777  static void getReferenceEdgeTangent(ArrayEdgeTangent & refEdgeTangent,
778  const int & edgeOrd,
779  const shards::CellTopology & parentCell);
780 
781 
782 
819  template<class ArrayFaceTangentU, class ArrayFaceTangentV>
820  static void getReferenceFaceTangents(ArrayFaceTangentU & refFaceTanU,
821  ArrayFaceTangentV & refFaceTanV,
822  const int & faceOrd,
823  const shards::CellTopology & parentCell);
824 
825 
826 
889  template<class ArraySideNormal>
890  static void getReferenceSideNormal(ArraySideNormal & refSideNormal,
891  const int & sideOrd,
892  const shards::CellTopology & parentCell);
893 
894 
895 
934  template<class ArrayFaceNormal>
935  static void getReferenceFaceNormal(ArrayFaceNormal & refFaceNormal,
936  const int & faceOrd,
937  const shards::CellTopology & parentCell);
938 
939 
940 
970  template<class ArrayEdgeTangent, class ArrayJac>
971  static void getPhysicalEdgeTangents(ArrayEdgeTangent & edgeTangents,
972  const ArrayJac & worksetJacobians,
973  const int & worksetEdgeOrd,
974  const shards::CellTopology & parentCell);
975 
976  /* template<class ArrayEdgeTangent, class ArrayJac>
977  static void getPhysicalEdgeTangentsTemp(ArrayEdgeTangent & edgeTangents,
978  const ArrayJac & worksetJacobians,
979  const int & worksetEdgeOrd,
980  const shards::CellTopology & parentCell); */
981 
1021  template<class ArrayFaceTangentU, class ArrayFaceTangentV, class ArrayJac>
1022  static void getPhysicalFaceTangents(ArrayFaceTangentU & faceTanU,
1023  ArrayFaceTangentV & faceTanV,
1024  const ArrayJac & worksetJacobians,
1025  const int & worksetFaceOrd,
1026  const shards::CellTopology & parentCell);
1027 
1028  /* template<class ArrayFaceTangentU, class ArrayFaceTangentV, class ArrayJac>
1029  static void getPhysicalFaceTangentsTemp(ArrayFaceTangentU & faceTanU,
1030  ArrayFaceTangentV & faceTanV,
1031  const ArrayJac & worksetJacobians,
1032  const int & worksetFaceOrd,
1033  const shards::CellTopology & parentCell);
1034  */
1095  template<class ArraySideNormal, class ArrayJac>
1096  static void getPhysicalSideNormals(ArraySideNormal & sideNormals,
1097  const ArrayJac & worksetJacobians,
1098  const int & worksetSideOrd,
1099  const shards::CellTopology & parentCell);
1100 
1101 
1102 
1141  template<class ArrayFaceNormal, class ArrayJac>
1142  static void getPhysicalFaceNormals(ArrayFaceNormal & faceNormals,
1143  const ArrayJac & worksetJacobians,
1144  const int & worksetFaceOrd,
1145  const shards::CellTopology & parentCell);
1146 
1147  /* template<class ArrayFaceNormal, class ArrayJac>
1148  static void getPhysicalFaceNormalsTemp(ArrayFaceNormal & faceNormals,
1149  const ArrayJac & worksetJacobians,
1150  const int & worksetFaceOrd,
1151  const shards::CellTopology & parentCell);
1152  */
1153  //============================================================================================//
1154  // //
1155  // Inclusion tests //
1156  // //
1157  //============================================================================================//
1158 
1169  static int checkPointInclusion(const Scalar* point,
1170  const int pointDim,
1171  const shards::CellTopology & cellTopo,
1172  const double & threshold = INTREPID_THRESHOLD);
1173 
1174 
1175 
1188  template<class ArrayPoint>
1189  static int checkPointsetInclusion(const ArrayPoint & points,
1190  const shards::CellTopology & cellTopo,
1191  const double & threshold = INTREPID_THRESHOLD);
1192 
1193 
1194 
1221  template<class ArrayIncl, class ArrayPoint>
1222  static void checkPointwiseInclusion(ArrayIncl & inRefCell,
1223  const ArrayPoint & points,
1224  const shards::CellTopology & cellTopo,
1225  const double & threshold = INTREPID_THRESHOLD);
1226 
1227 
1228 
1264  template<class ArrayIncl, class ArrayPoint, class ArrayCell>
1265  static void checkPointwiseInclusion(ArrayIncl & inCell,
1266  const ArrayPoint & points,
1267  const ArrayCell & cellWorkset,
1268  const shards::CellTopology & cell,
1269  const int & whichCell = -1,
1270  const double & threshold = INTREPID_THRESHOLD);
1271 
1272 
1273 
1284  static const double* getReferenceVertex(const shards::CellTopology& cell,
1285  const int vertexOrd);
1286 
1287 
1288 
1303  template<class ArraySubcellVert>
1304  static void getReferenceSubcellVertices(ArraySubcellVert & subcellVertices,
1305  const int subcellDim,
1306  const int subcellOrd,
1307  const shards::CellTopology& parentCell);
1308 
1309 
1310 
1326  static const double* getReferenceNode(const shards::CellTopology& cell,
1327  const int nodeOrd);
1328 
1329 
1330 
1344  template<class ArraySubcellNode>
1345  static void getReferenceSubcellNodes(ArraySubcellNode& subcellNodes,
1346  const int subcellDim,
1347  const int subcellOrd,
1348  const shards::CellTopology& parentCell);
1349 
1350 
1351 
1357  static int hasReferenceCell(const shards::CellTopology & cellTopo);
1358 
1359 
1360 
1361  //============================================================================================//
1362  // //
1363  // Debug //
1364  // //
1365  //============================================================================================//
1366 
1367 
1373  static void printSubcellVertices(const int subcellDim,
1374  const int subcellOrd,
1375  const shards::CellTopology & parentCell);
1376 
1377 
1378 
1382  template<class ArrayCell>
1383  static void printWorksetSubcell(const ArrayCell & cellWorkset,
1384  const shards::CellTopology & parentCell,
1385  const int& pCellOrd,
1386  const int& subcellDim,
1387  const int& subcellOrd,
1388  const int& fieldWidth = 3);
1389 
1390  //============================================================================================//
1391  // //
1392  // Control Volume Coordinates //
1393  // //
1394  //============================================================================================//
1395 
1469  template<class ArrayCVCoord, class ArrayCellCoord>
1470  static void getSubCVCoords(ArrayCVCoord & subCVcoords, const ArrayCellCoord & cellCoords,
1471  const shards::CellTopology& primaryCell);
1472 
1479  template<class ArrayCent, class ArrayCellCoord>
1480  static void getBarycenter(ArrayCent & barycenter, const ArrayCellCoord & cellCoords);
1481 
1482 
1483  }; // class CellTools
1484 
1485 } // namespace Intrepid
1486 
1487 // include templated function definitions
1488 #include "Intrepid_CellToolsDef.hpp"
1489 
1490 #endif
1491 
1492 /***************************************************************************************************
1493  ** **
1494  ** D O C U M E N T A T I O N P A G E S **
1495  ** **
1496  **************************************************************************************************/
1497 
1778 #endif
1779 
1780 #if defined(Intrepid_SHOW_DEPRECATED_WARNINGS)
1781 #ifdef __GNUC__
1782 #warning "The Intrepid package is deprecated"
1783 #endif
1784 #endif
1785 
static const double * getReferenceNode(const shards::CellTopology &cell, const int nodeOrd)
Retrieves the Cartesian coordinates of a reference cell node.
Header file for the Intrepid::HGRAD_WEDGE_I2_FEM class.
static void validateArguments_mapToReferenceFrame(const ArrayRefPoint &refPoints, const ArrayPhysPoint &physPoints, const ArrayCell &cellWorkset, const shards::CellTopology &cellTopo, const int &whichCell)
Validates arguments to Intrepid::CellTools::mapToReferenceFrame with default initial guess...
static void getReferenceEdgeTangent(ArrayEdgeTangent &refEdgeTangent, const int &edgeOrd, const shards::CellTopology &parentCell)
Computes constant tangent vectors to edges of 2D or 3D reference cells.
static void mapToReferenceSubcell(ArraySubcellPoint &refSubcellPoints, const ArrayParamPoint &paramPoints, const int subcellDim, const int subcellOrd, const shards::CellTopology &parentCell)
Computes parameterization maps of 1- and 2-subcells of reference cells.
Definition file for the Intrepid::CellTools class.
static void getPhysicalFaceNormals(ArrayFaceNormal &faceNormals, const ArrayJac &worksetJacobians, const int &worksetFaceOrd, const shards::CellTopology &parentCell)
Computes non-normalized normal vectors to physical faces in a face workset ; (see Subcell worksets fo...
static void setJacobianInv(ArrayJacInv &jacobianInv, const ArrayJac &jacobian)
Computes the inverse of the Jacobian matrix DF of the reference-to-physical frame map F...
static void getPhysicalEdgeTangents(ArrayEdgeTangent &edgeTangents, const ArrayJac &worksetJacobians, const int &worksetEdgeOrd, const shards::CellTopology &parentCell)
Computes non-normalized tangent vectors to physical edges in an edge workset ; (see Subcell worksets ...
static void mapToReferenceFrame(ArrayRefPoint &refPoints, const ArrayPhysPoint &physPoints, const ArrayCell &cellWorkset, const shards::CellTopology &cellTopo, const int &whichCell=-1)
Computes , the inverse of the reference-to-physical frame map using a default initial guess...
static void getReferenceSideNormal(ArraySideNormal &refSideNormal, const int &sideOrd, const shards::CellTopology &parentCell)
Computes constant normal vectors to sides of 2D or 3D reference cells.
static void getReferenceFaceTangents(ArrayFaceTangentU &refFaceTanU, ArrayFaceTangentV &refFaceTanV, const int &faceOrd, const shards::CellTopology &parentCell)
Computes pairs of constant tangent vectors to faces of a 3D reference cells.
static void printSubcellVertices(const int subcellDim, const int subcellOrd, const shards::CellTopology &parentCell)
Prints the reference space coordinates of the vertices of the specified subcell.
Header file for the Intrepid::HGRAD_TRI_C1_FEM class.
static void getReferenceFaceNormal(ArrayFaceNormal &refFaceNormal, const int &faceOrd, const shards::CellTopology &parentCell)
Computes constant normal vectors to faces of 3D reference cell.
Header file for utility class to provide multidimensional containers.
Contains definitions of custom data types in Intrepid.
static void mapToPhysicalFrame(ArrayPhysPoint &physPoints, const ArrayRefPoint &refPoints, const ArrayCell &cellWorkset, const shards::CellTopology &cellTopo, const int &whichCell=-1)
Computes F, the reference-to-physical frame map.
Header file for the Intrepid::HGRAD_TET_C2_FEM class.
static void getReferenceSubcellVertices(ArraySubcellVert &subcellVertices, const int subcellDim, const int subcellOrd, const shards::CellTopology &parentCell)
Retrieves the Cartesian coordinates of all vertices of a reference subcell.
static void getPhysicalSideNormals(ArraySideNormal &sideNormals, const ArrayJac &worksetJacobians, const int &worksetSideOrd, const shards::CellTopology &parentCell)
Computes non-normalized normal vectors to physical sides in a side workset .
static int checkPointsetInclusion(const ArrayPoint &points, const shards::CellTopology &cellTopo, const double &threshold=INTREPID_THRESHOLD)
Checks if a set of points belongs to a reference cell.
static int hasReferenceCell(const shards::CellTopology &cellTopo)
Checks if a cell topology has reference cell.
Header file for the Intrepid::HGRAD_TRI_C2_FEM class.
Intrepid utilities.
Header file for the abstract base class Intrepid::Basis.
Header file for the Intrepid::G_WEDGE_C2_FEM class.
static void mapToReferenceFrameInitGuess(ArrayRefPoint &refPoints, const ArrayInitGuess &initGuess, const ArrayPhysPoint &physPoints, const ArrayCell &cellWorkset, const shards::CellTopology &cellTopo, const int &whichCell=-1)
Computation of , the inverse of the reference-to-physical frame map using user-supplied initial guess...
static void setSubcellParametrization(FieldContainer< double > &subcellParametrization, const int subcellDim, const shards::CellTopology &parentCell)
Defines orientation-preserving parametrizations of reference edges and faces of cell topologies with ...
Header file for the Intrepid::HGRAD_HEX_I2_FEM class.
static void getReferenceSubcellNodes(ArraySubcellNode &subcellNodes, const int subcellDim, const int subcellOrd, const shards::CellTopology &parentCell)
Retrieves the Cartesian coordinates of all nodes of a reference subcell.
An abstract base class that defines interface for concrete basis implementations for Finite Element (...
static int checkPointInclusion(const Scalar *point, const int pointDim, const shards::CellTopology &cellTopo, const double &threshold=INTREPID_THRESHOLD)
Checks if a point belongs to a reference cell.
static const double * getReferenceVertex(const shards::CellTopology &cell, const int vertexOrd)
Retrieves the Cartesian coordinates of a reference cell vertex.
static void printWorksetSubcell(const ArrayCell &cellWorkset, const shards::CellTopology &parentCell, const int &pCellOrd, const int &subcellDim, const int &subcellOrd, const int &fieldWidth=3)
Prints the nodes of a subcell from a cell workset.
static void validateArguments_checkPointwiseInclusion(ArrayIncl &inCell, const ArrayPoint &physPoints, const ArrayCell &cellWorkset, const int &whichCell, const shards::CellTopology &cell)
Validates arguments to Intrepid::CellTools::checkPointwiseInclusion.
Computes the Jacobian matrix DF of the reference-to-physical frame map F.
static void validateArguments_setJacobianDetArgs(const ArrayJacDet &jacobianDet, const ArrayJac &jacobian)
Validates arguments to Intrepid::CellTools::setJacobianDet.
Header file for classes providing basic linear algebra functionality in 1D, 2D and 3D...
Header file for the Intrepid::HGRAD_HEX_C2_FEM class.
static void validateArguments_setJacobianInv(const ArrayJacInv &jacobianInv, const ArrayJac &jacobian)
Validates arguments to Intrepid::CellTools::setJacobianInv.
Computes F, the reference-to-physical frame map.
Header file for the Intrepid::HGRAD_TET_COMP12_FEM class.
static void validateArguments_mapToPhysicalFrame(const ArrayPhysPoint &physPoints, const ArrayRefPoint &refPoints, const ArrayCell &cellWorkset, const shards::CellTopology &cellTopo, const int &whichCell)
Validates arguments to Intrepid::CellTools::mapToPhysicalFrame.
static void checkPointwiseInclusion(ArrayIncl &inRefCell, const ArrayPoint &points, const shards::CellTopology &cellTopo, const double &threshold=INTREPID_THRESHOLD)
Checks every point in a set for inclusion in a reference cell.
static void validateArguments_setJacobian(const ArrayJac &jacobian, const ArrayPoint &points, const ArrayCell &cellWorkset, const int &whichCell, const shards::CellTopology &cellTopo)
Validates arguments to Intrepid::CellTools::setJacobian.
Header file for the Intrepid::HGRAD_PYR_I2_FEM class.
static const FieldContainer< double > & getSubcellParametrization(const int subcellDim, const shards::CellTopology &parentCell)
Returns array with the coefficients of the parametrization maps for the edges or faces of a reference...
static void getBarycenter(ArrayCent &barycenter, const ArrayCellCoord &cellCoords)
Compute cell barycenters.
static void getSubCVCoords(ArrayCVCoord &subCVcoords, const ArrayCellCoord &cellCoords, const shards::CellTopology &primaryCell)
Computes coordinates of sub-control volumes in each primary cell.
A stateless class for operations on cell data. Provides methods for:
static void setJacobianDet(ArrayJacDet &jacobianDet, const ArrayJac &jacobian)
Computes the determinant of the Jacobian matrix DF of the reference-to-physical frame map F...
CellTools()
Default constructor.
Header file for the Intrepid::HGRAD_LINE_C1_FEM class.
static void getPhysicalFaceTangents(ArrayFaceTangentU &faceTanU, ArrayFaceTangentV &faceTanV, const ArrayJac &worksetJacobians, const int &worksetFaceOrd, const shards::CellTopology &parentCell)
Computes non-normalized tangent vector pairs to physical faces in a face workset ; (see Subcell works...