1 /* 2 XLiFE++ is an extended library of finite elements written in C++ 3 Copyright (C) 2014 Lunéville, Eric; Kielbasiewicz, Nicolas; Lafranche, Yvon; Nguyen, Manh-Ha; Chambeyron, Colin 4 5 This program is free software: you can redistribute it and/or modify 6 it under the terms of the GNU General Public License as published by 7 the Free Software Foundation, either version 3 of the License, or 8 (at your option) any later version. 9 This program is distributed in the hope that it will be useful, 10 but WITHOUT ANY WARRANTY; without even the implied warranty of 11 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 12 GNU General Public License for more details. 13 You should have received a copy of the GNU General Public License 14 along with this program. If not, see <http://www.gnu.org/licenses/>. 15 */ 16 17 /*! 18 \file GeomMapData.hpp 19 \author D. Martin, E. Lunéville 20 \since 13 apr 2012 21 \date 23 may 2012 22 23 \brief Definition of the xlifepp::GeomMapData class 24 25 A xlifepp::GeomMapData is a data structure to store data related to the map from reference element 26 to geometric element (jacobian matrix, its determinant and the differential element,...) 27 */ 28 29 30 #ifndef GEOM_MAP_DATA_HPP 31 #define GEOM_MAP_DATA_HPP 32 33 #include "utils.h" 34 35 namespace xlifepp 36 { 37 //forward declaration 38 class MeshElement; 39 class Quadrature; 40 class ShapeValues; 41 class ExtensionData; 42 43 //------------------------------------------------------------------------------------ 44 /*! 45 \class GeomMapData 46 object handles useful computational data of a geometric element 47 It stores jacobian, its determinant, ... depending of a point 48 */ 49 //------------------------------------------------------------------------------------ 50 class GeomMapData 51 { 52 private: 53 const MeshElement* geomElement_p; //!< current geometric element 54 Point currentPoint; //!< point in reference element used for last computation 55 56 public: 57 Matrix<real_t> jacobianMatrix; //!< jacobian matrix of mapping from reference element to current one 58 Matrix<real_t> inverseJacobianMatrix; //!< inverse jacobian matrix (when needed) 59 real_t jacobianDeterminant; //!< jacobian determinant 60 real_t differentialElement; //!< differential element for elementary integrals (abs(jac. det.)) 61 Vector<real_t> normalVector; //!< unit outward normal vector at a given boundary point 62 Matrix<real_t> metricTensor; //!< symetric metric tensor (t_i|t_j) 63 real_t metricTensorDeterminant; //!< metric_tensor determinant for differential geometry 64 65 dimen_t elementDim; //!< dim of reference space 66 dimen_t spaceDim; //!< dim of physical space 67 68 std::map<Quadrature*,std::vector<Point> > phyPoints; //!< to store image of quadrature points in physical space 69 ExtensionData* extdata; //!< pointer to additionnal data used by extension 70 mutable std::vector< Vector<real_t> > sideNormalVectors; //!< normal to the sides of the element (build if required) 71 std::vector< Vector<real_t> >& sideNV() const; //!< get normals to the sides of the element 72 73 //constructor 74 GeomMapData(const MeshElement*, const Point&); //!< basic constructor 75 GeomMapData(const MeshElement*, std::vector<real_t>::const_iterator); //!< basic constructor 76 GeomMapData(const MeshElement*); //!< basic constructor (element centroid as point) 77 78 //point mapping 79 Point geomMap(number_t side = 0); //!< mapping from ref. elt onto geom. elt. Element (shape functions known) 80 Point geomMap(const ShapeValues& shv, number_t side); //!< mapping from ref. elt onto geom. elt. Element (shape functions given, thread safe) 81 Point geomMap(const std::vector<real_t>& p, number_t side = 0);//!< mapping from ref. elt onto geom. elt. Element of point p 82 Point geomMap(std::vector<real_t>::const_iterator p, number_t side = 0); //!< mapping from ref. elt onto geom. elt. Element of point p 83 Point geomMapInverse(const std::vector<real_t>&, real_t eps=theTolerance, 84 number_t maxIter = 50); //!< mapping from elt onto ref. elt. of point p 85 Point piolaMap(number_t side = 0); //!< Piola mapping from ref. elt onto geom. elt. Element (shape functions known) 86 Matrix<real_t> covariantPiolaMap(const Point& =Point()); //!< return the covariant Piola map matrix at current point 87 Matrix<real_t> contravariantPiolaMap(const Point& =Point()); //!< return the contravariant Piola map matrix at current point 88 89 //geometric computation functions 90 void computeJacobianMatrix(number_t side = 0); //!< compute jacobian matrix (shape functions known) 91 void computeJacobianMatrix(const ShapeValues& shv, number_t side =0); //!< compute jacobian matrix (shape functions given) 92 void computeJacobianMatrix(const std::vector<real_t>&, 93 number_t side = 0); //!< compute jacobian matrix at a point given as vector or Point 94 void computeJacobianMatrix(std::vector<real_t>::const_iterator, 95 number_t side=0); //!< compute jacobian matrix at a point given by iterator 96 97 real_t computeJacobianDeterminant(); //!< compute jacobian determinant (jacobian already computed) 98 void invertJacobianMatrix(); //!< compute inverse of jacobian matrix (jacobian already computed) 99 void computeDifferentialElement(); //!< compute differential element assuming jacobian matrix is update 100 void computeNormalVector(); //!< computes normal vector at a (boundary) point of element 101 real_t diffElement(); //!< update and return differential element 102 real_t diffElement(number_t); //!< update and return differential element on side 103 void normalize(); //!< normalize normal vector to element 104 void computeOutwardNormal(); //!< compute outward unit normal vector, assuming jacobian matrix is update 105 void computeMetricTensor(); //!< compute metric tensor for surface differential geometry 106 void computeSurfaceGradient(real_t, real_t, std::vector<real_t>&); //!< compute surface gradient gradS from 2D reference gradient (grad0, grad1) 107 void print(std::ostream&) const; //!< print GeoMapData print(PrintStream & os) const108 void print(PrintStream& os) const {print(os.currentStream());} 109 110 friend std::ostream& operator<<(std::ostream&, const GeomMapData&); 111 }; 112 113 std::ostream& operator<<(std::ostream&, const GeomMapData&); //!< print operator 114 115 } // end of namespace xlifepp 116 117 #endif // GEOM_MAP_DATA_HPP 118