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_TETRA_4_HPP_ 7 #define MINT_TETRA_4_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 Tetrahedron 23 * 24 * \verbatim 25 * 26 * tetra_4: 27 * 28 * 3 29 * x 30 * /|\ 31 * / | \ 32 * / | \ 33 * / | \ 34 * 0 x----|----x 2 35 * \ | / 36 * \ | / 37 * \ | / 38 * x 39 * 1 40 * 41 * \endverbatim 42 * 43 * \see ShapeFunction 44 */ 45 template <> 46 class Lagrange<mint::TET> : public ShapeFunction<Lagrange<mint::TET>> 47 { 48 public: getCellType()49 static CellType getCellType() { return mint::TET; } 50 getType()51 static int getType() { return MINT_LAGRANGE_BASIS; } 52 getNumDofs()53 static int getNumDofs() { return 4; } 54 getMaxNewtonIters()55 static int getMaxNewtonIters() { return 16; } 56 getDimension()57 static int getDimension() { return 3; } 58 getMin()59 static double getMin() { return 0; } 60 getMax()61 static double getMax() { return 1; } 62 getCenter(double * center)63 static void getCenter(double* center) 64 { 65 SLIC_ASSERT(center != nullptr); 66 center[0] = center[1] = center[2] = 0.25; 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] = 0.0; 85 coords[7] = 1.0; 86 coords[8] = 0.0; 87 88 // node 3 89 coords[9] = 0.0; 90 coords[10] = 0.0; 91 coords[11] = 1.0; 92 } 93 computeShape(const double * xr,double * phi)94 static void computeShape(const double* xr, double* phi) 95 { 96 SLIC_ASSERT(xr != nullptr); 97 SLIC_ASSERT(phi != nullptr); 98 99 const double r = xr[0]; 100 const double s = xr[1]; 101 const double t = xr[2]; 102 103 phi[0] = 1 - r - s - t; 104 phi[1] = r; 105 phi[2] = s; 106 phi[3] = t; 107 } 108 computeDerivatives(const double * AXOM_UNUSED_PARAM (xr),double * phidot)109 static void computeDerivatives(const double* AXOM_UNUSED_PARAM(xr), 110 double* phidot) 111 { 112 SLIC_ASSERT(phidot != nullptr); 113 114 // r derivatives 115 phidot[0] = -1; 116 phidot[1] = 1; 117 phidot[2] = 0; 118 phidot[3] = 0; 119 120 // s derivatives 121 phidot[4] = -1; 122 phidot[5] = 0; 123 phidot[6] = 1; 124 phidot[7] = 0; 125 126 // t derivatives 127 phidot[8] = -1; 128 phidot[9] = 0; 129 phidot[10] = 0; 130 phidot[11] = 1; 131 } 132 }; 133 134 } /* namespace mint */ 135 } /* namespace axom */ 136 #endif /* MINT_TETRA_4_HPP_ */ 137