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 LagrangeSegment.hpp 19 \authors D. Martin, N. Kielbasiewicz 20 \since 11 jan 2003 21 \date 6 aug 2012 22 23 \brief Definition of the xlifepp::LagrangeSegment class 24 25 xlifepp::LagrangeSegment defines Lagrange Reference Element interpolation data on a segment 26 27 child classes 28 - xlifepp::LagrangeStdSegment 29 Lagrange standard Reference Element (interpolation at equidistant nodes) 30 - xlifepp::LagrangeGLSegment 31 Lagrange Gauss-Lobatto Reference Element (interpolation at Gauss-Lobatto points) 32 33 member functions 34 - interpolationData interpolation data 35 - outputAsP1 36 - pointCoordinates point coordinates 37 38 Local node numbering of Lagrange Pk 1D elements, k=0,1, 2,.., 5 39 40 \verbatim 41 -o- o---o o---o---o o---o---o---o o---o---o---o---o o---o---o---o---o---o 42 1 2 1 2 3 1 2 4 3 1 2 5 4 3 1 2 6 5 4 3 1 43 k=0 k=1 k=2 k=3 k=4 k=5 44 \endverbatim 45 */ 46 47 #ifndef LAGRANGE_SEGMENT_HPP 48 #define LAGRANGE_SEGMENT_HPP 49 50 #include "config.h" 51 #include "RefSegment.hpp" 52 53 namespace xlifepp 54 { 55 56 class Interpolation; 57 58 /*! 59 \class LagrangeSegment 60 Grand-Parent class: ReferenceElement 61 Parent class: ReferenceSegment 62 Child classes: LagrangeStdSegment LagrangeGLSegment 63 */ 64 class LagrangeSegment : public RefSegment 65 { 66 public: 67 LagrangeSegment(const Interpolation* interp_p); //!< constructor by interpolation 68 virtual ~LagrangeSegment(); //!< destructor 69 virtual void pointCoordinates() = 0; //!< compute point coordinates, must exist il all subclasses 70 protected: 71 void interpolationData(); //!< build necessary interpolation data 72 void sideNumbering(); //!< numbering of side D.O.F's in standard Lagrange elements 73 std::vector<std::vector<number_t> > splitP1() const; //!< return nodes numbers of first order segment elements when splitting current element 74 splitvec_t splitO1() const; //!< return nodes numbers of first order segment elements when splitting current element 75 void outputAsP1(std::ofstream& os, const int, const std::vector<number_t>& ranks) const; //!< output of Lagrange D.o.f's numbers on underlyin P1 D.o.f's 76 77 }; // end of class LagrangeSegment 78 79 /*! 80 \class LagrangeStdSegment 81 child class of LagrangeHexahedron with standard interpolation 82 83 */ 84 class LagrangeStdSegment : public LagrangeSegment 85 { 86 public: 87 LagrangeStdSegment(const Interpolation* interp_p); //!< constructor by interpolation ~LagrangeStdSegment()88 ~LagrangeStdSegment() {} //!< destructor 89 90 void pointCoordinates(); //!< point coordinates (!) (equidistant on [0,1]) 91 void computeShapeFunctions(); //!< compute polynomials representation 92 using RefElement::computeShapeValues; //to unmask computeShapeValues(vector<real_t>::const_iterator, const bool) 93 void computeShapeValues(std::vector<real_t>::const_iterator it_pt, ShapeValues& shv, const bool withDeriv = true) const; //!< compute shape functions at given point getO1splitting() const94 const splitvec_t& getO1splitting() const // override 95 { return splitO1Scheme; } 96 private: 97 // storage of order 1 splitting scheme (intended to graphical softwares usage) 98 splitvec_t splitO1Scheme; 99 }; // end of class LagrangeStdSegment 100 101 /*! 102 \class LagrangeGLSegment 103 child class of LagrangeHexahedron with Gauss-Lobatto interpolation 104 */ 105 class LagrangeGLSegment : public LagrangeSegment 106 { 107 public: 108 LagrangeGLSegment(const Interpolation* interp_p); //!< constructor by interpolation ~LagrangeGLSegment()109 ~LagrangeGLSegment() {} //!< destructor 110 111 void pointCoordinates(); //!< compute point coordinates (!) (from Gauss-Lobatto points on [-1,+1]) 112 using RefElement::computeShapeValues; //to unmask computeShapeValues(vector<real_t>::const_iterator, const bool) 113 void computeShapeValues(std::vector<real_t>::const_iterator it_pt, ShapeValues& shv, const bool withDeriv = true) const; //!< compute shape functions at given point getO1splitting() const114 const splitvec_t& getO1splitting() const // override 115 { return splitO1Scheme; } 116 private: 117 // storage of order 1 splitting scheme (intended to graphical softwares usage) 118 splitvec_t splitO1Scheme; 119 }; // end of class LagrangeGLSegment 120 121 } // end of namespace xlifepp 122 123 #endif /* LAGRANGE_SEGMENT_HPP */ 124