Line data Source code
1 : // The libMesh Finite Element Library.
2 : // Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 :
4 : // This library is free software; you can redistribute it and/or
5 : // modify it under the terms of the GNU Lesser General Public
6 : // License as published by the Free Software Foundation; either
7 : // version 2.1 of the License, or (at your option) any later version.
8 :
9 : // This library 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 GNU
12 : // Lesser General Public License for more details.
13 :
14 : // You should have received a copy of the GNU Lesser General Public
15 : // License along with this library; if not, write to the Free Software
16 : // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 :
18 :
19 :
20 : // Local includes
21 : #include "libmesh/elem.h"
22 : #include "libmesh/enum_to_string.h"
23 : #include "libmesh/fe.h"
24 : #include "libmesh/fe_interface.h"
25 : #include "libmesh/fe_macro.h"
26 :
27 : namespace libMesh
28 : {
29 :
30 : // ------------------------------------------------------------
31 : // Hermite-specific implementations
32 :
33 : // Anonymous namespace for local helper functions
34 : namespace {
35 :
36 4653679 : void hermite_nodal_soln(const Elem * elem,
37 : const Order order,
38 : const std::vector<Number> & elem_soln,
39 : std::vector<Number> & nodal_soln,
40 : const bool add_p_level)
41 : {
42 4653679 : const unsigned int n_nodes = elem->n_nodes();
43 :
44 4653679 : const ElemType elem_type = elem->type();
45 :
46 4653679 : nodal_soln.resize(n_nodes);
47 :
48 : // FEType object to be passed to various FEInterface functions below.
49 422959 : FEType fe_type(order, HERMITE);
50 :
51 : const unsigned int n_sf =
52 4653679 : FEInterface::n_shape_functions(fe_type, elem, add_p_level);
53 :
54 845918 : std::vector<Point> refspace_nodes;
55 4653679 : FEBase::get_refspace_nodes(elem_type,refspace_nodes);
56 422959 : libmesh_assert_equal_to (refspace_nodes.size(), n_nodes);
57 422959 : libmesh_assert_equal_to (elem_soln.size(), n_sf);
58 :
59 : // Zero before summation
60 422959 : std::fill(nodal_soln.begin(), nodal_soln.end(), 0);
61 :
62 14011064 : for (unsigned int n=0; n<n_nodes; n++)
63 : // u_i = Sum (alpha_i phi_i)
64 47459570 : for (unsigned int i=0; i<n_sf; i++)
65 41557728 : nodal_soln[n] += elem_soln[i] *
66 41557737 : FEInterface::shape(fe_type, elem, i, refspace_nodes[n], add_p_level);
67 :
68 4653679 : } // hermite_nodal_soln()
69 :
70 :
71 :
72 12670752 : unsigned int HERMITE_n_dofs(const ElemType t, const Order o)
73 : {
74 : #ifdef DEBUG
75 5787202 : libmesh_error_msg_if(o < 3, "Error: Hermite elements require order>=3, but you asked for order=" << o);
76 : #endif
77 :
78 : // Piecewise (bi/tri)cubic C1 Hermite splines
79 12670752 : switch (t)
80 : {
81 0 : case NODEELEM:
82 0 : return 1;
83 5718494 : case EDGE2:
84 5718494 : libmesh_assert_less (o, 4);
85 : libmesh_fallthrough();
86 : case EDGE3:
87 12439275 : return (o+1);
88 :
89 0 : case QUAD4:
90 : case QUADSHELL4:
91 : case QUAD8:
92 : case QUADSHELL8:
93 0 : libmesh_assert_less (o, 4);
94 : libmesh_fallthrough();
95 : case QUAD9:
96 : case QUADSHELL9:
97 171231 : return ((o+1)*(o+1));
98 :
99 0 : case HEX8:
100 : case HEX20:
101 0 : libmesh_assert_less (o, 4);
102 : libmesh_fallthrough();
103 : case HEX27:
104 60246 : return ((o+1)*(o+1)*(o+1));
105 :
106 0 : case INVALID_ELEM:
107 0 : return 0;
108 :
109 0 : default:
110 0 : libmesh_error_msg("ERROR: Bad ElemType = " << Utility::enum_to_string(t) << " for " << Utility::enum_to_string(o) << " order approximation!");
111 : }
112 : } // HERMITE_n_dofs()
113 :
114 :
115 :
116 12008165 : unsigned int HERMITE_n_dofs(const Elem * e, const Order o)
117 : {
118 5787202 : libmesh_assert(e);
119 12670752 : return HERMITE_n_dofs(e->type(), o);
120 : }
121 :
122 :
123 :
124 154016342 : unsigned int HERMITE_n_dofs_at_node(const ElemType t,
125 : const Order o,
126 : const unsigned int n)
127 : {
128 10583568 : libmesh_assert_greater (o, 2);
129 : // Piecewise (bi/tri)cubic C1 Hermite splines
130 154016342 : switch (t)
131 : {
132 0 : case NODEELEM:
133 0 : return 1;
134 150640478 : case EDGE2:
135 : case EDGE3:
136 : {
137 150640478 : switch (n)
138 : {
139 10267916 : case 0:
140 : case 1:
141 10267916 : return 2;
142 61159 : case 2:
143 : // Interior DoFs are carried on Elems
144 : // return (o-3);
145 61159 : return 0;
146 :
147 0 : default:
148 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for EDGE2/3!");
149 : }
150 : }
151 :
152 3038152 : case QUAD4:
153 : case QUADSHELL4:
154 0 : libmesh_assert_less (o, 4);
155 : libmesh_fallthrough();
156 : case QUAD8:
157 : case QUADSHELL8:
158 : case QUAD9:
159 : case QUADSHELL9:
160 : {
161 3318165 : switch (n)
162 : {
163 : // Vertices
164 129802 : case 0:
165 : case 1:
166 : case 2:
167 : case 3:
168 129802 : return 4;
169 : // Edges
170 118781 : case 4:
171 : case 5:
172 : case 6:
173 : case 7:
174 1405506 : return (2*(o-3));
175 62808 : case 8:
176 : // Interior DoFs are carried on Elems
177 : // return ((o-3)*(o-3));
178 62808 : return 0;
179 :
180 0 : default:
181 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD4/8/9!");
182 : }
183 : }
184 :
185 52614 : case HEX8:
186 : case HEX20:
187 0 : libmesh_assert_less (o, 4);
188 : libmesh_fallthrough();
189 : case HEX27:
190 : {
191 57699 : switch (n)
192 : {
193 : // Vertices
194 1593 : case 0:
195 : case 1:
196 : case 2:
197 : case 3:
198 : case 4:
199 : case 5:
200 : case 6:
201 : case 7:
202 1593 : return 8;
203 : // Edges
204 2169 : case 8:
205 : case 9:
206 : case 10:
207 : case 11:
208 : case 12:
209 : case 13:
210 : case 14:
211 : case 15:
212 : case 16:
213 : case 17:
214 : case 18:
215 : case 19:
216 24840 : return (4*(o-3));
217 : // Faces
218 1134 : case 20:
219 : case 21:
220 : case 22:
221 : case 23:
222 : case 24:
223 : case 25:
224 13014 : return (2*(o-3)*(o-3));
225 378 : case 26:
226 : // Interior DoFs are carried on Elems
227 : // return ((o-3)*(o-3)*(o-3));
228 378 : return 0;
229 :
230 0 : default:
231 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for HEX8/20/27!");
232 : }
233 : }
234 :
235 0 : case INVALID_ELEM:
236 0 : return 0;
237 :
238 0 : default:
239 0 : libmesh_error_msg("ERROR: Bad ElemType = " << Utility::enum_to_string(t) << " for " << Utility::enum_to_string(o) << " order approximation!");
240 : }
241 : } // HERMITE_n_dofs_at_node()
242 :
243 :
244 :
245 1863669 : unsigned int HERMITE_n_dofs_at_node(const Elem & e,
246 : const Order o,
247 : const unsigned int n)
248 : {
249 1969328 : return HERMITE_n_dofs_at_node(e.type(), o, n);
250 : }
251 :
252 :
253 :
254 75214208 : unsigned int HERMITE_n_dofs_per_elem(const ElemType t,
255 : const Order o)
256 : {
257 5155262 : libmesh_assert_greater (o, 2);
258 :
259 75214208 : switch (t)
260 : {
261 0 : case NODEELEM:
262 0 : return 0;
263 5125925 : case EDGE2:
264 : case EDGE3:
265 74867082 : return (o-3);
266 0 : case QUAD4:
267 : case QUADSHELL4:
268 0 : libmesh_assert_less (o, 4);
269 : libmesh_fallthrough();
270 : case QUAD8:
271 : case QUADSHELL8:
272 : case QUAD9:
273 : case QUADSHELL9:
274 344957 : return ((o-3)*(o-3));
275 0 : case HEX8:
276 0 : libmesh_assert_less (o, 4);
277 : libmesh_fallthrough();
278 : case HEX20:
279 : case HEX27:
280 2169 : return ((o-3)*(o-3)*(o-3));
281 :
282 0 : case INVALID_ELEM:
283 0 : return 0;
284 :
285 0 : default:
286 0 : libmesh_error_msg("ERROR: Bad ElemType = " << Utility::enum_to_string(t) << " for " << Utility::enum_to_string(o) << " order approximation!");
287 : }
288 : } // HERMITE_n_dofs_per_elem()
289 :
290 :
291 :
292 70267841 : unsigned int HERMITE_n_dofs_per_elem(const Elem & e,
293 : const Order o)
294 : {
295 75214208 : return HERMITE_n_dofs_per_elem(e.type(), o);
296 : }
297 :
298 :
299 : } // anonymous namespace
300 :
301 :
302 : // Instantiate (side_) nodal_soln() function for every dimension
303 4653679 : LIBMESH_FE_NODAL_SOLN(HERMITE, hermite_nodal_soln)
304 0 : LIBMESH_FE_SIDE_NODAL_SOLN(HERMITE)
305 :
306 :
307 : // Instantiate n_dofs*() functions for every dimension
308 241901302 : LIBMESH_DEFAULT_NDOFS(HERMITE)
309 :
310 :
311 : // Hermite FEMs are C^1 continuous
312 0 : template <> FEContinuity FE<0,HERMITE>::get_continuity() const { return C_ONE; }
313 69333 : template <> FEContinuity FE<1,HERMITE>::get_continuity() const { return C_ONE; }
314 26027 : template <> FEContinuity FE<2,HERMITE>::get_continuity() const { return C_ONE; }
315 3342 : template <> FEContinuity FE<3,HERMITE>::get_continuity() const { return C_ONE; }
316 :
317 : // Hermite FEMs are hierarchic
318 0 : template <> bool FE<0,HERMITE>::is_hierarchic() const { return true; }
319 24 : template <> bool FE<1,HERMITE>::is_hierarchic() const { return true; }
320 24 : template <> bool FE<2,HERMITE>::is_hierarchic() const { return true; }
321 12 : template <> bool FE<3,HERMITE>::is_hierarchic() const { return true; }
322 :
323 :
324 : #ifdef LIBMESH_ENABLE_AMR
325 : // compute_constraints() specializations are only needed for 2 and 3D
326 : template <>
327 18144 : void FE<2,HERMITE>::compute_constraints (DofConstraints & constraints,
328 : DofMap & dof_map,
329 : const unsigned int variable_number,
330 : const Elem * elem)
331 18144 : { compute_proj_constraints(constraints, dof_map, variable_number, elem); }
332 :
333 : template <>
334 108 : void FE<3,HERMITE>::compute_constraints (DofConstraints & constraints,
335 : DofMap & dof_map,
336 : const unsigned int variable_number,
337 : const Elem * elem)
338 108 : { compute_proj_constraints(constraints, dof_map, variable_number, elem); }
339 : #endif // #ifdef LIBMESH_ENABLE_AMR
340 :
341 : // Hermite FEM shapes need reinit
342 0 : template <> bool FE<0,HERMITE>::shapes_need_reinit() const { return true; }
343 84792671 : template <> bool FE<1,HERMITE>::shapes_need_reinit() const { return true; }
344 175314 : template <> bool FE<2,HERMITE>::shapes_need_reinit() const { return true; }
345 36 : template <> bool FE<3,HERMITE>::shapes_need_reinit() const { return true; }
346 :
347 : } // namespace libMesh
|
