libMesh
fe_hermite_shape_3D.C
Go to the documentation of this file.
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 // Local includes
20 #include "libmesh/fe.h"
21 #include "libmesh/elem.h"
22 #include "libmesh/fe_interface.h"
23 #include "libmesh/number_lookups.h"
24 #include "libmesh/enum_to_string.h"
25 
26 namespace
27 {
28 using namespace libMesh;
29 
30 // Compute the static coefficients for an element
31 void hermite_compute_coefs(const Elem * elem, std::vector<std::vector<Real>> & dxdxi
32 
33 #ifdef DEBUG
34  , std::vector<Real> & dydxi, std::vector<Real> & dzdeta, std::vector<Real> & dxdzeta,
35  std::vector<Real> & dzdxi, std::vector<Real> & dxdeta, std::vector<Real> & dydzeta
36 #endif
37  )
38 {
39  const FEFamily mapping_family = FEMap::map_fe_type(*elem);
40  const FEType map_fe_type(elem->default_order(), mapping_family);
41 
42  // Note: we explicitly don't consider the elem->p_level() when
43  // computing the number of mapping shape functions.
44  const int n_mapping_shape_functions =
45  FEInterface::n_shape_functions (map_fe_type, /*extra_order=*/0, elem);
46 
47  static const Point dofpt[2] = {Point(-1,-1,-1), Point(1,1,1)};
48 
49  FEInterface::shape_deriv_ptr shape_deriv_ptr =
50  FEInterface::shape_deriv_function(map_fe_type, elem);
51 
52  for (int p = 0; p != 2; ++p)
53  {
54  dxdxi[0][p] = 0;
55  dxdxi[1][p] = 0;
56  dxdxi[2][p] = 0;
57 #ifdef DEBUG
58  dydxi[p] = 0;
59  dzdeta[p] = 0;
60  dxdzeta[p] = 0;
61  dzdxi[p] = 0;
62  dxdeta[p] = 0;
63  dydzeta[p] = 0;
64 #endif
65  for (int i = 0; i != n_mapping_shape_functions; ++i)
66  {
67  const Real ddxi = shape_deriv_ptr
68  (map_fe_type, elem, i, 0, dofpt[p], /*add_p_level=*/false);
69  const Real ddeta = shape_deriv_ptr
70  (map_fe_type, elem, i, 1, dofpt[p], /*add_p_level=*/false);
71  const Real ddzeta = shape_deriv_ptr
72  (map_fe_type, elem, i, 2, dofpt[p], /*add_p_level=*/false);
73 
74  // dxdeta, dxdzeta, dydxi, dydzeta, dzdxi, dzdeta should all
75  // be 0!
76  const Point & point_i = elem->point(i);
77  dxdxi[0][p] += point_i(0) * ddxi;
78  dxdxi[1][p] += point_i(1) * ddeta;
79  dxdxi[2][p] += point_i(2) * ddzeta;
80 #ifdef DEBUG
81  dydxi[p] += point_i(1) * ddxi;
82  dzdeta[p] += point_i(2) * ddeta;
83  dxdzeta[p] += point_i(0) * ddzeta;
84  dzdxi[p] += point_i(2) * ddxi;
85  dxdeta[p] += point_i(0) * ddeta;
86  dydzeta[p] += point_i(1) * ddzeta;
87 #endif
88  }
89 
90  // No singular elements!
91  libmesh_assert(dxdxi[0][p]);
92  libmesh_assert(dxdxi[1][p]);
93  libmesh_assert(dxdxi[2][p]);
94  // No non-rectilinear or non-axis-aligned elements!
95 #ifdef DEBUG
96  libmesh_assert_less (std::abs(dydxi[p]), TOLERANCE);
97  libmesh_assert_less (std::abs(dzdeta[p]), TOLERANCE);
98  libmesh_assert_less (std::abs(dxdzeta[p]), TOLERANCE);
99  libmesh_assert_less (std::abs(dzdxi[p]), TOLERANCE);
100  libmesh_assert_less (std::abs(dxdeta[p]), TOLERANCE);
101  libmesh_assert_less (std::abs(dydzeta[p]), TOLERANCE);
102 #endif
103  }
104 }
105 
106 
107 
108 Real hermite_bases_3D (std::vector<unsigned int> & bases1D,
109  const std::vector<std::vector<Real>> & dxdxi,
110  const Order & o,
111  unsigned int i)
112 {
113  bases1D.clear();
114  bases1D.resize(3,0);
115  Real coef = 1.0;
116 
117  unsigned int e = o-2;
118 
119  // Nodes
120  if (i < 64)
121  {
122  switch (i / 8)
123  {
124  case 0:
125  break;
126  case 1:
127  bases1D[0] = 1;
128  break;
129  case 2:
130  bases1D[0] = 1;
131  bases1D[1] = 1;
132  break;
133  case 3:
134  bases1D[1] = 1;
135  break;
136  case 4:
137  bases1D[2] = 1;
138  break;
139  case 5:
140  bases1D[0] = 1;
141  bases1D[2] = 1;
142  break;
143  case 6:
144  bases1D[0] = 1;
145  bases1D[1] = 1;
146  bases1D[2] = 1;
147  break;
148  case 7:
149  bases1D[1] = 1;
150  bases1D[2] = 1;
151  break;
152  default:
153  libmesh_error_msg("Invalid basis node = " << i/8);
154  }
155 
156  unsigned int basisnum = i%8;
157  switch (basisnum) // DoF type
158  {
159  case 0: // DoF = value at node
160  coef = 1.0;
161  break;
162  case 1: // DoF = x derivative at node
163  coef = dxdxi[0][bases1D[0]];
164  bases1D[0] += 2; break;
165  case 2: // DoF = y derivative at node
166  coef = dxdxi[1][bases1D[1]];
167  bases1D[1] += 2; break;
168  case 3: // DoF = xy derivative at node
169  coef = dxdxi[0][bases1D[0]] * dxdxi[1][bases1D[1]];
170  bases1D[0] += 2; bases1D[1] += 2; break;
171  case 4: // DoF = z derivative at node
172  coef = dxdxi[2][bases1D[2]];
173  bases1D[2] += 2; break;
174  case 5: // DoF = xz derivative at node
175  coef = dxdxi[0][bases1D[0]] * dxdxi[2][bases1D[2]];
176  bases1D[0] += 2; bases1D[2] += 2; break;
177  case 6: // DoF = yz derivative at node
178  coef = dxdxi[1][bases1D[1]] * dxdxi[2][bases1D[2]];
179  bases1D[1] += 2; bases1D[2] += 2; break;
180  case 7: // DoF = xyz derivative at node
181  coef = dxdxi[0][bases1D[0]] * dxdxi[1][bases1D[1]] * dxdxi[2][bases1D[2]];
182  bases1D[0] += 2; bases1D[1] += 2; bases1D[2] += 2; break;
183  default:
184  libmesh_error_msg("Invalid basisnum = " << basisnum);
185  }
186  }
187  // Edges
188  else if (i < 64 + 12*4*e)
189  {
190  unsigned int basisnum = (i - 64) % (4*e);
191  switch ((i - 64) / (4*e))
192  {
193  case 0:
194  bases1D[0] = basisnum / 4 + 4;
195  bases1D[1] = basisnum % 4 / 2 * 2;
196  bases1D[2] = basisnum % 2 * 2;
197  if (basisnum % 4 / 2)
198  coef *= dxdxi[1][0];
199  if (basisnum % 2)
200  coef *= dxdxi[2][0];
201  break;
202  case 1:
203  bases1D[0] = basisnum % 4 / 2 * 2 + 1;
204  bases1D[1] = basisnum / 4 + 4;
205  bases1D[2] = basisnum % 2 * 2;
206  if (basisnum % 4 / 2)
207  coef *= dxdxi[0][1];
208  if (basisnum % 2)
209  coef *= dxdxi[2][0];
210  break;
211  case 2:
212  bases1D[0] = basisnum / 4 + 4;
213  bases1D[1] = basisnum % 4 / 2 * 2 + 1;
214  bases1D[2] = basisnum % 2 * 2;
215  if (basisnum % 4 / 2)
216  coef *= dxdxi[1][1];
217  if (basisnum % 2)
218  coef *= dxdxi[2][0];
219  break;
220  case 3:
221  bases1D[0] = basisnum % 4 / 2 * 2;
222  bases1D[1] = basisnum / 4 + 4;
223  bases1D[2] = basisnum % 2 * 2;
224  if (basisnum % 4 / 2)
225  coef *= dxdxi[0][0];
226  if (basisnum % 2)
227  coef *= dxdxi[2][0];
228  break;
229  case 4:
230  bases1D[0] = basisnum % 4 / 2 * 2;
231  bases1D[1] = basisnum % 2 * 2;
232  bases1D[2] = basisnum / 4 + 4;
233  if (basisnum % 4 / 2)
234  coef *= dxdxi[0][0];
235  if (basisnum % 2)
236  coef *= dxdxi[1][0];
237  break;
238  case 5:
239  bases1D[0] = basisnum % 4 / 2 * 2 + 1;
240  bases1D[1] = basisnum % 2 * 2;
241  bases1D[2] = basisnum / 4 + 4;
242  if (basisnum % 4 / 2)
243  coef *= dxdxi[0][1];
244  if (basisnum % 2)
245  coef *= dxdxi[1][0];
246  break;
247  case 6:
248  bases1D[0] = basisnum % 4 / 2 * 2 + 1;
249  bases1D[1] = basisnum % 2 * 2 + 1;
250  bases1D[2] = basisnum / 4 + 4;
251  if (basisnum % 4 / 2)
252  coef *= dxdxi[0][1];
253  if (basisnum % 2)
254  coef *= dxdxi[1][1];
255  break;
256  case 7:
257  bases1D[0] = basisnum % 4 / 2 * 2;
258  bases1D[1] = basisnum % 2 * 2 + 1;
259  bases1D[2] = basisnum / 4 + 4;
260  if (basisnum % 4 / 2)
261  coef *= dxdxi[0][0];
262  if (basisnum % 2)
263  coef *= dxdxi[1][1];
264  break;
265  case 8:
266  bases1D[0] = basisnum / 4 + 4;
267  bases1D[1] = basisnum % 4 / 2 * 2;
268  bases1D[2] = basisnum % 2 * 2 + 1;
269  if (basisnum % 4 / 2)
270  coef *= dxdxi[1][0];
271  if (basisnum % 2)
272  coef *= dxdxi[2][1];
273  break;
274  case 9:
275  bases1D[0] = basisnum % 4 / 2 * 2 + 1;
276  bases1D[1] = basisnum / 4 + 4;
277  bases1D[2] = basisnum % 2 * 2;
278  if (basisnum % 4 / 2)
279  coef *= dxdxi[0][1];
280  if (basisnum % 2)
281  coef *= dxdxi[2][1];
282  break;
283  case 10:
284  bases1D[0] = basisnum / 4 + 4;
285  bases1D[1] = basisnum % 4 / 2 * 2 + 1;
286  bases1D[2] = basisnum % 2 * 2 + 1;
287  if (basisnum % 4 / 2)
288  coef *= dxdxi[1][1];
289  if (basisnum % 2)
290  coef *= dxdxi[2][1];
291  break;
292  case 11:
293  bases1D[0] = basisnum % 4 / 2 * 2;
294  bases1D[1] = basisnum / 4 + 4;
295  bases1D[2] = basisnum % 2 * 2 + 1;
296  if (basisnum % 4 / 2)
297  coef *= dxdxi[0][0];
298  if (basisnum % 2)
299  coef *= dxdxi[2][1];
300  break;
301  default:
302  libmesh_error_msg("Invalid basis node = " << (i - 64) / (4*e));
303  }
304  }
305  // Faces
306  else if (i < 64 + 12*4*e + 6*2*e*e)
307  {
308  unsigned int basisnum = (i - 64 - 12*4*e) % (2*e*e);
309  switch ((i - 64 - 12*4*e) / (2*e*e))
310  {
311  case 0:
312  bases1D[0] = square_number_column[basisnum / 2];
313  bases1D[1] = square_number_row[basisnum / 2];
314  bases1D[2] = basisnum % 2 * 2;
315  if (basisnum % 2)
316  coef *= dxdxi[2][0];
317  break;
318  case 1:
319  bases1D[0] = square_number_column[basisnum / 2];
320  bases1D[1] = basisnum % 2 * 2;
321  bases1D[2] = square_number_row[basisnum / 2];
322  if (basisnum % 2)
323  coef *= dxdxi[1][0];
324  break;
325  case 2:
326  bases1D[0] = basisnum % 2 * 2 + 1;
327  bases1D[1] = square_number_column[basisnum / 2];
328  bases1D[2] = square_number_row[basisnum / 2];
329  if (basisnum % 2)
330  coef *= dxdxi[0][1];
331  break;
332  case 3:
333  bases1D[0] = square_number_column[basisnum / 2];
334  bases1D[1] = basisnum % 2 * 2 + 1;
335  bases1D[2] = square_number_row[basisnum / 2];
336  if (basisnum % 2)
337  coef *= dxdxi[1][1];
338  break;
339  case 4:
340  bases1D[0] = basisnum % 2 * 2;
341  bases1D[1] = square_number_column[basisnum / 2];
342  bases1D[2] = square_number_row[basisnum / 2];
343  if (basisnum % 2)
344  coef *= dxdxi[0][0];
345  break;
346  case 5:
347  bases1D[0] = square_number_column[basisnum / 2];
348  bases1D[1] = square_number_row[basisnum / 2];
349  bases1D[2] = basisnum % 2 * 2 + 1;
350  if (basisnum % 2)
351  coef *= dxdxi[2][1];
352  break;
353  default:
354  libmesh_error_msg("Invalid basis node = " << (i - 64 - 12*4*e) / (2*e*e));
355  }
356  }
357  // Interior
358  else
359  {
360  unsigned int basisnum = i - 64 - 12*4*e;
361  bases1D[0] = cube_number_column[basisnum] + 4;
362  bases1D[1] = cube_number_row[basisnum] + 4;
363  bases1D[2] = cube_number_page[basisnum] + 4;
364  }
365 
366  // No singular elements
367  libmesh_assert(coef);
368  return coef;
369 }
370 
371 
372 } // end anonymous namespace
373 
374 
375 namespace libMesh
376 {
377 
378 
379 LIBMESH_DEFAULT_VECTORIZED_FE(3,HERMITE)
380 
381 
382 template <>
383 Real FE<3,HERMITE>::shape(const Elem * elem,
384  const Order order,
385  const unsigned int i,
386  const Point & p,
387  const bool add_p_level)
388 {
389  libmesh_assert(elem);
390 
391  std::vector<std::vector<Real>> dxdxi(3, std::vector<Real>(2, 0));
392 
393 #ifdef DEBUG
394  std::vector<Real> dydxi(2), dzdeta(2), dxdzeta(2);
395  std::vector<Real> dzdxi(2), dxdeta(2), dydzeta(2);
396 #endif //DEBUG
397 
398  hermite_compute_coefs(elem, dxdxi
399 #ifdef DEBUG
400  , dydxi, dzdeta, dxdzeta, dzdxi, dxdeta, dydzeta
401 #endif
402  );
403 
404  const ElemType type = elem->type();
405 
406  const Order totalorder =
407  order + add_p_level*elem->p_level();
408 
409  switch (totalorder)
410  {
411  // 3rd-order tricubic Hermite functions
412  case THIRD:
413  {
414  switch (type)
415  {
416  case HEX8:
417  case HEX20:
418  case HEX27:
419  {
420  libmesh_assert_less (i, 64);
421 
422  std::vector<unsigned int> bases1D;
423 
424  Real coef = hermite_bases_3D(bases1D, dxdxi, totalorder, i);
425 
426  return coef *
427  FEHermite<1>::hermite_raw_shape(bases1D[0],p(0)) *
428  FEHermite<1>::hermite_raw_shape(bases1D[1],p(1)) *
429  FEHermite<1>::hermite_raw_shape(bases1D[2],p(2));
430  }
431  default:
432  libmesh_error_msg("ERROR: Unsupported element type " << Utility::enum_to_string(type));
433  }
434  }
435  // by default throw an error
436  default:
437  libmesh_error_msg("ERROR: Unsupported polynomial order " << totalorder);
438  }
439 }
440 
441 
442 template <>
443 Real FE<3,HERMITE>::shape(const ElemType,
444  const Order,
445  const unsigned int,
446  const Point &)
447 {
448  libmesh_error_msg("Hermite elements require the real element \nto construct gradient-based degrees of freedom.");
449  return 0.;
450 }
451 
452 template <>
453 Real FE<3,HERMITE>::shape(const FEType fet,
454  const Elem * elem,
455  const unsigned int i,
456  const Point & p,
457  const bool add_p_level)
458 {
459  return FE<3,HERMITE>::shape(elem, fet.order, i, p, add_p_level);
460 }
461 
462 
463 
464 template <>
465 Real FE<3,HERMITE>::shape_deriv(const Elem * elem,
466  const Order order,
467  const unsigned int i,
468  const unsigned int j,
469  const Point & p,
470  const bool add_p_level)
471 {
472  libmesh_assert(elem);
473  libmesh_assert (j == 0 || j == 1 || j == 2);
474 
475  std::vector<std::vector<Real>> dxdxi(3, std::vector<Real>(2, 0));
476 
477 #ifdef DEBUG
478  std::vector<Real> dydxi(2), dzdeta(2), dxdzeta(2);
479  std::vector<Real> dzdxi(2), dxdeta(2), dydzeta(2);
480 #endif //DEBUG
481 
482  hermite_compute_coefs(elem, dxdxi
483 #ifdef DEBUG
484  , dydxi, dzdeta, dxdzeta, dzdxi, dxdeta, dydzeta
485 #endif
486  );
487 
488  const ElemType type = elem->type();
489 
490  const Order totalorder =
491  order + add_p_level*elem->p_level();
492 
493  switch (totalorder)
494  {
495  // 3rd-order tricubic Hermite functions
496  case THIRD:
497  {
498  switch (type)
499  {
500  case HEX8:
501  case HEX20:
502  case HEX27:
503  {
504  libmesh_assert_less (i, 64);
505 
506  std::vector<unsigned int> bases1D;
507 
508  Real coef = hermite_bases_3D(bases1D, dxdxi, totalorder, i);
509 
510  switch (j) // Derivative type
511  {
512  case 0:
513  return coef *
514  FEHermite<1>::hermite_raw_shape_deriv(bases1D[0],p(0)) *
515  FEHermite<1>::hermite_raw_shape(bases1D[1],p(1)) *
516  FEHermite<1>::hermite_raw_shape(bases1D[2],p(2));
517  break;
518  case 1:
519  return coef *
520  FEHermite<1>::hermite_raw_shape(bases1D[0],p(0)) *
521  FEHermite<1>::hermite_raw_shape_deriv(bases1D[1],p(1)) *
522  FEHermite<1>::hermite_raw_shape(bases1D[2],p(2));
523  break;
524  case 2:
525  return coef *
526  FEHermite<1>::hermite_raw_shape(bases1D[0],p(0)) *
527  FEHermite<1>::hermite_raw_shape(bases1D[1],p(1)) *
528  FEHermite<1>::hermite_raw_shape_deriv(bases1D[2],p(2));
529  break;
530  default:
531  libmesh_error_msg("Invalid shape function derivative j = " << j);
532  }
533 
534  }
535  default:
536  libmesh_error_msg("ERROR: Unsupported element type " << Utility::enum_to_string(type));
537  }
538  }
539  // by default throw an error
540  default:
541  libmesh_error_msg("ERROR: Unsupported polynomial order " << totalorder);
542  }
543 }
544 
545 
546 template <>
547 Real FE<3,HERMITE>::shape_deriv(const ElemType,
548  const Order,
549  const unsigned int,
550  const unsigned int,
551  const Point &)
552 {
553  libmesh_error_msg("Hermite elements require the real element \nto construct gradient-based degrees of freedom.");
554  return 0.;
555 }
556 
557 
558 template <>
559 Real FE<3,HERMITE>::shape_deriv(const FEType fet,
560  const Elem * elem,
561  const unsigned int i,
562  const unsigned int j,
563  const Point & p,
564  const bool add_p_level)
565 {
566  return FE<3,HERMITE>::shape_deriv(elem, fet.order, i, j, p, add_p_level);
567 }
568 
569 
570 
571 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
572 
573 
574 template <>
575 Real FE<3,HERMITE>::shape_second_deriv(const Elem * elem,
576  const Order order,
577  const unsigned int i,
578  const unsigned int j,
579  const Point & p,
580  const bool add_p_level)
581 {
582  libmesh_assert(elem);
583 
584  std::vector<std::vector<Real>> dxdxi(3, std::vector<Real>(2, 0));
585 
586 #ifdef DEBUG
587  std::vector<Real> dydxi(2), dzdeta(2), dxdzeta(2);
588  std::vector<Real> dzdxi(2), dxdeta(2), dydzeta(2);
589 #endif //DEBUG
590 
591  hermite_compute_coefs(elem, dxdxi
592 #ifdef DEBUG
593  , dydxi, dzdeta, dxdzeta, dzdxi, dxdeta, dydzeta
594 #endif
595  );
596 
597  const ElemType type = elem->type();
598 
599  const Order totalorder =
600  order + add_p_level*elem->p_level();
601 
602  switch (totalorder)
603  {
604  // 3rd-order tricubic Hermite functions
605  case THIRD:
606  {
607  switch (type)
608  {
609  case HEX8:
610  case HEX20:
611  case HEX27:
612  {
613  libmesh_assert_less (i, 64);
614 
615  std::vector<unsigned int> bases1D;
616 
617  Real coef = hermite_bases_3D(bases1D, dxdxi, totalorder, i);
618 
619  switch (j) // Derivative type
620  {
621  case 0:
622  return coef *
623  FEHermite<1>::hermite_raw_shape_second_deriv(bases1D[0],p(0)) *
624  FEHermite<1>::hermite_raw_shape(bases1D[1],p(1)) *
625  FEHermite<1>::hermite_raw_shape(bases1D[2],p(2));
626  break;
627  case 1:
628  return coef *
629  FEHermite<1>::hermite_raw_shape_deriv(bases1D[0],p(0)) *
630  FEHermite<1>::hermite_raw_shape_deriv(bases1D[1],p(1)) *
631  FEHermite<1>::hermite_raw_shape(bases1D[2],p(2));
632  break;
633  case 2:
634  return coef *
635  FEHermite<1>::hermite_raw_shape(bases1D[0],p(0)) *
636  FEHermite<1>::hermite_raw_shape_second_deriv(bases1D[1],p(1)) *
637  FEHermite<1>::hermite_raw_shape(bases1D[2],p(2));
638  break;
639  case 3:
640  return coef *
641  FEHermite<1>::hermite_raw_shape_deriv(bases1D[0],p(0)) *
642  FEHermite<1>::hermite_raw_shape(bases1D[1],p(1)) *
643  FEHermite<1>::hermite_raw_shape_deriv(bases1D[2],p(2));
644  break;
645  case 4:
646  return coef *
647  FEHermite<1>::hermite_raw_shape(bases1D[0],p(0)) *
648  FEHermite<1>::hermite_raw_shape_deriv(bases1D[1],p(1)) *
649  FEHermite<1>::hermite_raw_shape_deriv(bases1D[2],p(2));
650  break;
651  case 5:
652  return coef *
653  FEHermite<1>::hermite_raw_shape(bases1D[0],p(0)) *
654  FEHermite<1>::hermite_raw_shape(bases1D[1],p(1)) *
655  FEHermite<1>::hermite_raw_shape_second_deriv(bases1D[2],p(2));
656  break;
657  default:
658  libmesh_error_msg("Invalid shape function derivative j = " << j);
659  }
660 
661  }
662  default:
663  libmesh_error_msg("ERROR: Unsupported element type " << Utility::enum_to_string(type));
664  }
665  }
666  // by default throw an error
667  default:
668  libmesh_error_msg("ERROR: Unsupported polynomial order " << totalorder);
669  }
670 }
671 
672 
673 template <>
674 Real FE<3,HERMITE>::shape_second_deriv(const ElemType,
675  const Order,
676  const unsigned int,
677  const unsigned int,
678  const Point &)
679 {
680  libmesh_error_msg("Hermite elements require the real element \nto construct gradient-based degrees of freedom.");
681  return 0.;
682 }
683 
684 template <>
685 Real FE<3,HERMITE>::shape_second_deriv(const FEType fet,
686  const Elem * elem,
687  const unsigned int i,
688  const unsigned int j,
689  const Point & p,
690  const bool add_p_level)
691 {
692  return FE<3,HERMITE>::shape_second_deriv(elem, fet.order, i, j, p, add_p_level);
693 }
694 
695 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
696 
697 } // namespace libMesh
libMesh::FEInterface::shape_deriv_ptr
Real(* shape_deriv_ptr)(const FEType fet, const Elem *elem, const unsigned int i, const unsigned int j, const Point &p, const bool add_p_level)
Typedef for pointer to a function that returns FE shape function derivative values.
Definition: fe_interface.h:648
libMesh::FEType
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:196
libMesh::ElemType
ElemType
Defines an enum for geometric element types.
Definition: enum_elem_type.h:33
libMesh::Order
Order
defines an enum for polynomial orders.
Definition: enum_order.h:40
libMesh::FE::shape
static OutputShape shape(const ElemType t, const Order o, const unsigned int i, const Point &p)
libMesh::TOLERANCE
static constexpr Real TOLERANCE
Definition: libmesh_common.h:151
libMesh::Elem
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
libMesh::cube_number_column
const unsigned char cube_number_column[]
Definition: number_lookups.C:84
libMesh::FE::shape_deriv
static OutputShape shape_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
libMesh::Elem::p_level
unsigned int p_level() const
Definition: elem.h:3108
libMesh::FEType::order
OrderWrapper order
The approximation order of the element.
Definition: fe_type.h:215
libMesh
The libMesh namespace provides an interface to certain functionality in the library.
libMesh::FEHermite::hermite_raw_shape_second_deriv
static Real hermite_raw_shape_second_deriv(const unsigned int basis_num, const Real xi)
1D hermite functions on unit interval
libMesh::FEHermite::hermite_raw_shape
static Real hermite_raw_shape(const unsigned int basis_num, const Real xi)
libMesh::cube_number_row
const unsigned char cube_number_row[]
Definition: number_lookups.C:196
libMesh::square_number_column
const unsigned char square_number_column[]
Definition: number_lookups.C:56
libMesh::cube_number_page
const unsigned char cube_number_page[]
Definition: number_lookups.C:308
libMesh::LIBMESH_DEFAULT_VECTORIZED_FE
LIBMESH_DEFAULT_VECTORIZED_FE(template<>Real FE< 0, BERNSTEIN)
Definition: fe_bernstein_shape_0D.C:30
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::libmesh_assert
libmesh_assert(ctx)
libMesh::FEInterface::shape_deriv_function
static shape_deriv_ptr shape_deriv_function(const unsigned int dim, const FEType &fe_t, const ElemType t)
libMesh::Utility::enum_to_string
std::string enum_to_string(const T e)
libMesh::Real
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Definition: libmesh_common.h:144
libMesh::square_number_row
const unsigned char square_number_row[]
Definition: number_lookups.C:69
libMesh::FEHermite::hermite_raw_shape_deriv
static Real hermite_raw_shape_deriv(const unsigned int basis_num, const Real xi)
libMesh::FEFamily
FEFamily
defines an enum for finite element families.
Definition: enum_fe_family.h:34
libMesh::Elem::default_order
virtual Order default_order() const =0
libMesh::Elem::type
virtual ElemType type() const =0
libMesh::Point
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39
libMesh::Elem::point
const Point & point(const unsigned int i) const
Definition: elem.h:2453
libMesh::FE::shape_second_deriv
static OutputShape shape_second_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
libMesh::FEMap::map_fe_type
static FEFamily map_fe_type(const Elem &elem)
Definition: fe_map.C:46
Generated on Thu Jul 17 2025 01:29:05 for libMesh by   doxygen 1.8.14