libMesh
Loading...
Searching...
No Matches
fe_hermite.C
Go to the documentation of this file.
1// The libMesh Finite Element Library.
2// Copyright (C) 2002-2026 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
27namespace libMesh
28{
29
30// ------------------------------------------------------------
31// Hermite-specific implementations
32
33// Anonymous namespace for local helper functions
34namespace {
35
36void 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 const unsigned int n_nodes = elem->n_nodes();
43
44 const ElemType elem_type = elem->type();
45
46 nodal_soln.resize(n_nodes);
47
48 // FEType object to be passed to various FEInterface functions below.
49 FEType fe_type(order, HERMITE);
50
51 const unsigned int n_sf =
52 FEInterface::n_shape_functions(fe_type, elem, add_p_level);
53
54 std::vector<Point> refspace_nodes;
55 FEBase::get_refspace_nodes(elem_type,refspace_nodes);
56 libmesh_assert_equal_to (refspace_nodes.size(), n_nodes);
57 libmesh_assert_equal_to (elem_soln.size(), n_sf);
58
59 // Zero before summation
60 std::fill(nodal_soln.begin(), nodal_soln.end(), 0);
61
62 for (unsigned int n=0; n<n_nodes; n++)
63 // u_i = Sum (alpha_i phi_i)
64 for (unsigned int i=0; i<n_sf; i++)
65 nodal_soln[n] += elem_soln[i] *
66 FEInterface::shape(fe_type, elem, i, refspace_nodes[n], add_p_level);
67
68} // hermite_nodal_soln()
69
70
71
72unsigned int HERMITE_n_dofs(const ElemType t, const Order o)
73{
74#ifdef DEBUG
75 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 switch (t)
80 {
81 case NODEELEM:
82 return 1;
83 case EDGE2:
84 libmesh_assert_less (o, 4);
85 libmesh_fallthrough();
86 case EDGE3:
87 return (o+1);
88
89 case QUAD4:
90 case QUADSHELL4:
91 case QUAD8:
92 case QUADSHELL8:
93 libmesh_assert_less (o, 4);
94 libmesh_fallthrough();
95 case QUAD9:
96 case QUADSHELL9:
97 return ((o+1)*(o+1));
98
99 case HEX8:
100 case HEX20:
101 libmesh_assert_less (o, 4);
102 libmesh_fallthrough();
103 case HEX27:
104 return ((o+1)*(o+1)*(o+1));
105
106 case INVALID_ELEM:
107 return 0;
108
109 default:
110 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
116unsigned int HERMITE_n_dofs(const Elem * e, const Order o)
117{
118 libmesh_assert(e);
119 return HERMITE_n_dofs(e->type(), o);
120}
121
122
123
124unsigned int HERMITE_n_dofs_at_node(const ElemType t,
125 const Order o,
126 const unsigned int n)
127{
128 libmesh_assert_greater (o, 2);
129 // Piecewise (bi/tri)cubic C1 Hermite splines
130 switch (t)
131 {
132 case NODEELEM:
133 return 1;
134 case EDGE2:
135 case EDGE3:
136 {
137 switch (n)
138 {
139 case 0:
140 case 1:
141 return 2;
142 case 2:
143 // Interior DoFs are carried on Elems
144 // return (o-3);
145 return 0;
146
147 default:
148 libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for EDGE2/3!");
149 }
150 }
151
152 case QUAD4:
153 case QUADSHELL4:
154 libmesh_assert_less (o, 4);
155 libmesh_fallthrough();
156 case QUAD8:
157 case QUADSHELL8:
158 case QUAD9:
159 case QUADSHELL9:
160 {
161 switch (n)
162 {
163 // Vertices
164 case 0:
165 case 1:
166 case 2:
167 case 3:
168 return 4;
169 // Edges
170 case 4:
171 case 5:
172 case 6:
173 case 7:
174 return (2*(o-3));
175 case 8:
176 // Interior DoFs are carried on Elems
177 // return ((o-3)*(o-3));
178 return 0;
179
180 default:
181 libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD4/8/9!");
182 }
183 }
184
185 case HEX8:
186 case HEX20:
187 libmesh_assert_less (o, 4);
188 libmesh_fallthrough();
189 case HEX27:
190 {
191 switch (n)
192 {
193 // Vertices
194 case 0:
195 case 1:
196 case 2:
197 case 3:
198 case 4:
199 case 5:
200 case 6:
201 case 7:
202 return 8;
203 // Edges
204 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 return (4*(o-3));
217 // Faces
218 case 20:
219 case 21:
220 case 22:
221 case 23:
222 case 24:
223 case 25:
224 return (2*(o-3)*(o-3));
225 case 26:
226 // Interior DoFs are carried on Elems
227 // return ((o-3)*(o-3)*(o-3));
228 return 0;
229
230 default:
231 libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for HEX8/20/27!");
232 }
233 }
234
235 case INVALID_ELEM:
236 return 0;
237
238 default:
239 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
245unsigned int HERMITE_n_dofs_at_node(const Elem & e,
246 const Order o,
247 const unsigned int n)
248{
249 return HERMITE_n_dofs_at_node(e.type(), o, n);
250}
251
252
253
254unsigned int HERMITE_n_dofs_per_elem(const ElemType t,
255 const Order o)
256{
257 libmesh_assert_greater (o, 2);
258
259 switch (t)
260 {
261 case NODEELEM:
262 return 0;
263 case EDGE2:
264 case EDGE3:
265 return (o-3);
266 case QUAD4:
267 case QUADSHELL4:
268 libmesh_assert_less (o, 4);
269 libmesh_fallthrough();
270 case QUAD8:
271 case QUADSHELL8:
272 case QUAD9:
273 case QUADSHELL9:
274 return ((o-3)*(o-3));
275 case HEX8:
276 libmesh_assert_less (o, 4);
277 libmesh_fallthrough();
278 case HEX20:
279 case HEX27:
280 return ((o-3)*(o-3)*(o-3));
281
282 case INVALID_ELEM:
283 return 0;
284
285 default:
286 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
292unsigned int HERMITE_n_dofs_per_elem(const Elem & e,
293 const Order o)
294{
295 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
303LIBMESH_FE_NODAL_SOLN(HERMITE, hermite_nodal_soln)
304LIBMESH_FE_SIDE_NODAL_SOLN(HERMITE)
305
306
307// Instantiate n_dofs*() functions for every dimension
308LIBMESH_DEFAULT_NDOFS(HERMITE)
309
310
311// Hermite FEMs are C^1 continuous
312template <> FEContinuity FE<0,HERMITE>::get_continuity() const { return C_ONE; }
313template <> FEContinuity FE<1,HERMITE>::get_continuity() const { return C_ONE; }
314template <> FEContinuity FE<2,HERMITE>::get_continuity() const { return C_ONE; }
315template <> FEContinuity FE<3,HERMITE>::get_continuity() const { return C_ONE; }
316
317// Hermite FEMs are hierarchic
318template <> bool FE<0,HERMITE>::is_hierarchic() const { return true; }
319template <> bool FE<1,HERMITE>::is_hierarchic() const { return true; }
320template <> bool FE<2,HERMITE>::is_hierarchic() const { return true; }
321template <> 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
326template <>
327void FE<2,HERMITE>::compute_constraints (DofConstraints & constraints,
328 DofMap & dof_map,
329 const unsigned int variable_number,
330 const Elem * elem)
331{ compute_proj_constraints(constraints, dof_map, variable_number, elem); }
332
333template <>
334void FE<3,HERMITE>::compute_constraints (DofConstraints & constraints,
335 DofMap & dof_map,
336 const unsigned int variable_number,
337 const Elem * elem)
338{ compute_proj_constraints(constraints, dof_map, variable_number, elem); }
339#endif // #ifdef LIBMESH_ENABLE_AMR
340
341// Hermite FEM shapes need reinit
342template <> bool FE<0,HERMITE>::shapes_need_reinit() const { return true; }
343template <> bool FE<1,HERMITE>::shapes_need_reinit() const { return true; }
344template <> bool FE<2,HERMITE>::shapes_need_reinit() const { return true; }
345template <> bool FE<3,HERMITE>::shapes_need_reinit() const { return true; }
346
347} // namespace libMesh
libMesh::FEAbstract::get_refspace_nodes
static void get_refspace_nodes(const ElemType t, std::vector< Point > &nodes)
Definition fe_abstract.C:400
libMesh::FEInterface::shape
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
Definition fe_interface.C:761
libMesh::FEInterface::n_shape_functions
static unsigned int n_shape_functions(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition fe_interface.C:272
libMesh
The libMesh namespace provides an interface to certain functionality in the library.
libMesh::ElemType
ElemType
Defines an enum for geometric element types.
Definition enum_elem_type.h:33
n_nodes
const dof_id_type n_nodes
Definition tecplot_io.C:67
Generated on Fri Aug 21 2026 13:32:28 for libMesh by doxygen 1.9.8