1 // Copyright (c) 2017-2021, Lawrence Livermore National Security, LLC and 2 // other Axom Project Developers. See the top-level LICENSE file for details. 3 // 4 // SPDX-License-Identifier: (BSD-3-Clause) 5 6 #ifndef MINT_LAGRANGE_PYRA_5_HPP_ 7 #define MINT_LAGRANGE_PYRA_5_HPP_ 8 9 // Mint includes 10 #include "axom/mint/mesh/CellTypes.hpp" 11 #include "axom/mint/fem/FEBasisTypes.hpp" 12 #include "axom/mint/fem/shape_functions/ShapeFunction.hpp" 13 14 // Slic includes 15 #include "axom/slic/interface/slic.hpp" 16 17 namespace axom 18 { 19 namespace mint 20 { 21 /*! 22 * \brief Lagrange Finite Element definition for the Linear Pyramnid 23 * 24 * \verbatim 25 * 26 * pyra_5: 27 * 28 * 4 29 * / |\ 30 * / | \ 31 * / | \ 32 * 0/_ __|_ 3\ 33 * \ | \ 34 * \ _ |_ _ \ 35 * 1 2 36 * 37 * \endverbatim 38 * 39 * \warning The jacobian for the pyramid elements may become singular near the 40 * apex in which case the isoparametric mapping will fail. 41 * 42 * \see ShapeFunction 43 */ 44 template <> 45 class Lagrange<mint::PYRAMID> : public ShapeFunction<Lagrange<mint::PYRAMID>> 46 { 47 public: getCellType()48 static CellType getCellType() { return mint::PYRAMID; } 49 getType()50 static int getType() { return MINT_LAGRANGE_BASIS; } 51 getNumDofs()52 static int getNumDofs() { return 5; } 53 getMaxNewtonIters()54 static int getMaxNewtonIters() { return 16; } 55 getDimension()56 static int getDimension() { return 3; } 57 getMin()58 static double getMin() { return 0; } 59 getMax()60 static double getMax() { return 1; } 61 getCenter(double * center)62 static void getCenter(double* center) 63 { 64 SLIC_ASSERT(center != nullptr); 65 center[0] = center[1] = 0.4; 66 center[2] = 0.2; 67 } 68 getCoords(double * coords)69 static void getCoords(double* coords) 70 { 71 SLIC_ASSERT(coords != nullptr); 72 73 // node 0 74 coords[0] = 0.0; 75 coords[1] = 0.0; 76 coords[2] = 0.0; 77 78 // node 1 79 coords[3] = 1.0; 80 coords[4] = 0.0; 81 coords[5] = 0.0; 82 83 // node 2 84 coords[6] = 1.0; 85 coords[7] = 1.0; 86 coords[8] = 0.0; 87 88 // node 3 89 coords[9] = 0.0; 90 coords[10] = 1.0; 91 coords[11] = 0.0; 92 93 // node 4 94 coords[12] = 0.0; 95 coords[13] = 0.0; 96 coords[14] = 1.0; 97 } 98 computeShape(const double * xr,double * phi)99 static void computeShape(const double* xr, double* phi) 100 { 101 SLIC_ASSERT(xr != nullptr); 102 SLIC_ASSERT(phi != nullptr); 103 104 const double r = xr[0]; 105 const double s = xr[1]; 106 const double t = xr[2]; 107 const double rm = 1. - r; 108 const double sm = 1. - s; 109 const double tm = 1. - t; 110 111 phi[0] = rm * sm * tm; 112 phi[1] = r * sm * tm; 113 phi[2] = r * s * tm; 114 phi[3] = rm * s * tm; 115 phi[4] = t; 116 } 117 computeDerivatives(const double * xr,double * phidot)118 static void computeDerivatives(const double* xr, double* phidot) 119 { 120 SLIC_ASSERT(xr != nullptr); 121 SLIC_ASSERT(phidot != nullptr); 122 123 const double r = xr[0]; 124 const double s = xr[1]; 125 const double t = xr[2]; 126 const double rm = 1. - r; 127 const double sm = 1. - s; 128 const double tm = 1. - t; 129 130 // r-derivatives 131 phidot[0] = -sm * tm; 132 phidot[1] = sm * tm; 133 phidot[2] = s * tm; 134 phidot[3] = -s * tm; 135 phidot[4] = 0.0; 136 137 // s-derivatives 138 phidot[5] = -rm * tm; 139 phidot[6] = -r * tm; 140 phidot[7] = r * tm; 141 phidot[8] = rm * tm; 142 phidot[9] = 0.0; 143 144 // t-derivatives 145 phidot[10] = -rm * sm; 146 phidot[11] = -r * sm; 147 phidot[12] = -r * s; 148 phidot[13] = -rm * s; 149 phidot[14] = 1.0; 150 } 151 }; 152 153 } /* namespace mint */ 154 } /* namespace axom */ 155 156 #endif /* MINT_LAGRANGE_PYRA_5_HPP_ */ 157