|
Intrepid
|
00001 // @HEADER 00002 // ************************************************************************ 00003 // 00004 // Intrepid Package 00005 // Copyright (2007) Sandia Corporation 00006 // 00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive 00008 // license for use of this work by or on behalf of the U.S. Government. 00009 // 00010 // This library is free software; you can redistribute it and/or modify 00011 // it under the terms of the GNU Lesser General Public License as 00012 // published by the Free Software Foundation; either version 2.1 of the 00013 // License, or (at your option) any later version. 00014 // 00015 // This library is distributed in the hope that it will be useful, but 00016 // WITHOUT ANY WARRANTY; without even the implied warranty of 00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00018 // Lesser General Public License for more details. 00019 // 00020 // You should have received a copy of the GNU Lesser General Public 00021 // License along with this library; if not, write to the Free Software 00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 00023 // USA 00024 // Questions? Contact Pavel Bochev (pbboche@sandia.gov) 00025 // Denis Ridzal (dridzal@sandia.gov), or 00026 // Kara Peterson (kjpeter@sandia.gov) 00027 // 00028 // ************************************************************************ 00029 // @HEADER 00030 00036 #ifndef INTREPID_CELLTOOLS_HPP 00037 #define INTREPID_CELLTOOLS_HPP 00038 00039 00040 #include "Intrepid_FieldContainer.hpp" 00041 #include "Intrepid_RealSpaceTools.hpp" 00042 #include "Intrepid_ConfigDefs.hpp" 00043 #include "Intrepid_Types.hpp" 00044 #include "Intrepid_Utils.hpp" 00045 #include "Intrepid_Basis.hpp" 00046 #include "Intrepid_HGRAD_TRI_C1_FEM.hpp" 00047 #include "Intrepid_HGRAD_QUAD_C1_FEM.hpp" 00048 #include "Intrepid_HGRAD_TET_C1_FEM.hpp" 00049 #include "Intrepid_HGRAD_WEDGE_C1_FEM.hpp" 00050 #include "Intrepid_HGRAD_HEX_C1_FEM.hpp" 00051 00052 #include "Intrepid_HGRAD_LINE_C1_FEM.hpp" 00053 00054 #include "Intrepid_HGRAD_TRI_C2_FEM.hpp" 00055 #include "Intrepid_HGRAD_QUAD_C2_FEM.hpp" 00056 #include "Intrepid_HGRAD_TET_C2_FEM.hpp" 00057 #include "Intrepid_HGRAD_WEDGE_C2_FEM.hpp" 00058 #include "Intrepid_HGRAD_HEX_C2_FEM.hpp" 00059 00060 #include "Shards_CellTopology.hpp" 00061 #include "Shards_BasicTopologies.hpp" 00062 00063 #include "Teuchos_TestForException.hpp" 00064 #include "Teuchos_RCP.hpp" 00065 00066 namespace Intrepid { 00067 00068 00069 //============================================================================================// 00070 // // 00071 // CellTools // 00072 // // 00073 //============================================================================================// 00074 00083 template<class Scalar> 00084 class CellTools { 00085 private: 00086 00087 //============================================================================================// 00088 // // 00089 // Parametrization coefficients of edges and faces of reference cells // 00090 // // 00091 //============================================================================================// 00092 00093 00106 static const FieldContainer<double>& getSubcellParametrization(const int subcellDim, 00107 const shards::CellTopology& parentCell); 00108 00109 00110 00150 static void setSubcellParametrization(FieldContainer<double>& subcellParametrization, 00151 const int subcellDim, 00152 const shards::CellTopology& parentCell); 00153 00154 //============================================================================================// 00155 // // 00156 // Validation of input/output arguments for CellTools methods // 00157 // // 00158 //============================================================================================// 00159 00167 template<class ArrayJac, class ArrayPoint, class ArrayCell> 00168 static void validateArguments_setJacobian(const ArrayJac & jacobian, 00169 const ArrayPoint & points, 00170 const ArrayCell & cellWorkset, 00171 const int & whichCell, 00172 const shards::CellTopology & cellTopo); 00173 00174 00175 00180 template<class ArrayJacInv, class ArrayJac> 00181 static void validateArguments_setJacobianInv(const ArrayJacInv & jacobianInv, 00182 const ArrayJac & jacobian); 00183 00184 00185 00190 template<class ArrayJacDet, class ArrayJac> 00191 static void validateArguments_setJacobianDetArgs(const ArrayJacDet & jacobianDet, 00192 const ArrayJac & jacobian); 00193 00194 00195 00203 template<class ArrayPhysPoint, class ArrayRefPoint, class ArrayCell> 00204 static void validateArguments_mapToPhysicalFrame(const ArrayPhysPoint & physPoints, 00205 const ArrayRefPoint & refPoints, 00206 const ArrayCell & cellWorkset, 00207 const shards::CellTopology & cellTopo, 00208 const int& whichCell); 00209 00210 00211 00219 template<class ArrayRefPoint, class ArrayPhysPoint, class ArrayCell> 00220 static void validateArguments_mapToReferenceFrame(const ArrayRefPoint & refPoints, 00221 const ArrayPhysPoint & physPoints, 00222 const ArrayCell & cellWorkset, 00223 const shards::CellTopology & cellTopo, 00224 const int& whichCell); 00225 00226 00227 00236 template<class ArrayRefPoint, class ArrayInitGuess, class ArrayPhysPoint, class ArrayCell> 00237 static void validateArguments_mapToReferenceFrame(const ArrayRefPoint & refPoints, 00238 const ArrayInitGuess & initGuess, 00239 const ArrayPhysPoint & physPoints, 00240 const ArrayCell & cellWorkset, 00241 const shards::CellTopology & cellTopo, 00242 const int& whichCell); 00243 00244 00245 00253 template<class ArrayIncl, class ArrayPoint, class ArrayCell> 00254 static void validateArguments_checkPointwiseInclusion(ArrayIncl & inCell, 00255 const ArrayPoint & physPoints, 00256 const ArrayCell & cellWorkset, 00257 const int & whichCell, 00258 const shards::CellTopology & cell); 00259 public: 00260 00263 CellTools(){ }; 00264 00265 00268 ~CellTools(){ }; 00269 00270 //============================================================================================// 00271 // // 00272 // Jacobian, inverse Jacobian and Jacobian determinant // 00273 // // 00274 //============================================================================================// 00275 00323 template<class ArrayJac, class ArrayPoint, class ArrayCell> 00324 static void setJacobian(ArrayJac & jacobian, 00325 const ArrayPoint & points, 00326 const ArrayCell & cellWorkset, 00327 const shards::CellTopology & cellTopo, 00328 const int & whichCell = -1); 00329 00330 00331 00344 template<class ArrayJacInv, class ArrayJac> 00345 static void setJacobianInv(ArrayJacInv & jacobianInv, 00346 const ArrayJac & jacobian); 00347 00348 00349 00362 template<class ArrayJacDet, class ArrayJac> 00363 static void setJacobianDet(ArrayJacDet & jacobianDet, 00364 const ArrayJac & jacobian); 00365 00366 //============================================================================================// 00367 // // 00368 // Reference-to-physical frame mapping and its inverse // 00369 // // 00370 //============================================================================================// 00371 00427 template<class ArrayPhysPoint, class ArrayRefPoint, class ArrayCell> 00428 static void mapToPhysicalFrame(ArrayPhysPoint & physPoints, 00429 const ArrayRefPoint & refPoints, 00430 const ArrayCell & cellWorkset, 00431 const shards::CellTopology & cellTopo, 00432 const int & whichCell = -1); 00433 00434 00435 00494 template<class ArrayRefPoint, class ArrayPhysPoint, class ArrayCell> 00495 static void mapToReferenceFrame(ArrayRefPoint & refPoints, 00496 const ArrayPhysPoint & physPoints, 00497 const ArrayCell & cellWorkset, 00498 const shards::CellTopology & cellTopo, 00499 const int & whichCell = -1); 00500 00501 00502 00548 template<class ArrayRefPoint, class ArrayInitGuess, class ArrayPhysPoint, class ArrayCell> 00549 static void mapToReferenceFrameInitGuess(ArrayRefPoint & refPoints, 00550 const ArrayInitGuess & initGuess, 00551 const ArrayPhysPoint & physPoints, 00552 const ArrayCell & cellWorkset, 00553 const shards::CellTopology & cellTopo, 00554 const int & whichCell = -1); 00555 00556 00557 00608 template<class ArraySubcellPoint, class ArrayParamPoint> 00609 static void mapToReferenceSubcell(ArraySubcellPoint & refSubcellPoints, 00610 const ArrayParamPoint & paramPoints, 00611 const int subcellDim, 00612 const int subcellOrd, 00613 const shards::CellTopology & parentCell); 00614 00615 00616 00642 template<class ArrayEdgeTangent> 00643 static void getReferenceEdgeTangent(ArrayEdgeTangent & refEdgeTangent, 00644 const int & edgeOrd, 00645 const shards::CellTopology & parentCell); 00646 00647 00648 00685 template<class ArrayFaceTangentU, class ArrayFaceTangentV> 00686 static void getReferenceFaceTangents(ArrayFaceTangentU & refFaceTanU, 00687 ArrayFaceTangentV & refFaceTanV, 00688 const int & faceOrd, 00689 const shards::CellTopology & parentCell); 00690 00691 00692 00755 template<class ArraySideNormal> 00756 static void getReferenceSideNormal(ArraySideNormal & refSideNormal, 00757 const int & sideOrd, 00758 const shards::CellTopology & parentCell); 00759 00760 00761 00800 template<class ArrayFaceNormal> 00801 static void getReferenceFaceNormal(ArrayFaceNormal & refFaceNormal, 00802 const int & faceOrd, 00803 const shards::CellTopology & parentCell); 00804 00805 00806 00836 template<class ArrayEdgeTangent, class ArrayJac> 00837 static void getPhysicalEdgeTangents(ArrayEdgeTangent & edgeTangents, 00838 const ArrayJac & worksetJacobians, 00839 const int & worksetEdgeOrd, 00840 const shards::CellTopology & parentCell); 00841 00842 00843 00883 template<class ArrayFaceTangentU, class ArrayFaceTangentV, class ArrayJac> 00884 static void getPhysicalFaceTangents(ArrayFaceTangentU & faceTanU, 00885 ArrayFaceTangentV & faceTanV, 00886 const ArrayJac & worksetJacobians, 00887 const int & worksetFaceOrd, 00888 const shards::CellTopology & parentCell); 00889 00890 00891 00952 template<class ArraySideNormal, class ArrayJac> 00953 static void getPhysicalSideNormals(ArraySideNormal & sideNormals, 00954 const ArrayJac & worksetJacobians, 00955 const int & worksetSideOrd, 00956 const shards::CellTopology & parentCell); 00957 00958 00959 00998 template<class ArrayFaceNormal, class ArrayJac> 00999 static void getPhysicalFaceNormals(ArrayFaceNormal & faceNormals, 01000 const ArrayJac & worksetJacobians, 01001 const int & worksetFaceOrd, 01002 const shards::CellTopology & parentCell); 01003 01004 01005 //============================================================================================// 01006 // // 01007 // Inclusion tests // 01008 // // 01009 //============================================================================================// 01010 01021 static int checkPointInclusion(const Scalar* point, 01022 const int pointDim, 01023 const shards::CellTopology & cellTopo, 01024 const double & threshold = INTREPID_THRESHOLD); 01025 01026 01027 01040 template<class ArrayPoint> 01041 static int checkPointsetInclusion(const ArrayPoint & points, 01042 const shards::CellTopology & cellTopo, 01043 const double & threshold = INTREPID_THRESHOLD); 01044 01045 01046 01073 template<class ArrayIncl, class ArrayPoint> 01074 static void checkPointwiseInclusion(ArrayIncl & inRefCell, 01075 const ArrayPoint & points, 01076 const shards::CellTopology & cellTopo, 01077 const double & threshold = INTREPID_THRESHOLD); 01078 01079 01080 01116 template<class ArrayIncl, class ArrayPoint, class ArrayCell> 01117 static void checkPointwiseInclusion(ArrayIncl & inCell, 01118 const ArrayPoint & points, 01119 const ArrayCell & cellWorkset, 01120 const shards::CellTopology & cell, 01121 const int & whichCell = -1, 01122 const double & threshold = INTREPID_THRESHOLD); 01123 01124 01125 01136 static const double* getReferenceVertex(const shards::CellTopology& cell, 01137 const int vertexOrd); 01138 01139 01140 01155 template<class ArraySubcellVert> 01156 static void getReferenceSubcellVertices(ArraySubcellVert & subcellVertices, 01157 const int subcellDim, 01158 const int subcellOrd, 01159 const shards::CellTopology& parentCell); 01160 01161 01162 01178 static const double* getReferenceNode(const shards::CellTopology& cell, 01179 const int nodeOrd); 01180 01181 01182 01196 template<class ArraySubcellNode> 01197 static void getReferenceSubcellNodes(ArraySubcellNode& subcellNodes, 01198 const int subcellDim, 01199 const int subcellOrd, 01200 const shards::CellTopology& parentCell); 01201 01202 01203 01209 static int hasReferenceCell(const shards::CellTopology & cellTopo); 01210 01211 01212 01213 //============================================================================================// 01214 // // 01215 // Debug // 01216 // // 01217 //============================================================================================// 01218 01219 01225 static void printSubcellVertices(const int subcellDim, 01226 const int subcellOrd, 01227 const shards::CellTopology & parentCell); 01228 01229 01230 01234 template<class ArrayCell> 01235 static void printWorksetSubcell(const ArrayCell & cellWorkset, 01236 const shards::CellTopology & parentCell, 01237 const int& pCellOrd, 01238 const int& subcellDim, 01239 const int& subcellOrd, 01240 const int& fieldWidth = 3); 01241 01242 }; // class CellTools 01243 01244 } // namespace Intrepid 01245 01246 // include templated function definitions 01247 01248 #include "Intrepid_CellToolsDef.hpp" 01249 01250 #endif 01251 01252 /*************************************************************************************************** 01253 ** ** 01254 ** D O C U M E N T A T I O N P A G E S ** 01255 ** ** 01256 **************************************************************************************************/ 01257
1.7.4