libMesh
fe_bernstein_shape_2D.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2019 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/libmesh_config.h"
21 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
22 
23 #include "libmesh/fe.h"
24 #include "libmesh/elem.h"
25 #include "libmesh/number_lookups.h"
26 #include "libmesh/utility.h"
27 
28 
29 namespace libMesh
30 {
31 
32 
33 template <>
35  const Order,
36  const unsigned int,
37  const Point &)
38 {
39  libmesh_error_msg("Bernstein polynomials require the element type \nbecause edge orientation is needed.");
40  return 0.;
41 }
42 
43 
44 
45 template <>
47  const Order order,
48  const unsigned int i,
49  const Point & p,
50  const bool add_p_level)
51 {
52  libmesh_assert(elem);
53 
54  const ElemType type = elem->type();
55 
56  const Order totalorder =
57  static_cast<Order>(order + add_p_level * elem->p_level());
58 
59  // Declare that we are using our own special power function
60  // from the Utility namespace. This saves typing later.
61  using Utility::pow;
62 
63  switch (type)
64  {
65  // Hierarchic shape functions on the quadrilateral.
66  case QUAD4:
67  case QUADSHELL4:
68  case QUAD9:
69  {
70  // Compute quad shape functions as a tensor-product
71  const Real xi = p(0);
72  const Real eta = p(1);
73 
74  libmesh_assert_less (i, (totalorder+1u)*(totalorder+1u));
75 
76  // Example i, i0, i1 values for totalorder = 5:
77  // 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35
78  // static const unsigned int i0[] = {0, 1, 1, 0, 2, 3, 4, 5, 1, 1, 1, 1, 2, 3, 4, 5, 0, 0, 0, 0, 2, 3, 3, 2, 4, 4, 4, 3, 2, 5, 5, 5, 5, 4, 3, 2};
79  // static const unsigned int i1[] = {0, 0, 1, 1, 0, 0, 0, 0, 2, 3, 4, 5, 1, 1, 1, 1, 2, 3, 4, 5, 2, 2, 3, 3, 2, 3, 4, 4, 4, 2, 3, 4, 5, 5, 5, 5};
80 
81  unsigned int i0, i1;
82 
83  // Vertex DoFs
84  if (i == 0)
85  { i0 = 0; i1 = 0; }
86  else if (i == 1)
87  { i0 = 1; i1 = 0; }
88  else if (i == 2)
89  { i0 = 1; i1 = 1; }
90  else if (i == 3)
91  { i0 = 0; i1 = 1; }
92 
93 
94  // Edge DoFs
95  else if (i < totalorder + 3u)
96  { i0 = i - 2; i1 = 0; }
97  else if (i < 2u*totalorder + 2)
98  { i0 = 1; i1 = i - totalorder - 1; }
99  else if (i < 3u*totalorder + 1)
100  { i0 = i - 2u*totalorder; i1 = 1; }
101  else if (i < 4u*totalorder)
102  { i0 = 0; i1 = i - 3u*totalorder + 1; }
103  // Interior DoFs. Use Roy's number look up
104  else
105  {
106  unsigned int basisnum = i - 4*totalorder;
107  i0 = square_number_column[basisnum] + 2;
108  i1 = square_number_row[basisnum] + 2;
109  }
110 
111 
112  // Flip odd degree of freedom values if necessary
113  // to keep continuity on sides.
114  if ((i>= 4 && i<= 4+ totalorder-2u) && elem->point(0) > elem->point(1)) i0=totalorder+2-i0;//
115  else if ((i>= 4+ totalorder-1u && i<= 4+2*totalorder-3u) && elem->point(1) > elem->point(2)) i1=totalorder+2-i1;
116  else if ((i>= 4+2*totalorder-2u && i<= 4+3*totalorder-4u) && elem->point(3) > elem->point(2)) i0=totalorder+2-i0;
117  else if ((i>= 4+3*totalorder-3u && i<= 4+4*totalorder-5u) && elem->point(0) > elem->point(3)) i1=totalorder+2-i1;
118 
119 
120  return (FE<1,BERNSTEIN>::shape(EDGE3, totalorder, i0, xi)*
121  FE<1,BERNSTEIN>::shape(EDGE3, totalorder, i1, eta));
122  }
123  // handle serendipity QUAD8 element separately
124  case QUAD8:
125  case QUADSHELL8:
126  {
127  libmesh_assert_less (totalorder, 3);
128 
129  const Real xi = p(0);
130  const Real eta = p(1);
131 
132  libmesh_assert_less (i, 8);
133 
134  // 0 1 2 3 4 5 6 7 8
135  static const unsigned int i0[] = {0, 1, 1, 0, 2, 1, 2, 0, 2};
136  static const unsigned int i1[] = {0, 0, 1, 1, 0, 2, 1, 2, 2};
137  static const Real scal[] = {-0.25, -0.25, -0.25, -0.25, 0.5, 0.5, 0.5, 0.5};
138 
139  //B_t,i0(i)|xi * B_s,i1(i)|eta + scal(i) * B_t,2|xi * B_t,2|eta
140  return (FE<1,BERNSTEIN>::shape(EDGE3, totalorder, i0[i], xi)*
141  FE<1,BERNSTEIN>::shape(EDGE3, totalorder, i1[i], eta)
142  +scal[i]*
143  FE<1,BERNSTEIN>::shape(EDGE3, totalorder, i0[8], xi)*
144  FE<1,BERNSTEIN>::shape(EDGE3, totalorder, i1[8], eta));
145 
146  }
147 
148  case TRI3:
149  case TRISHELL3:
150  libmesh_assert_less (totalorder, 2);
151  libmesh_fallthrough();
152  case TRI6:
153  switch (totalorder)
154  {
155  case FIRST:
156  {
157  const Real x=p(0);
158  const Real y=p(1);
159  const Real r=1.-x-y;
160 
161  libmesh_assert_less (i, 3);
162 
163  switch(i)
164  {
165  case 0: return r; //f0,0,1
166  case 1: return x; //f0,1,1
167  case 2: return y; //f1,0,1
168 
169  default:
170  libmesh_error_msg("Invalid shape function index i = " << i);
171  }
172  }
173  case SECOND:
174  {
175  const Real x=p(0);
176  const Real y=p(1);
177  const Real r=1.-x-y;
178 
179  libmesh_assert_less (i, 6);
180 
181  switch(i)
182  {
183  case 0: return r*r;
184  case 1: return x*x;
185  case 2: return y*y;
186 
187  case 3: return 2.*x*r;
188  case 4: return 2.*x*y;
189  case 5: return 2.*r*y;
190 
191  default:
192  libmesh_error_msg("Invalid shape function index i = " << i);
193  }
194  }
195  case THIRD:
196  {
197  const Real x=p(0);
198  const Real y=p(1);
199  const Real r=1.-x-y;
200  libmesh_assert_less (i, 10);
201 
202  unsigned int shape=i;
203 
204 
205  if ((i==3||i==4) && elem->point(0) > elem->point(1)) shape=7-i;
206  if ((i==5||i==6) && elem->point(1) > elem->point(2)) shape=11-i;
207  if ((i==7||i==8) && elem->point(0) > elem->point(2)) shape=15-i;
208 
209  switch(shape)
210  {
211  case 0: return r*r*r;
212  case 1: return x*x*x;
213  case 2: return y*y*y;
214 
215  case 3: return 3.*x*r*r;
216  case 4: return 3.*x*x*r;
217 
218  case 5: return 3.*y*x*x;
219  case 6: return 3.*y*y*x;
220 
221  case 7: return 3.*y*r*r;
222  case 8: return 3.*y*y*r;
223 
224  case 9: return 6.*x*y*r;
225 
226  default:
227  libmesh_error_msg("Invalid shape function index shape = " << shape);
228  }
229  }
230  case FOURTH:
231  {
232  const Real x=p(0);
233  const Real y=p(1);
234  const Real r=1-x-y;
235  unsigned int shape=i;
236 
237  libmesh_assert_less (i, 15);
238 
239  if ((i==3||i== 5) && elem->point(0) > elem->point(1)) shape=8-i;
240  if ((i==6||i== 8) && elem->point(1) > elem->point(2)) shape=14-i;
241  if ((i==9||i==11) && elem->point(0) > elem->point(2)) shape=20-i;
242 
243 
244  switch(shape)
245  {
246  // point functions
247  case 0: return r*r*r*r;
248  case 1: return x*x*x*x;
249  case 2: return y*y*y*y;
250 
251  // edge functions
252  case 3: return 4.*x*r*r*r;
253  case 4: return 6.*x*x*r*r;
254  case 5: return 4.*x*x*x*r;
255 
256  case 6: return 4.*y*x*x*x;
257  case 7: return 6.*y*y*x*x;
258  case 8: return 4.*y*y*y*x;
259 
260  case 9: return 4.*y*r*r*r;
261  case 10: return 6.*y*y*r*r;
262  case 11: return 4.*y*y*y*r;
263 
264  // inner functions
265  case 12: return 12.*x*y*r*r;
266  case 13: return 12.*x*x*y*r;
267  case 14: return 12.*x*y*y*r;
268 
269  default:
270  libmesh_error_msg("Invalid shape function index shape = " << shape);
271  }
272  }
273  case FIFTH:
274  {
275  const Real x=p(0);
276  const Real y=p(1);
277  const Real r=1-x-y;
278  unsigned int shape=i;
279 
280  libmesh_assert_less (i, 21);
281 
282  if ((i>= 3&&i<= 6) && elem->point(0) > elem->point(1)) shape=9-i;
283  if ((i>= 7&&i<=10) && elem->point(1) > elem->point(2)) shape=17-i;
284  if ((i>=11&&i<=14) && elem->point(0) > elem->point(2)) shape=25-i;
285 
286  switch(shape)
287  {
288  //point functions
289  case 0: return pow<5>(r);
290  case 1: return pow<5>(x);
291  case 2: return pow<5>(y);
292 
293  //edge functions
294  case 3: return 5.*x *pow<4>(r);
295  case 4: return 10.*pow<2>(x)*pow<3>(r);
296  case 5: return 10.*pow<3>(x)*pow<2>(r);
297  case 6: return 5.*pow<4>(x)*r;
298 
299  case 7: return 5.*y *pow<4>(x);
300  case 8: return 10.*pow<2>(y)*pow<3>(x);
301  case 9: return 10.*pow<3>(y)*pow<2>(x);
302  case 10: return 5.*pow<4>(y)*x;
303 
304  case 11: return 5.*y *pow<4>(r);
305  case 12: return 10.*pow<2>(y)*pow<3>(r);
306  case 13: return 10.*pow<3>(y)*pow<2>(r);
307  case 14: return 5.*pow<4>(y)*r;
308 
309  //inner functions
310  case 15: return 20.*x*y*pow<3>(r);
311  case 16: return 30.*x*pow<2>(y)*pow<2>(r);
312  case 17: return 30.*pow<2>(x)*y*pow<2>(r);
313  case 18: return 20.*x*pow<3>(y)*r;
314  case 19: return 20.*pow<3>(x)*y*r;
315  case 20: return 30.*pow<2>(x)*pow<2>(y)*r;
316 
317  default:
318  libmesh_error_msg("Invalid shape function index shape = " << shape);
319  }
320  }
321  case SIXTH:
322  {
323  const Real x=p(0);
324  const Real y=p(1);
325  const Real r=1-x-y;
326  unsigned int shape=i;
327 
328  libmesh_assert_less (i, 28);
329 
330  if ((i>= 3&&i<= 7) && elem->point(0) > elem->point(1)) shape=10-i;
331  if ((i>= 8&&i<=12) && elem->point(1) > elem->point(2)) shape=20-i;
332  if ((i>=13&&i<=17) && elem->point(0) > elem->point(2)) shape=30-i;
333 
334  switch(shape)
335  {
336  //point functions
337  case 0: return pow<6>(r);
338  case 1: return pow<6>(x);
339  case 2: return pow<6>(y);
340 
341  //edge functions
342  case 3: return 6.*x *pow<5>(r);
343  case 4: return 15.*pow<2>(x)*pow<4>(r);
344  case 5: return 20.*pow<3>(x)*pow<3>(r);
345  case 6: return 15.*pow<4>(x)*pow<2>(r);
346  case 7: return 6.*pow<5>(x)*r;
347 
348  case 8: return 6.*y *pow<5>(x);
349  case 9: return 15.*pow<2>(y)*pow<4>(x);
350  case 10: return 20.*pow<3>(y)*pow<3>(x);
351  case 11: return 15.*pow<4>(y)*pow<2>(x);
352  case 12: return 6.*pow<5>(y)*x;
353 
354  case 13: return 6.*y *pow<5>(r);
355  case 14: return 15.*pow<2>(y)*pow<4>(r);
356  case 15: return 20.*pow<3>(y)*pow<3>(r);
357  case 16: return 15.*pow<4>(y)*pow<2>(r);
358  case 17: return 6.*pow<5>(y)*r;
359 
360  //inner functions
361  case 18: return 30.*x*y*pow<4>(r);
362  case 19: return 60.*x*pow<2>(y)*pow<3>(r);
363  case 20: return 60.* pow<2>(x)*y*pow<3>(r);
364  case 21: return 60.*x*pow<3>(y)*pow<2>(r);
365  case 22: return 60.*pow<3>(x)*y*pow<2>(r);
366  case 23: return 90.*pow<2>(x)*pow<2>(y)*pow<2>(r);
367  case 24: return 30.*x*pow<4>(y)*r;
368  case 25: return 60.*pow<2>(x)*pow<3>(y)*r;
369  case 26: return 60.*pow<3>(x)*pow<2>(y)*r;
370  case 27: return 30.*pow<4>(x)*y*r;
371 
372  default:
373  libmesh_error_msg("Invalid shape function index shape = " << shape);
374  } // switch shape
375  } // case TRI6
376  default:
377  libmesh_error_msg("Invalid totalorder = " << totalorder);
378  } // switch order
379 
380  default:
381  libmesh_error_msg("ERROR: Unsupported element type = " << type);
382  } // switch type
383 }
384 
385 
386 
387 template <>
389  const Order,
390  const unsigned int,
391  const unsigned int,
392  const Point &)
393 {
394  libmesh_error_msg("Bernstein polynomials require the element type \nbecause edge orientation is needed.");
395  return 0.;
396 }
397 
398 
399 
400 template <>
402  const Order order,
403  const unsigned int i,
404  const unsigned int j,
405  const Point & p,
406  const bool add_p_level)
407 {
408  libmesh_assert(elem);
409 
410  const ElemType type = elem->type();
411 
412  const Order totalorder =
413  static_cast<Order>(order + add_p_level * elem->p_level());
414 
415  switch (type)
416  {
417  // Hierarchic shape functions on the quadrilateral.
418  case QUAD4:
419  case QUAD9:
420  {
421  // Compute quad shape functions as a tensor-product
422  const Real xi = p(0);
423  const Real eta = p(1);
424 
425  libmesh_assert_less (i, (totalorder+1u)*(totalorder+1u));
426 
427  unsigned int i0, i1;
428 
429  // Vertex DoFs
430  if (i == 0)
431  { i0 = 0; i1 = 0; }
432  else if (i == 1)
433  { i0 = 1; i1 = 0; }
434  else if (i == 2)
435  { i0 = 1; i1 = 1; }
436  else if (i == 3)
437  { i0 = 0; i1 = 1; }
438 
439 
440  // Edge DoFs
441  else if (i < totalorder + 3u)
442  { i0 = i - 2; i1 = 0; }
443  else if (i < 2u*totalorder + 2)
444  { i0 = 1; i1 = i - totalorder - 1; }
445  else if (i < 3u*totalorder + 1)
446  { i0 = i - 2u*totalorder; i1 = 1; }
447  else if (i < 4u*totalorder)
448  { i0 = 0; i1 = i - 3u*totalorder + 1; }
449  // Interior DoFs
450  else
451  {
452  unsigned int basisnum = i - 4*totalorder;
453  i0 = square_number_column[basisnum] + 2;
454  i1 = square_number_row[basisnum] + 2;
455  }
456 
457 
458  // Flip odd degree of freedom values if necessary
459  // to keep continuity on sides
460  if ((i>= 4 && i<= 4+ totalorder-2u) && elem->point(0) > elem->point(1)) i0=totalorder+2-i0;
461  else if ((i>= 4+ totalorder-1u && i<= 4+2*totalorder-3u) && elem->point(1) > elem->point(2)) i1=totalorder+2-i1;
462  else if ((i>= 4+2*totalorder-2u && i<= 4+3*totalorder-4u) && elem->point(3) > elem->point(2)) i0=totalorder+2-i0;
463  else if ((i>= 4+3*totalorder-3u && i<= 4+4*totalorder-5u) && elem->point(0) > elem->point(3)) i1=totalorder+2-i1;
464 
465  switch (j)
466  {
467  // d()/dxi
468  case 0:
469  return (FE<1,BERNSTEIN>::shape_deriv(EDGE3, totalorder, i0, 0, xi)*
470  FE<1,BERNSTEIN>::shape (EDGE3, totalorder, i1, eta));
471 
472  // d()/deta
473  case 1:
474  return (FE<1,BERNSTEIN>::shape (EDGE3, totalorder, i0, xi)*
475  FE<1,BERNSTEIN>::shape_deriv(EDGE3, totalorder, i1, 0, eta));
476 
477  default:
478  libmesh_error_msg("Invalid shape function derivative j = " << j);
479  }
480  }
481 
482  // Bernstein shape functions on the 8-noded quadrilateral
483  // is handled separately.
484  case QUAD8:
485  case QUADSHELL8:
486  {
487  libmesh_assert_less (totalorder, 3);
488 
489  const Real xi = p(0);
490  const Real eta = p(1);
491 
492  libmesh_assert_less (i, 8);
493 
494  // 0 1 2 3 4 5 6 7 8
495  static const unsigned int i0[] = {0, 1, 1, 0, 2, 1, 2, 0, 2};
496  static const unsigned int i1[] = {0, 0, 1, 1, 0, 2, 1, 2, 2};
497  static const Real scal[] = {-0.25, -0.25, -0.25, -0.25, 0.5, 0.5, 0.5, 0.5};
498  switch (j)
499  {
500  // d()/dxi
501  case 0:
502  return (FE<1,BERNSTEIN>::shape_deriv(EDGE3, totalorder, i0[i], 0, xi)*
503  FE<1,BERNSTEIN>::shape (EDGE3, totalorder, i1[i], eta)
504  +scal[i]*
505  FE<1,BERNSTEIN>::shape_deriv(EDGE3, totalorder, i0[8], 0, xi)*
506  FE<1,BERNSTEIN>::shape (EDGE3, totalorder, i1[8], eta));
507 
508  // d()/deta
509  case 1:
510  return (FE<1,BERNSTEIN>::shape (EDGE3, totalorder, i0[i], xi)*
511  FE<1,BERNSTEIN>::shape_deriv(EDGE3, totalorder, i1[i], 0, eta)
512  +scal[i]*
513  FE<1,BERNSTEIN>::shape (EDGE3, totalorder, i0[8], xi)*
514  FE<1,BERNSTEIN>::shape_deriv(EDGE3, totalorder, i1[8], 0, eta));
515 
516  default:
517  libmesh_error_msg("Invalid shape function derivative j = " << j);
518  }
519  }
520 
521  case TRI3:
522  case TRISHELL3:
523  libmesh_assert_less (totalorder, 2);
524  libmesh_fallthrough();
525  case TRI6:
526  {
527  // I have been lazy here and am using finite differences
528  // to compute the derivatives!
529  const Real eps = 1.e-6;
530 
531  switch (j)
532  {
533  // d()/dxi
534  case 0:
535  {
536  const Point pp(p(0)+eps, p(1));
537  const Point pm(p(0)-eps, p(1));
538 
539  return (FE<2,BERNSTEIN>::shape(elem, totalorder, i, pp) -
540  FE<2,BERNSTEIN>::shape(elem, totalorder, i, pm))/2./eps;
541  }
542 
543  // d()/deta
544  case 1:
545  {
546  const Point pp(p(0), p(1)+eps);
547  const Point pm(p(0), p(1)-eps);
548 
549  return (FE<2,BERNSTEIN>::shape(elem, totalorder, i, pp) -
550  FE<2,BERNSTEIN>::shape(elem, totalorder, i, pm))/2./eps;
551  }
552 
553 
554  default:
555  libmesh_error_msg("Invalid shape function derivative j = " << j);
556  }
557  }
558 
559  default:
560  libmesh_error_msg("ERROR: Unsupported element type = " << type);
561  }
562 }
563 
564 
565 
566 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
567 template <>
569  const Order,
570  const unsigned int,
571  const unsigned int,
572  const Point &)
573 {
574  static bool warning_given = false;
575 
576  if (!warning_given)
577  libMesh::err << "Second derivatives for Bernstein elements "
578  << "are not yet implemented!"
579  << std::endl;
580 
581  warning_given = true;
582  return 0.;
583 }
584 
585 
586 
587 template <>
589  const Order,
590  const unsigned int,
591  const unsigned int,
592  const Point &,
593  const bool)
594 {
595  static bool warning_given = false;
596 
597  if (!warning_given)
598  libMesh::err << "Second derivatives for Bernstein elements "
599  << "are not yet implemented!"
600  << std::endl;
601 
602  warning_given = true;
603  return 0.;
604 }
605 
606 #endif
607 
608 } // namespace libMesh
609 
610 
611 #endif// LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
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::SIXTH
Definition: enum_order.h:47
libMesh
The libMesh namespace provides an interface to certain functionality in the library.
Definition: factoryfunction.C:55
libMesh::Order
Order
Definition: enum_order.h:40
libMesh::Elem::p_level
unsigned int p_level() const
Definition: elem.h:2512
libMesh::FIFTH
Definition: enum_order.h:46
libMesh::SECOND
Definition: enum_order.h:43
libMesh::Elem::point
const Point & point(const unsigned int i) const
Definition: elem.h:1955
libMesh::square_number_column
const unsigned char square_number_column[]
Definition: number_lookups.C:56
libMesh::libmesh_assert
libmesh_assert(ctx)
libMesh::FE::shape
static OutputShape shape(const ElemType t, const Order o, const unsigned int i, const Point &p)
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::FE
A specific instantiation of the FEBase class.
Definition: fe.h:95
libMesh::QUAD4
Definition: enum_elem_type.h:41
libMesh::Point
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
libMesh::TRI3
Definition: enum_elem_type.h:39
libMesh::FOURTH
Definition: enum_order.h:45
libMesh::TRI6
Definition: enum_elem_type.h:40
libMesh::Utility::pow
T pow(const T &x)
Definition: utility.h:246
libMesh::QUADSHELL8
Definition: enum_elem_type.h:73
libMesh::EDGE3
Definition: enum_elem_type.h:36
libMesh::square_number_row
const unsigned char square_number_row[]
Definition: number_lookups.C:69
libMesh::QUADSHELL4
Definition: enum_elem_type.h:72
libMesh::Elem
This is the base class from which all geometric element types are derived.
Definition: elem.h:100
libMesh::THIRD
Definition: enum_order.h:44
libMesh::QUAD9
Definition: enum_elem_type.h:43
libMesh::err
OStreamProxy err
libMesh::Real
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Definition: libmesh_common.h:121
libMesh::TRISHELL3
Definition: enum_elem_type.h:71
libMesh::Elem::type
virtual ElemType type() const =0
libMesh::FIRST
Definition: enum_order.h:42
libMesh::QUAD8
Definition: enum_elem_type.h:42
libMesh::ElemType
ElemType
Defines an enum for geometric element types.
Definition: enum_elem_type.h:33