BELFEM 0.9.0
Berkeley Lab Finite Element Framework
Loading...
Searching...
No Matches
cl_IFG_TET10C.hpp
Go to the documentation of this file.
1/*
2 * BELFEM -- The Berkeley Lab Finite Element Framework
3 * Copyright (c) 2026, The Regents of the University of California,
4 * through Lawrence Berkeley National Laboratory (subject to receipt of any required
5 * approvals from the U.S. Dept. of Energy). All rights reserved.
6 *
7 * Developers: Christian Messe, Gregory Giard
8 *
9 * See the top-level LICENSE file for the complete license and disclaimer.
10 */
11
12#ifndef BELFEM_CL_IFG_TET10C_HPP
13#define BELFEM_CL_IFG_TET10C_HPP
14
15#include "cl_IF_InterpolationFunctionTemplate.hpp"
16
17namespace belfem
18{
19 namespace fem
20 {
21//------------------------------------------------------------------------------
22
30
31//------------------------------------------------------------------------------
32
33 template<>
34 InterpolationOrder
35 InterpolationFunctionTemplate<
36 GeometryType::TET, InterpolationType::BubbleFace2, 3, 3 >
41
42//------------------------------------------------------------------------------
43
44 template<>
45 ElementType
46 InterpolationFunctionTemplate<
47 GeometryType::TET, InterpolationType::BubbleFace2, 3, 3 >
48 ::element_type() const
49 {
50 return ElementType::TET10;
51 }
52
53//------------------------------------------------------------------------------
54
55
56 template<>
57 void
58 InterpolationFunctionTemplate<
59 GeometryType::TET, InterpolationType::BubbleFace2, 3, 3 >
60 ::param_coords( Matrix <real> & aXiHat ) const
61 {
62 aXiHat.set_size( 3, 3 );
63
64 // point 0
65 aXiHat( 0, 0 ) = 2./3. ;
66 aXiHat( 0, 1 ) = 0.0 ;
67 aXiHat( 0, 2 ) = 1./6. ;
68
69 // point 1
70 aXiHat( 1, 0 ) = 1./6. ;
71 aXiHat( 1, 1 ) = 0.0 ;
72 aXiHat( 1, 2 ) = 1./6. ;
73
74 // point 2
75 aXiHat( 2, 0 ) = 1./6. ;
76 aXiHat( 2, 1 ) = 0.0 ;
77 aXiHat( 2, 2 ) = 2./3. ;
78 }
79
80//------------------------------------------------------------------------------
81
82 template<>
83 void
84 InterpolationFunctionTemplate<
85 GeometryType::TET, InterpolationType::BubbleFace2, 3, 3 >::N(
86 const Vector <real> & aXi,
87 Matrix <real> & aN ) const
88 {
89 const real xi = aXi( 0 );
90 const real eta = aXi( 1 );
91 const real zeta = aXi( 2 );
92 const real tau = 1. - xi - eta - zeta ;
93
94 real K = 64.8 ;
95
96 real xi2 = xi*xi ;
97 real zeta2 = zeta*zeta ;
98 real tau2 = tau*tau ;
99
100 real a = 2.*tau2 - xi2 - zeta2 ;
101 real b = 2.*xi2 - zeta2 - tau2 ;
102 real c = 2.*zeta2 - xi2 - tau2 ;
103 real d = xi*zeta*tau ;
104
105 aN.set_size( 1, 3 );
106 aN( 0, 0 ) = K * d * a ;
107 aN( 0, 1 ) = K * d * b ;
108 aN( 0, 2 ) = K * d * c ;
109 }
110
111
112//------------------------------------------------------------------------------
113
114 template<>
115 void
116 InterpolationFunctionTemplate<
117 GeometryType::TET, InterpolationType::BubbleFace2, 3, 3 >::dNdXi(
118 const Vector <real> & aXi,
119 Matrix <real> & adNdXi ) const
120 {
121 const real xi = aXi( 0 );
122 const real eta = aXi( 1 );
123 const real zeta = aXi( 2 );
124 const real tau = 1. - xi - eta - zeta ;
125
126 real K = 64.8 ;
127
128 real xi2 = xi*xi ;
129 real zeta2 = zeta*zeta ;
130 real tau2 = tau*tau ;
131
132 real a = 2.*tau2 - xi2 - zeta2 ;
133 real b = 2.*xi2 - zeta2 - tau2 ;
134 real c = 2.*zeta2 - xi2 - tau2 ;
135 real d = xi*zeta*tau ;
136
137 real a_xi = 2.* ( -xi - 2.* tau );
138 real a_eta = -4. * tau ;
139 real a_zeta = 2.* ( -zeta -2*tau );
140
141 real b_xi = 2. * ( tau + 2. * xi );
142 real b_eta = 2. * tau ;
143 real b_zeta = 2.* ( tau - zeta ) ;
144
145 real c_xi = 2. * ( tau - xi ) ;
146 real c_eta = 2. * tau ;
147 real c_zeta = 2. * ( tau + 2.*zeta );
148
149 real d_xi = zeta * ( tau - xi );
150 real d_eta = -xi*zeta ;
151 real d_zeta = xi*( tau - zeta );
152
153 adNdXi.set_size( 3, 3 );
154
155 adNdXi( 0, 0 ) = K * ( d_xi * a + d * a_xi ); // dN1/dxi
156 adNdXi( 0, 1 ) = K * ( d_eta * a + d * a_eta ); // dN1/deta
157 adNdXi( 0, 2 ) = K * ( d_zeta * a + d * a_zeta ); // dN1/dzeta
158
159 adNdXi( 1, 0 ) = K * ( d_xi * b + d * b_xi ); // dN2/dxi
160 adNdXi( 1, 1 ) = K * ( d_eta * b + d * b_eta ); // dN2/deta
161 adNdXi( 1, 2 ) = K * ( d_zeta * b + d * b_zeta ); // dN2/dzeta
162
163 adNdXi( 2, 0 ) = K * ( d_xi * c + d * c_xi ); // dN3/dxi
164 adNdXi( 2, 1 ) = K * ( d_eta * c + d * c_eta ); // dN3/deta
165 adNdXi( 2, 2 ) = K * ( d_zeta * c + d * c_zeta ); // dN3/dzeta
166
167 }
168
169//------------------------------------------------------------------------------
170
171
172 template<>
173 void
174 InterpolationFunctionTemplate<
175 GeometryType::TET, InterpolationType::BubbleFace2, 3, 3 >::d2NdXi2(
176 const Vector <real> & aXi,
177 Matrix <real> & ad2NdXi2 ) const
178 {
179 const real xi = aXi( 0 );
180 const real eta = aXi( 1 );
181 const real zeta = aXi( 2 );
182 const real tau = 1. - xi - eta - zeta ;
183
184 real K = 64.8 ;
185
186 real xi2 = xi*xi ;
187 real zeta2 = zeta*zeta ;
188 real tau2 = tau*tau ;
189
190 real a = 2.*tau2 - xi2 - zeta2 ;
191 real b = 2.*xi2 - zeta2 - tau2 ;
192 real c = 2.*zeta2 - xi2 - tau2 ;
193 real d = xi*zeta*tau ;
194
195 real a_xi = 2.* ( -xi - 2.* tau );
196 real a_eta = -4. * tau ;
197 real a_zeta = 2.* ( -zeta -2*tau );
198
199 real b_xi = 2. * ( tau + 2. * xi );
200 real b_eta = 2. * tau ;
201 real b_zeta = 2.* ( tau - zeta ) ;
202
203 real c_xi = 2. * ( tau - xi ) ;
204 real c_eta = 2. * tau ;
205 real c_zeta = 2. * ( tau + 2.*zeta );
206
207 real d_xi = zeta * ( tau - xi );
208 real d_eta = -xi*zeta ;
209 real d_zeta = xi*( tau - zeta );
210
211 ad2NdXi2.set_size( 6, 3 );
212
213 // d²N1/dxi²
214 ad2NdXi2( 0, 0 ) = 2.*K*(a_xi*d_xi - a*zeta+d);
215
216 // d²N1/deta²
217 ad2NdXi2( 1, 0 ) = 2.*K*(a_eta*d_eta+2.*d);
218
219 // d²N1/dzeta²
220 ad2NdXi2( 2, 0 ) = 2.*K*(a_zeta*d_zeta - a*xi+d);
221
222
223 // d²N1/(deta*dzeta)
224 ad2NdXi2( 3, 0 ) = K*(a_eta*d_zeta + a_zeta*d_eta - a*xi + 4.*d);
225
226 // d²N1/(dxi*dzeta)
227 ad2NdXi2( 4, 0 ) = K*(a_xi*d_zeta + a_zeta*d_xi - a*(xi+zeta-tau) + 4.*d);
228
229 // d²N1/(dxi*deta)
230 ad2NdXi2( 5, 0 ) = K*(a_eta*d_xi + a_xi*d_eta - a*zeta + 4.*d);
231
232 // d²N2/dxi²
233 ad2NdXi2( 0, 1 ) = 2.*K*(b_xi*d_xi - b*zeta+d) ;
234
235 // d²N2/deta²
236 ad2NdXi2( 1, 1 ) = 2.*K*(b_eta*d_eta-d);
237
238 // d²N2/dzeta²
239 ad2NdXi2( 2, 1 ) = 2.*K*(b_zeta*d_zeta - b*xi-2.*d);
240
241 // d²N2/(deta*dzeta)
242 ad2NdXi2( 3, 1 ) = K*( b_eta*d_zeta + b_zeta*d_eta - b*xi - 2.*d);
243
244 // d²N2/(dxi*dzeta)
245 ad2NdXi2( 4, 1 ) = K*(b_xi*d_zeta + b_zeta*d_xi - b*(xi+zeta-tau) - 2.*d);
246
247 // d²N2/(dxi*deta)
248 ad2NdXi2( 5, 1 ) = K*(b_eta*d_xi + b_xi*d_eta - b*zeta - 2.*d);
249
250 // d²N3/dxi²
251 ad2NdXi2( 0, 2 ) = 2.*K*( c_xi*d_xi - c*zeta-2.*d) ;
252
253 // d²N3/deta²
254 ad2NdXi2( 1, 2 ) = 2.*K*( c_eta*d_eta-d);
255
256 // d²N3/dzeta²
257 ad2NdXi2( 2, 2 ) = 2.*K*( c_zeta*d_zeta - c*xi+d);
258
259 // d²N3/(deta*dzeta)
260 ad2NdXi2( 3, 2 ) = K*(c_eta*d_zeta + c_zeta*d_eta - c*xi-2.*d);
261
262 // d²N3/(dxi*dzeta)
263 ad2NdXi2( 4, 2 ) = K*(c_xi*d_zeta + c_zeta*d_xi - c*(xi+zeta-tau)-2.*d);
264
265 // d²N3/(dxi*deta)
266 ad2NdXi2( 5, 2 ) = K*(c_eta*d_xi + c_xi*d_eta - c*zeta-2.*d);
267
268 }
269
270//------------------------------------------------------------------------------
271 }
272}
273#endif //BELFEM_CL_IFG_TET10C_HPP
cl_IF_InterpolationFunctionTemplate.hpp
belfem::Matrix::set_size
void set_size(const size_t aNumRows, const size_t aNumCols)
Definition cl_AR_Matrix.hpp:186
belfem::fem::InterpolationFunctionTemplate
shape function templated class G : Geometry T : Type D : Dimension B : Number of Basis
Definition cl_IF_InterpolationFunctionTemplate.hpp:25
belfem::fem::InterpolationFunctionTemplate::param_coords
void param_coords(Matrix< real > &aXiHat) const override
returns a matrix containing the parameter coordinates of the nodes < number of dimensions x number of...
Definition cl_IF_InterpolationFunctionTemplate.hpp:49
belfem::fem::InterpolationFunctionTemplate::interpolation_order
InterpolationOrder interpolation_order() const override
returns the interpolation order
Definition cl_IF_InterpolationFunctionTemplate.hpp:145
belfem::fem::InterpolationFunctionTemplate::d2NdXi2
void d2NdXi2(const Vector< real > &aXi, Matrix< real > &ad2NdXi2) const override
calculates the second derivative of the shape function in parameter space
Definition cl_IF_InterpolationFunctionTemplate.hpp:110
belfem::fem::InterpolationFunctionTemplate::dNdXi
void dNdXi(const Vector< real > &aXi, Matrix< real > &adNdXi) const override
calculates the first derivative of the shape function in parameter space
Definition cl_IF_InterpolationFunctionTemplate.hpp:89
belfem::fem::InterpolationFunctionTemplate::N
void N(const Vector< real > &aXi, Matrix< real > &aN) const override
evaluates the shape function at a given point
Definition cl_IF_InterpolationFunctionTemplate.hpp:68
belfem::fem
Definition cl_IFB_LINE3.hpp:21
belfem::fem::MatrixFlag::K
@ K
Stiffness matrix.
Definition cl_TimestepMatrices.hpp:35
belfem
USER GUIDES:
Definition cl_Capacitor.cpp:16
belfem::InterpolationType::BubbleFace2
@ BubbleFace2
Definition Mesh_Enums.hpp:108
belfem::element_type
ElementType element_type(const std::string &aStr)
Definition Mesh_Enums.hpp:370
belfem::ElementType
ElementType
Element types.
Definition Mesh_Enums.hpp:27
belfem::ElementType::TET10
@ TET10
Definition Mesh_Enums.hpp:39
belfem::InterpolationOrder
InterpolationOrder
Definition Mesh_Enums.hpp:85
belfem::InterpolationOrder::CUBIC
@ CUBIC
Definition Mesh_Enums.hpp:90
belfem::real
double real
Definition typedefs.hpp:36
belfem::GeometryType::TET
@ TET
Definition Mesh_Enums.hpp:75