libMesh
cell_prism18.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/side.h"
21 #include "libmesh/cell_prism18.h"
22 #include "libmesh/edge_edge3.h"
23 #include "libmesh/face_quad9.h"
24 #include "libmesh/face_tri6.h"
25 #include "libmesh/enum_io_package.h"
26 #include "libmesh/enum_order.h"
27 #include "libmesh/int_range.h"
28 
29 namespace libMesh
30 {
31 
32 
33 
34 // ------------------------------------------------------------
35 // Prism18 class static member initializations
36 const int Prism18::num_nodes;
37 const int Prism18::num_sides;
38 const int Prism18::num_edges;
39 const int Prism18::num_children;
40 const int Prism18::nodes_per_side;
41 const int Prism18::nodes_per_edge;
42 
43 const unsigned int Prism18::side_nodes_map[Prism18::num_sides][Prism18::nodes_per_side] =
44  {
45  {0, 2, 1, 8, 7, 6, 99, 99, 99}, // Side 0
46  {0, 1, 4, 3, 6, 10, 12, 9, 15}, // Side 1
47  {1, 2, 5, 4, 7, 11, 13, 10, 16}, // Side 2
48  {2, 0, 3, 5, 8, 9, 14, 11, 17}, // Side 3
49  {3, 4, 5, 12, 13, 14, 99, 99, 99} // Side 4
50  };
51 
52 const unsigned int Prism18::edge_nodes_map[Prism18::num_edges][Prism18::nodes_per_edge] =
53  {
54  {0, 1, 6}, // Edge 0
55  {1, 2, 7}, // Edge 1
56  {0, 2, 8}, // Edge 2
57  {0, 3, 9}, // Edge 3
58  {1, 4, 10}, // Edge 4
59  {2, 5, 11}, // Edge 5
60  {3, 4, 12}, // Edge 6
61  {4, 5, 13}, // Edge 7
62  {3, 5, 14} // Edge 8
63  };
64 
65 
66 // ------------------------------------------------------------
67 // Prism18 class member functions
68 
69 bool Prism18::is_vertex(const unsigned int i) const
70 {
71  if (i < 6)
72  return true;
73  return false;
74 }
75 
76 bool Prism18::is_edge(const unsigned int i) const
77 {
78  if (i < 6)
79  return false;
80  if (i > 14)
81  return false;
82  return true;
83 }
84 
85 bool Prism18::is_face(const unsigned int i) const
86 {
87  if (i > 14)
88  return true;
89  return false;
90 }
91 
92 bool Prism18::is_node_on_side(const unsigned int n,
93  const unsigned int s) const
94 {
95  libmesh_assert_less (s, n_sides());
96  return std::find(std::begin(side_nodes_map[s]),
97  std::end(side_nodes_map[s]),
98  n) != std::end(side_nodes_map[s]);
99 }
100 
101 std::vector<unsigned>
102 Prism18::nodes_on_side(const unsigned int s) const
103 {
104  libmesh_assert_less(s, n_sides());
105  auto trim = (s > 0 && s < 4) ? 0 : 3;
106  return {std::begin(side_nodes_map[s]), std::end(side_nodes_map[s]) - trim};
107 }
108 
109 bool Prism18::is_node_on_edge(const unsigned int n,
110  const unsigned int e) const
111 {
112  libmesh_assert_less (e, n_edges());
113  return std::find(std::begin(edge_nodes_map[e]),
114  std::end(edge_nodes_map[e]),
115  n) != std::end(edge_nodes_map[e]);
116 }
117 
118 
119 
120 bool Prism18::has_affine_map() const
121 {
122  // Make sure z edges are affine
123  Point v = this->point(3) - this->point(0);
124  if (!v.relative_fuzzy_equals(this->point(4) - this->point(1)) ||
125  !v.relative_fuzzy_equals(this->point(5) - this->point(2)))
126  return false;
127  // Make sure edges are straight
128  v /= 2;
129  if (!v.relative_fuzzy_equals(this->point(9) - this->point(0)) ||
130  !v.relative_fuzzy_equals(this->point(10) - this->point(1)) ||
131  !v.relative_fuzzy_equals(this->point(11) - this->point(2)) ||
132  !v.relative_fuzzy_equals(this->point(15) - this->point(6)) ||
133  !v.relative_fuzzy_equals(this->point(16) - this->point(7)) ||
134  !v.relative_fuzzy_equals(this->point(17) - this->point(8)))
135  return false;
136  v = (this->point(1) - this->point(0))/2;
137  if (!v.relative_fuzzy_equals(this->point(6) - this->point(0)) ||
138  !v.relative_fuzzy_equals(this->point(12) - this->point(3)))
139  return false;
140  v = (this->point(2) - this->point(0))/2;
141  if (!v.relative_fuzzy_equals(this->point(8) - this->point(0)) ||
142  !v.relative_fuzzy_equals(this->point(14) - this->point(3)))
143  return false;
144  v = (this->point(2) - this->point(1))/2;
145  if (!v.relative_fuzzy_equals(this->point(7) - this->point(1)) ||
146  !v.relative_fuzzy_equals(this->point(13) - this->point(4)))
147  return false;
148  return true;
149 }
150 
151 
152 
153 Order Prism18::default_order() const
154 {
155  return SECOND;
156 }
157 
158 dof_id_type Prism18::key (const unsigned int s) const
159 {
160  libmesh_assert_less (s, this->n_sides());
161 
162  switch (s)
163  {
164  case 0: // the triangular face at z=0
165  {
166  return Prism::key(0);
167  }
168  case 1: // the quad face at y=0
169  {
170  return Elem::compute_key (this->node_id(15));
171  }
172  case 2: // the other quad face
173  {
174  return Elem::compute_key (this->node_id(16));
175  }
176  case 3: // the quad face at x=0
177  {
178  return Elem::compute_key (this->node_id(17));
179  }
180  case 4: // the triangular face at z=1
181  {
182  return Prism::key(4);
183  }
184  default:
185  libmesh_error_msg("Invalid side " << s);
186  }
187 }
188 
189 
190 
191 unsigned int Prism18::which_node_am_i(unsigned int side,
192  unsigned int side_node) const
193 {
194  libmesh_assert_less (side, this->n_sides());
195 
196  // Never more than 9 nodes per side.
197  libmesh_assert_less(side_node, 9);
198 
199  // Some sides have 6 nodes.
200  libmesh_assert(!(side==0 || side==4) || side_node < 6);
201 
202  return Prism18::side_nodes_map[side][side_node];
203 }
204 
205 
206 
207 std::unique_ptr<Elem> Prism18::build_side_ptr (const unsigned int i,
208  bool proxy)
209 {
210  libmesh_assert_less (i, this->n_sides());
211 
212  if (proxy)
213  {
214  switch(i)
215  {
216  case 0:
217  case 4:
218  return libmesh_make_unique<Side<Tri6,Prism18>>(this,i);
219 
220  case 1:
221  case 2:
222  case 3:
223  return libmesh_make_unique<Side<Quad9,Prism18>>(this,i);
224 
225  default:
226  libmesh_error_msg("Invalid side i = " << i);
227  }
228  }
229 
230  else
231  {
232  // Return value
233  std::unique_ptr<Elem> face;
234 
235  switch (i)
236  {
237  case 0: // the triangular face at z=-1
238  case 4: // the triangular face at z=1
239  {
240  face = libmesh_make_unique<Tri6>();
241  break;
242  }
243  case 1: // the quad face at y=0
244  case 2: // the other quad face
245  case 3: // the quad face at x=0
246  {
247  face = libmesh_make_unique<Quad9>();
248  break;
249  }
250  default:
251  libmesh_error_msg("Invalid side i = " << i);
252  }
253 
254  face->subdomain_id() = this->subdomain_id();
255 
256  // Set the nodes
257  for (auto n : face->node_index_range())
258  face->set_node(n) = this->node_ptr(Prism18::side_nodes_map[i][n]);
259 
260  return face;
261  }
262 }
263 
264 
265 
266 void Prism18::build_side_ptr (std::unique_ptr<Elem> & side,
267  const unsigned int i)
268 {
269  libmesh_assert_less (i, this->n_sides());
270 
271  switch (i)
272  {
273  case 0: // the triangular face at z=-1
274  case 4: // the triangular face at z=1
275  {
276  if (!side.get() || side->type() != TRI6)
277  {
278  side = this->build_side_ptr(i, false);
279  return;
280  }
281  break;
282  }
283 
284  case 1: // the quad face at y=0
285  case 2: // the other quad face
286  case 3: // the quad face at x=0
287  {
288  if (!side.get() || side->type() != QUAD9)
289  {
290  side = this->build_side_ptr(i, false);
291  return;
292  }
293  break;
294  }
295 
296  default:
297  libmesh_error_msg("Invalid side i = " << i);
298  }
299 
300  side->subdomain_id() = this->subdomain_id();
301 
302  // Set the nodes
303  for (auto n : side->node_index_range())
304  side->set_node(n) = this->node_ptr(Prism18::side_nodes_map[i][n]);
305 }
306 
307 
308 
309 std::unique_ptr<Elem> Prism18::build_edge_ptr (const unsigned int i)
310 {
311  libmesh_assert_less (i, this->n_edges());
312 
313  return libmesh_make_unique<SideEdge<Edge3,Prism18>>(this,i);
314 }
315 
316 
317 
318 void Prism18::connectivity(const unsigned int sc,
319  const IOPackage iop,
320  std::vector<dof_id_type> & conn) const
321 {
322  libmesh_assert(_nodes);
323  libmesh_assert_less (sc, this->n_sub_elem());
324  libmesh_assert_not_equal_to (iop, INVALID_IO_PACKAGE);
325 
326  switch (iop)
327  {
328  case TECPLOT:
329  {
330  conn.resize(8);
331  switch (sc)
332  {
333 
334  case 0:
335  {
336  conn[0] = this->node_id(0)+1;
337  conn[1] = this->node_id(6)+1;
338  conn[2] = this->node_id(8)+1;
339  conn[3] = this->node_id(8)+1;
340  conn[4] = this->node_id(9)+1;
341  conn[5] = this->node_id(15)+1;
342  conn[6] = this->node_id(17)+1;
343  conn[7] = this->node_id(17)+1;
344 
345  return;
346  }
347 
348  case 1:
349  {
350  conn[0] = this->node_id(6)+1;
351  conn[1] = this->node_id(1)+1;
352  conn[2] = this->node_id(7)+1;
353  conn[3] = this->node_id(7)+1;
354  conn[4] = this->node_id(15)+1;
355  conn[5] = this->node_id(10)+1;
356  conn[6] = this->node_id(16)+1;
357  conn[7] = this->node_id(16)+1;
358 
359  return;
360  }
361 
362  case 2:
363  {
364  conn[0] = this->node_id(8)+1;
365  conn[1] = this->node_id(7)+1;
366  conn[2] = this->node_id(2)+1;
367  conn[3] = this->node_id(2)+1;
368  conn[4] = this->node_id(17)+1;
369  conn[5] = this->node_id(16)+1;
370  conn[6] = this->node_id(11)+1;
371  conn[7] = this->node_id(11)+1;
372 
373  return;
374  }
375 
376  case 3:
377  {
378  conn[0] = this->node_id(6)+1;
379  conn[1] = this->node_id(7)+1;
380  conn[2] = this->node_id(8)+1;
381  conn[3] = this->node_id(8)+1;
382  conn[4] = this->node_id(15)+1;
383  conn[5] = this->node_id(16)+1;
384  conn[6] = this->node_id(17)+1;
385  conn[7] = this->node_id(17)+1;
386 
387  return;
388  }
389 
390  case 4:
391  {
392  conn[0] = this->node_id(9)+1;
393  conn[1] = this->node_id(15)+1;
394  conn[2] = this->node_id(17)+1;
395  conn[3] = this->node_id(17)+1;
396  conn[4] = this->node_id(3)+1;
397  conn[5] = this->node_id(12)+1;
398  conn[6] = this->node_id(14)+1;
399  conn[7] = this->node_id(14)+1;
400 
401  return;
402  }
403 
404  case 5:
405  {
406  conn[0] = this->node_id(15)+1;
407  conn[1] = this->node_id(10)+1;
408  conn[2] = this->node_id(16)+1;
409  conn[3] = this->node_id(16)+1;
410  conn[4] = this->node_id(12)+1;
411  conn[5] = this->node_id(4)+1;
412  conn[6] = this->node_id(13)+1;
413  conn[7] = this->node_id(13)+1;
414 
415  return;
416  }
417 
418  case 6:
419  {
420  conn[0] = this->node_id(17)+1;
421  conn[1] = this->node_id(16)+1;
422  conn[2] = this->node_id(11)+1;
423  conn[3] = this->node_id(11)+1;
424  conn[4] = this->node_id(14)+1;
425  conn[5] = this->node_id(13)+1;
426  conn[6] = this->node_id(5)+1;
427  conn[7] = this->node_id(5)+1;
428 
429  return;
430  }
431 
432  case 7:
433  {
434  conn[0] = this->node_id(15)+1;
435  conn[1] = this->node_id(16)+1;
436  conn[2] = this->node_id(17)+1;
437  conn[3] = this->node_id(17)+1;
438  conn[4] = this->node_id(12)+1;
439  conn[5] = this->node_id(13)+1;
440  conn[6] = this->node_id(14)+1;
441  conn[7] = this->node_id(14)+1;
442 
443  return;
444  }
445 
446  default:
447  libmesh_error_msg("Invalid sc = " << sc);
448  }
449 
450  }
451 
452  case VTK:
453  {
454  // VTK now supports VTK_BIQUADRATIC_QUADRATIC_WEDGE directly
455  const unsigned int conn_size = 18;
456  conn.resize(conn_size);
457 
458  // VTK's VTK_BIQUADRATIC_QUADRATIC_WEDGE first 9 (vertex) and
459  // last 3 (mid-face) nodes match. The middle and top layers
460  // of mid-edge nodes are reversed from LibMesh's.
461  for (auto i : index_range(conn))
462  conn[i] = this->node_id(i);
463 
464  // top "ring" of mid-edge nodes
465  conn[9] = this->node_id(12);
466  conn[10] = this->node_id(13);
467  conn[11] = this->node_id(14);
468 
469  // middle "ring" of mid-edge nodes
470  conn[12] = this->node_id(9);
471  conn[13] = this->node_id(10);
472  conn[14] = this->node_id(11);
473 
474  return;
475  }
476 
477  default:
478  libmesh_error_msg("Unsupported IO package " << iop);
479  }
480 }
481 
482 
483 
484 
485 unsigned int Prism18::n_second_order_adjacent_vertices (const unsigned int n) const
486 {
487  switch (n)
488  {
489  case 6:
490  case 7:
491  case 8:
492  case 9:
493  case 10:
494  case 11:
495  case 12:
496  case 13:
497  case 14:
498  return 2;
499 
500  case 15:
501  case 16:
502  case 17:
503  return 4;
504 
505  default:
506  libmesh_error_msg("Invalid node n = " << n);
507  }
508 }
509 
510 
511 
512 
513 
514 unsigned short int Prism18::second_order_adjacent_vertex (const unsigned int n,
515  const unsigned int v) const
516 {
517  libmesh_assert_greater_equal (n, this->n_vertices());
518  libmesh_assert_less (n, this->n_nodes());
519 
520  switch (n)
521  {
522  /*
523  * These nodes are unique to \p Prism18,
524  * let our _remaining_... matrix handle
525  * this.
526  */
527  case 15:
528  case 16:
529  case 17:
530  {
531  libmesh_assert_less (v, 4);
532  return _remaining_second_order_adjacent_vertices[n-15][v];
533  }
534 
535  /*
536  * All other second-order nodes (6,...,14) are
537  * identical with Prism15 and are therefore
538  * delegated to the _second_order matrix of
539  * \p Prism
540  */
541  default:
542  {
543  libmesh_assert_less (v, 2);
544  return _second_order_adjacent_vertices[n-this->n_vertices()][v];
545  }
546 
547  }
548 
549  libmesh_error_msg("We'll never ge here!");
550  return static_cast<unsigned short int>(-1);
551 }
552 
553 
554 
555 const unsigned short int Prism18::_remaining_second_order_adjacent_vertices[3][4] =
556  {
557  { 0, 1, 3, 4}, // vertices adjacent to node 15
558  { 1, 2, 4, 5}, // vertices adjacent to node 16
559  { 0, 2, 3, 5} // vertices adjacent to node 17
560  };
561 
562 
563 
564 std::pair<unsigned short int, unsigned short int>
565 Prism18::second_order_child_vertex (const unsigned int n) const
566 {
567  libmesh_assert_greater_equal (n, this->n_vertices());
568  libmesh_assert_less (n, this->n_nodes());
569 
570  return std::pair<unsigned short int, unsigned short int>
571  (_second_order_vertex_child_number[n],
572  _second_order_vertex_child_index[n]);
573 }
574 
575 
576 
577 Real Prism18::volume () const
578 {
579  // This specialization is good for Lagrange mappings only
580  if (this->mapping_type() != LAGRANGE_MAP)
581  return this->Elem::volume();
582 
583  // Make copies of our points. It makes the subsequent calculations a bit
584  // shorter and avoids dereferencing the same pointer multiple times.
585  Point
586  x0 = point(0), x1 = point(1), x2 = point(2), x3 = point(3), x4 = point(4), x5 = point(5),
587  x6 = point(6), x7 = point(7), x8 = point(8), x9 = point(9), x10 = point(10), x11 = point(11),
588  x12 = point(12), x13 = point(13), x14 = point(14), x15 = point(15), x16 = point(16), x17 = point(17);
589 
590  // The number of components in the dx_dxi, dx_deta, and dx_dzeta arrays.
591  const int n_components = 16;
592 
593  // Terms are copied directly from a Python script.
594  Point dx_dxi[n_components] =
595  {
596  -x10 + 4*x15 - 3*x9,
597  3*x0/2 + x1/2 + 2*x12 - 3*x3/2 - x4/2 - 2*x6,
598  -3*x0/2 - x1/2 + x10 + 2*x12 - 4*x15 - 3*x3/2 - x4/2 + 2*x6 + 3*x9,
599  -4*x15 + 4*x16 - 4*x17 + 4*x9,
600  -2*x0 - 2*x12 + 2*x13 - 2*x14 + 2*x3 + 2*x6 - 2*x7 + 2*x8,
601  2*x0 - 2*x12 + 2*x13 - 2*x14 + 4*x15 - 4*x16 + 4*x17 + 2*x3 - 2*x6 + 2*x7 - 2*x8 - 4*x9,
602  Point(0,0,0),
603  Point(0,0,0),
604  4*x10 - 8*x15 + 4*x9,
605  -2*x0 - 2*x1 - 4*x12 + 2*x3 + 2*x4 + 4*x6,
606  2*x0 + 2*x1 - 4*x10 - 4*x12 + 8*x15 + 2*x3 + 2*x4 - 4*x6 - 4*x9,
607  Point(0,0,0),
608  Point(0,0,0),
609  Point(0,0,0),
610  Point(0,0,0),
611  Point(0,0,0)
612  };
613 
614  Point dx_deta[n_components] =
615  {
616  -x11 + 4*x17 - 3*x9,
617  3*x0/2 + 2*x14 + x2/2 - 3*x3/2 - x5/2 - 2*x8,
618  -3*x0/2 + x11 + 2*x14 - 4*x17 - x2/2 - 3*x3/2 - x5/2 + 2*x8 + 3*x9,
619  4*x11 - 8*x17 + 4*x9,
620  -2*x0 - 4*x14 - 2*x2 + 2*x3 + 2*x5 + 4*x8,
621  2*x0 - 4*x11 - 4*x14 + 8*x17 + 2*x2 + 2*x3 + 2*x5 - 4*x8 - 4*x9,
622  Point(0,0,0),
623  Point(0,0,0),
624  -4*x15 + 4*x16 - 4*x17 + 4*x9,
625  -2*x0 - 2*x12 + 2*x13 - 2*x14 + 2*x3 + 2*x6 - 2*x7 + 2*x8,
626  2*x0 - 2*x12 + 2*x13 - 2*x14 + 4*x15 - 4*x16 + 4*x17 + 2*x3 - 2*x6 + 2*x7 - 2*x8 - 4*x9,
627  Point(0,0,0),
628  Point(0,0,0),
629  Point(0,0,0),
630  Point(0,0,0),
631  Point(0,0,0)
632  };
633 
634  Point dx_dzeta[n_components] =
635  {
636  -x0/2 + x3/2,
637  x0 + x3 - 2*x9,
638  Point(0,0,0),
639  3*x0/2 + 2*x14 + x2/2 - 3*x3/2 - x5/2 - 2*x8,
640  -3*x0 + 2*x11 + 4*x14 - 8*x17 - x2 - 3*x3 - x5 + 4*x8 + 6*x9,
641  Point(0,0,0),
642  -x0 - 2*x14 - x2 + x3 + x5 + 2*x8,
643  2*x0 - 4*x11 - 4*x14 + 8*x17 + 2*x2 + 2*x3 + 2*x5 - 4*x8 - 4*x9,
644  3*x0/2 + x1/2 + 2*x12 - 3*x3/2 - x4/2 - 2*x6,
645  -3*x0 - x1 + 2*x10 + 4*x12 - 8*x15 - 3*x3 - x4 + 4*x6 + 6*x9,
646  Point(0,0,0),
647  -2*x0 - 2*x12 + 2*x13 - 2*x14 + 2*x3 + 2*x6 - 2*x7 + 2*x8,
648  4*x0 - 4*x12 + 4*x13 - 4*x14 + 8*x15 - 8*x16 + 8*x17 + 4*x3 - 4*x6 + 4*x7 - 4*x8 - 8*x9,
649  Point(0,0,0),
650  -x0 - x1 - 2*x12 + x3 + x4 + 2*x6,
651  2*x0 + 2*x1 - 4*x10 - 4*x12 + 8*x15 + 2*x3 + 2*x4 - 4*x6 - 4*x9
652  };
653 
654  // The quadrature rule for the Prism18 is a tensor product between a
655  // FOURTH-order TRI rule (in xi, eta) and a FIFTH-order EDGE rule
656  // in zeta.
657 
658  // Number of points in the 2D quadrature rule.
659  const int N2D = 6;
660 
661  // Parameters of the 2D rule
662  static const Real
663  w1 = Real(1.1169079483900573284750350421656140e-01L),
664  w2 = Real(5.4975871827660933819163162450105264e-02L),
665  a1 = Real(4.4594849091596488631832925388305199e-01L),
666  a2 = Real(9.1576213509770743459571463402201508e-02L);
667 
668  // Points and weights of the 2D rule
669  static const Real w2D[N2D] = {w1, w1, w1, w2, w2, w2};
670 
671  // Quadrature point locations raised to powers. xi[0][2] is
672  // quadrature point 0, squared, xi[1][1] is quadrature point 1 to the
673  // first power, etc. This lets us avoid calling std::pow inside the
674  // loops below.
675  static const Real xi[N2D][3] =
676  {
677  // ^0 ^1 ^2
678  { 1., a1, a1*a1},
679  { 1., 1-2*a1, (1-2*a1)*(1-2*a1)},
680  { 1., a1, a1*a1},
681  { 1., a2, a2*a2},
682  { 1., 1-2*a2, (1-2*a2)*(1-2*a2)},
683  { 1., a2, a2*a2}
684  };
685 
686  static const Real eta[N2D][3] =
687  {
688  // ^0 ^1 ^2
689  { 1., a1, a1*a1},
690  { 1., a1, a1*a1},
691  { 1., 1-2*a1, (1-2*a1)*(1-2*a1)},
692  { 1., a2, a2*a2},
693  { 1., a2, a2*a2},
694  { 1., 1-2*a2, (1-2*a2)*(1-2*a2)}
695  };
696 
697  // Number of points in the 1D quadrature rule.
698  const int N1D = 3;
699 
700  // Points and weights of the 1D quadrature rule.
701  static const Real w1D[N1D] = {5./9, 8./9, 5./9};
702 
703  const Real zeta[N1D][3] =
704  {
705  //^0 ^1 ^2
706  { 1., -std::sqrt(15)/5., 15./25},
707  { 1., 0., 0.},
708  { 1., std::sqrt(15)/5., 15./25}
709  };
710 
711  // The integer exponents for each term.
712  static const int exponents[n_components][3] =
713  {
714  {0, 0, 0},
715  {0, 0, 1},
716  {0, 0, 2},
717  {0, 1, 0},
718  {0, 1, 1},
719  {0, 1, 2},
720  {0, 2, 0},
721  {0, 2, 1},
722  {1, 0, 0},
723  {1, 0, 1},
724  {1, 0, 2},
725  {1, 1, 0},
726  {1, 1, 1},
727  {1, 2, 0},
728  {2, 0, 0},
729  {2, 0, 1}
730  };
731 
732  Real vol = 0.;
733  for (int i=0; i<N2D; ++i)
734  for (int j=0; j<N1D; ++j)
735  {
736  // Compute dx_dxi, dx_deta, dx_dzeta at the current quadrature point.
737  Point dx_dxi_q, dx_deta_q, dx_dzeta_q;
738  for (int c=0; c<n_components; ++c)
739  {
740  Real coeff =
741  xi[i][exponents[c][0]]*
742  eta[i][exponents[c][1]]*
743  zeta[j][exponents[c][2]];
744 
745  dx_dxi_q += coeff * dx_dxi[c];
746  dx_deta_q += coeff * dx_deta[c];
747  dx_dzeta_q += coeff * dx_dzeta[c];
748  }
749 
750  // Compute scalar triple product, multiply by weight, and accumulate volume.
751  vol += w2D[i] * w1D[j] * triple_product(dx_dxi_q, dx_deta_q, dx_dzeta_q);
752  }
753 
754  return vol;
755 }
756 
757 
758 
759 
760 #ifdef LIBMESH_ENABLE_AMR
761 
762 const float Prism18::_embedding_matrix[Prism18::num_children][Prism18::num_nodes][Prism18::num_nodes] =
763  {
764  // embedding matrix for child 0
765  {
766  // 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17
767  { 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 0
768  { 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 1
769  { 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 2
770  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0.}, // 3
771  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0.}, // 4
772  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1.}, // 5
773  { 0.375, -0.125, 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 6
774  { 0., -0.125, -0.125, 0., 0., 0., 0.5, 0.25, 0.5, 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 7
775  { 0.375, 0., -0.125, 0., 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 8
776  { 0.375, 0., 0., -0.125, 0., 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0., 0., 0., 0.}, // 9
777  { 0., 0., 0., 0., 0., 0., 0.375, 0., 0., 0., 0., 0., -0.125, 0., 0., 0.75, 0., 0.}, // 10
778  { 0., 0., 0., 0., 0., 0., 0., 0., 0.375, 0., 0., 0., 0., 0., -0.125, 0., 0., 0.75}, // 11
779  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.375, -0.125, 0., 0., 0., 0., 0.75, 0., 0.}, // 12
780  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, -0.125, 0., 0., 0., 0.5, 0.25, 0.5}, // 13
781  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.375, 0., -0.125, 0., 0., 0., 0., 0., 0.75}, // 14
782  { 0.140625,-0.046875, 0.,-0.046875, 0.015625, 0., 0.28125, 0., 0., 0.28125, -0.09375, 0., -0.09375, 0., 0., 0.5625, 0., 0.}, // 15
783  { 0.,-0.046875,-0.046875, 0., 0.015625, 0.015625, 0.1875, 0.09375, 0.1875, 0., -0.09375, -0.09375, -0.0625, -0.03125, -0.0625, 0.375, 0.1875, 0.375}, // 16
784  { 0.140625, 0.,-0.046875,-0.046875, 0., 0.015625, 0., 0., 0.28125, 0.28125, 0., -0.09375, 0., 0., -0.09375, 0., 0., 0.5625} // 17
785  },
786 
787  // embedding matrix for child 1
788  {
789  // 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17
790  { 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 0
791  { 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 1
792  { 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 2
793  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0.}, // 3
794  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0.}, // 4
795  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0.}, // 5
796  { -0.125, 0.375, 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 6
797  { 0., 0.375, -0.125, 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 7
798  { -0.125, 0., -0.125, 0., 0., 0., 0.5, 0.5, 0.25, 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 8
799  { 0., 0., 0., 0., 0., 0., 0.375, 0., 0., 0., 0., 0., -0.125, 0., 0., 0.75, 0., 0.}, // 9
800  { 0., 0.375, 0., 0., -0.125, 0., 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0., 0., 0.}, // 10
801  { 0., 0., 0., 0., 0., 0., 0., 0.375, 0., 0., 0., 0., 0., -0.125, 0., 0., 0.75, 0.}, // 11
802  { 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, 0.375, 0., 0., 0., 0., 0.75, 0., 0.}, // 12
803  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.375, -0.125, 0., 0., 0., 0., 0.75, 0.}, // 13
804  { 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, 0., -0.125, 0., 0., 0., 0.5, 0.5, 0.25}, // 14
805  {-0.046875, 0.140625, 0., 0.015625,-0.046875, 0., 0.28125, 0., 0., -0.09375, 0.28125, 0., -0.09375, 0., 0., 0.5625, 0., 0.}, // 15
806  { 0., 0.140625,-0.046875, 0.,-0.046875, 0.015625, 0., 0.28125, 0., 0., 0.28125, -0.09375, 0., -0.09375, 0., 0., 0.5625, 0.}, // 16
807  {-0.046875, 0.,-0.046875, 0.015625, 0., 0.015625, 0.1875, 0.1875, 0.09375, -0.09375, 0., -0.09375, -0.0625, -0.0625, -0.03125, 0.375, 0.375, 0.1875} // 17
808  },
809 
810  // embedding matrix for child 2
811  {
812  // 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17
813  { 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 0
814  { 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 1
815  { 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 2
816  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1.}, // 3
817  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0.}, // 4
818  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0.}, // 5
819  { -0.125, -0.125, 0., 0., 0., 0., 0.25, 0.5, 0.5, 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 6
820  { 0., -0.125, 0.375, 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 7
821  { -0.125, 0., 0.375, 0., 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 8
822  { 0., 0., 0., 0., 0., 0., 0., 0., 0.375, 0., 0., 0., 0., 0., -0.125, 0., 0., 0.75}, // 9
823  { 0., 0., 0., 0., 0., 0., 0., 0.375, 0., 0., 0., 0., 0., -0.125, 0., 0., 0.75, 0.}, // 10
824  { 0., 0., 0.375, 0., 0., -0.125, 0., 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0., 0.}, // 11
825  { 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, -0.125, 0., 0., 0., 0., 0.25, 0.5, 0.5}, // 12
826  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, 0.375, 0., 0., 0., 0., 0.75, 0.}, // 13
827  { 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, 0., 0.375, 0., 0., 0., 0., 0., 0.75}, // 14
828  {-0.046875,-0.046875, 0., 0.015625, 0.015625, 0., 0.09375, 0.1875, 0.1875, -0.09375, -0.09375, 0., -0.03125, -0.0625, -0.0625, 0.1875, 0.375, 0.375}, // 15
829  { 0.,-0.046875, 0.140625, 0., 0.015625,-0.046875, 0., 0.28125, 0., 0., -0.09375, 0.28125, 0., -0.09375, 0., 0., 0.5625, 0.}, // 16
830  {-0.046875, 0., 0.140625, 0.015625, 0.,-0.046875, 0., 0., 0.28125, -0.09375, 0., 0.28125, 0., 0., -0.09375, 0., 0., 0.5625} // 17
831  },
832 
833  // embedding matrix for child 3
834  {
835  // 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17
836  { 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 0
837  { 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 1
838  { 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 2
839  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0.}, // 3
840  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0.}, // 4
841  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1.}, // 5
842  { -0.125, 0., -0.125, 0., 0., 0., 0.5, 0.5, 0.25, 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 6
843  { -0.125, -0.125, 0., 0., 0., 0., 0.25, 0.5, 0.5, 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 7
844  { 0., -0.125, -0.125, 0., 0., 0., 0.5, 0.25, 0.5, 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 8
845  { 0., 0., 0., 0., 0., 0., 0.375, 0., 0., 0., 0., 0., -0.125, 0., 0., 0.75, 0., 0.}, // 9
846  { 0., 0., 0., 0., 0., 0., 0., 0.375, 0., 0., 0., 0., 0., -0.125, 0., 0., 0.75, 0.}, // 10
847  { 0., 0., 0., 0., 0., 0., 0., 0., 0.375, 0., 0., 0., 0., 0., -0.125, 0., 0., 0.75}, // 11
848  { 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, 0., -0.125, 0., 0., 0., 0.5, 0.5, 0.25}, // 12
849  { 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, -0.125, 0., 0., 0., 0., 0.25, 0.5, 0.5}, // 13
850  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, -0.125, 0., 0., 0., 0.5, 0.25, 0.5}, // 14
851  {-0.046875, 0.,-0.046875, 0.015625, 0., 0.015625, 0.1875, 0.1875, 0.09375, -0.09375, 0., -0.09375, -0.0625, -0.0625, -0.03125, 0.375, 0.375, 0.1875}, // 15
852  {-0.046875,-0.046875, 0., 0.015625, 0.015625, 0., 0.09375, 0.1875, 0.1875, -0.09375, -0.09375, 0., -0.03125, -0.0625, -0.0625, 0.1875, 0.375, 0.375}, // 16
853  { 0.,-0.046875,-0.046875, 0., 0.015625, 0.015625, 0.1875, 0.09375, 0.1875, 0., -0.09375, -0.09375, -0.0625, -0.03125, -0.0625, 0.375, 0.1875, 0.375} // 17
854  },
855 
856  // embedding matrix for child 4
857  {
858  // 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17
859  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0.}, // 0
860  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0.}, // 1
861  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1.}, // 2
862  { 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 3
863  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0.}, // 4
864  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0.}, // 5
865  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.375, -0.125, 0., 0., 0., 0., 0.75, 0., 0.}, // 6
866  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, -0.125, 0., 0., 0., 0.5, 0.25, 0.5}, // 7
867  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.375, 0., -0.125, 0., 0., 0., 0., 0., 0.75}, // 8
868  { -0.125, 0., 0., 0.375, 0., 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0., 0., 0., 0.}, // 9
869  { 0., 0., 0., 0., 0., 0., -0.125, 0., 0., 0., 0., 0., 0.375, 0., 0., 0.75, 0., 0.}, // 10
870  { 0., 0., 0., 0., 0., 0., 0., 0., -0.125, 0., 0., 0., 0., 0., 0.375, 0., 0., 0.75}, // 11
871  { 0., 0., 0., 0.375, -0.125, 0., 0., 0., 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0.}, // 12
872  { 0., 0., 0., 0., -0.125, -0.125, 0., 0., 0., 0., 0., 0., 0.5, 0.25, 0.5, 0., 0., 0.}, // 13
873  { 0., 0., 0., 0.375, 0., -0.125, 0., 0., 0., 0., 0., 0., 0., 0., 0.75, 0., 0., 0.}, // 14
874  {-0.046875, 0.015625, 0., 0.140625,-0.046875, 0., -0.09375, 0., 0., 0.28125, -0.09375, 0., 0.28125, 0., 0., 0.5625, 0., 0.}, // 15
875  { 0., 0.015625, 0.015625, 0.,-0.046875,-0.046875, -0.0625, -0.03125, -0.0625, 0., -0.09375, -0.09375, 0.1875, 0.09375, 0.1875, 0.375, 0.1875, 0.375}, // 16
876  {-0.046875, 0., 0.015625, 0.140625, 0.,-0.046875, 0., 0., -0.09375, 0.28125, 0., -0.09375, 0., 0., 0.28125, 0., 0., 0.5625} // 17
877  },
878 
879  // embedding matrix for child 5
880  {
881  // 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17
882  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0.}, // 0
883  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0.}, // 1
884  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0.}, // 2
885  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0.}, // 3
886  { 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 4
887  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0.}, // 5
888  { 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, 0.375, 0., 0., 0., 0., 0.75, 0., 0.}, // 6
889  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.375, -0.125, 0., 0., 0., 0., 0.75, 0.}, // 7
890  { 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, 0., -0.125, 0., 0., 0., 0.5, 0.5, 0.25}, // 8
891  { 0., 0., 0., 0., 0., 0., -0.125, 0., 0., 0., 0., 0., 0.375, 0., 0., 0.75, 0., 0.}, // 9
892  { 0., -0.125, 0., 0., 0.375, 0., 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0., 0., 0.}, // 10
893  { 0., 0., 0., 0., 0., 0., 0., -0.125, 0., 0., 0., 0., 0., 0.375, 0., 0., 0.75, 0.}, // 11
894  { 0., 0., 0., -0.125, 0.375, 0., 0., 0., 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0.}, // 12
895  { 0., 0., 0., 0., 0.375, -0.125, 0., 0., 0., 0., 0., 0., 0., 0.75, 0., 0., 0., 0.}, // 13
896  { 0., 0., 0., -0.125, 0., -0.125, 0., 0., 0., 0., 0., 0., 0.5, 0.5, 0.25, 0., 0., 0.}, // 14
897  { 0.015625,-0.046875, 0.,-0.046875, 0.140625, 0., -0.09375, 0., 0., -0.09375, 0.28125, 0., 0.28125, 0., 0., 0.5625, 0., 0.}, // 15
898  { 0.,-0.046875, 0.015625, 0., 0.140625,-0.046875, 0., -0.09375, 0., 0., 0.28125, -0.09375, 0., 0.28125, 0., 0., 0.5625, 0.}, // 16
899  { 0.015625, 0., 0.015625,-0.046875, 0.,-0.046875, -0.0625, -0.0625, -0.03125, -0.09375, 0., -0.09375, 0.1875, 0.1875, 0.09375, 0.375, 0.375, 0.1875} // 17
900  },
901 
902  // embedding matrix for child 6
903  {
904  // 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17
905  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1.}, // 0
906  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0.}, // 1
907  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0.}, // 2
908  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0.}, // 3
909  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0.}, // 4
910  { 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.}, // 5
911  { 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, -0.125, 0., 0., 0., 0., 0.25, 0.5, 0.5}, // 6
912  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, 0.375, 0., 0., 0., 0., 0.75, 0.}, // 7
913  { 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, 0., 0.375, 0., 0., 0., 0., 0., 0.75}, // 8
914  { 0., 0., 0., 0., 0., 0., 0., 0., -0.125, 0., 0., 0., 0., 0., 0.375, 0., 0., 0.75}, // 9
915  { 0., 0., 0., 0., 0., 0., 0., -0.125, 0., 0., 0., 0., 0., 0.375, 0., 0., 0.75, 0.}, // 10
916  { 0., 0., -0.125, 0., 0., 0.375, 0., 0., 0., 0., 0., 0.75, 0., 0., 0., 0., 0., 0.}, // 11
917  { 0., 0., 0., -0.125, -0.125, 0., 0., 0., 0., 0., 0., 0., 0.25, 0.5, 0.5, 0., 0., 0.}, // 12
918  { 0., 0., 0., 0., -0.125, 0.375, 0., 0., 0., 0., 0., 0., 0., 0.75, 0., 0., 0., 0.}, // 13
919  { 0., 0., 0., -0.125, 0., 0.375, 0., 0., 0., 0., 0., 0., 0., 0., 0.75, 0., 0., 0.}, // 14
920  { 0.015625, 0.015625, 0.,-0.046875,-0.046875, 0., -0.03125, -0.0625, -0.0625, -0.09375, -0.09375, 0., 0.09375, 0.1875, 0.1875, 0.1875, 0.375, 0.375}, // 15
921  { 0., 0.015625,-0.046875, 0.,-0.046875, 0.140625, 0., -0.09375, 0., 0., -0.09375, 0.28125, 0., 0.28125, 0., 0., 0.5625, 0.}, // 16
922  { 0.015625, 0.,-0.046875,-0.046875, 0., 0.140625, 0., 0., -0.09375, -0.09375, 0., 0.28125, 0., 0., 0.28125, 0., 0., 0.5625} // 17
923  },
924 
925  // embedding matrix for child 7
926  {
927  // 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17
928  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0.}, // 0
929  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0.}, // 1
930  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1.}, // 2
931  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0.}, // 3
932  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0.}, // 4
933  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0.}, // 5
934  { 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, 0., -0.125, 0., 0., 0., 0.5, 0.5, 0.25}, // 6
935  { 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, -0.125, 0., 0., 0., 0., 0.25, 0.5, 0.5}, // 7
936  { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., -0.125, -0.125, 0., 0., 0., 0.5, 0.25, 0.5}, // 8
937  { 0., 0., 0., 0., 0., 0., -0.125, 0., 0., 0., 0., 0., 0.375, 0., 0., 0.75, 0., 0.}, // 9
938  { 0., 0., 0., 0., 0., 0., 0., -0.125, 0., 0., 0., 0., 0., 0.375, 0., 0., 0.75, 0.}, // 10
939  { 0., 0., 0., 0., 0., 0., 0., 0., -0.125, 0., 0., 0., 0., 0., 0.375, 0., 0., 0.75}, // 11
940  { 0., 0., 0., -0.125, 0., -0.125, 0., 0., 0., 0., 0., 0., 0.5, 0.5, 0.25, 0., 0., 0.}, // 12
941  { 0., 0., 0., -0.125, -0.125, 0., 0., 0., 0., 0., 0., 0., 0.25, 0.5, 0.5, 0., 0., 0.}, // 13
942  { 0., 0., 0., 0., -0.125, -0.125, 0., 0., 0., 0., 0., 0., 0.5, 0.25, 0.5, 0., 0., 0.}, // 14
943  { 0.015625, 0., 0.015625,-0.046875, 0.,-0.046875, -0.0625, -0.0625, -0.03125, -0.09375, 0., -0.09375, 0.1875, 0.1875, 0.09375, 0.375, 0.375, 0.1875}, // 15
944  { 0.015625, 0.015625, 0.,-0.046875,-0.046875, 0., -0.03125, -0.0625, -0.0625, -0.09375, -0.09375, 0., 0.09375, 0.1875, 0.1875, 0.1875, 0.375, 0.375}, // 16
945  { 0., 0.015625, 0.015625, 0.,-0.046875,-0.046875, -0.0625, -0.03125, -0.0625, 0., -0.09375, -0.09375, 0.1875, 0.09375, 0.1875, 0.375, 0.1875, 0.375} // 17
946  }
947  };
948 
949 #endif
950 
951 } // namespace libMesh
libMesh::Prism18::is_node_on_edge
virtual bool is_node_on_edge(const unsigned int n, const unsigned int e) const override
Definition: cell_prism18.C:109
libMesh::Prism18::nodes_on_side
virtual std::vector< unsigned int > nodes_on_side(const unsigned int s) const override
Definition: cell_prism18.C:102
libMesh::dof_id_type
uint8_t dof_id_type
Definition: id_types.h:67
libMesh::Prism::n_vertices
virtual unsigned int n_vertices() const override final
Definition: cell_prism.h:84
libMesh::Elem::compute_key
static dof_id_type compute_key(dof_id_type n0)
Definition: elem.h:2683
libMesh::Prism18::num_edges
static const int num_edges
Definition: cell_prism18.h:216
libMesh::Prism18::which_node_am_i
virtual unsigned int which_node_am_i(unsigned int side, unsigned int side_node) const override
Definition: cell_prism18.C:191
libMesh::Prism18::is_node_on_side
virtual bool is_node_on_side(const unsigned int n, const unsigned int s) const override
Definition: cell_prism18.C:92
libMesh::Prism18::nodes_per_edge
static const int nodes_per_edge
Definition: cell_prism18.h:219
libMesh::Prism18::second_order_child_vertex
virtual std::pair< unsigned short int, unsigned short int > second_order_child_vertex(const unsigned int n) const override
Definition: cell_prism18.C:565
libMesh::Prism18::is_face
virtual bool is_face(const unsigned int i) const override
Definition: cell_prism18.C:85
libMesh::IOPackage
IOPackage
libMesh interfaces with several different software packages for the purposes of creating,...
Definition: enum_io_package.h:37
libMesh::Prism18::build_edge_ptr
virtual std::unique_ptr< Elem > build_edge_ptr(const unsigned int i) override
Builds a EDGE3 or INFEDGE2 built coincident with edge i.
Definition: cell_prism18.C:309
libMesh::Prism18::num_sides
static const int num_sides
Definition: cell_prism18.h:215
libMesh::Prism::_second_order_vertex_child_index
static const unsigned short int _second_order_vertex_child_index[18]
Vector that names the child vertex index for each second order node.
Definition: cell_prism.h:164
libMesh::index_range
IntRange< std::size_t > index_range(const std::vector< T > &vec)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:106
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::INVALID_IO_PACKAGE
Definition: enum_io_package.h:48
end
IterBase * end
Also have a polymorphic pointer to the end object, this prevents iterating past the end.
Definition: variant_filter_iterator.h:343
std::sqrt
MetaPhysicL::DualNumber< T, D > sqrt(const MetaPhysicL::DualNumber< T, D > &in)
libMesh::Prism18::build_side_ptr
virtual std::unique_ptr< Elem > build_side_ptr(const unsigned int i, bool proxy=true) override
Builds a QUAD9 or TRI6 built coincident with face i.
Definition: cell_prism18.C:207
libMesh::Prism::n_edges
virtual unsigned int n_edges() const override final
Definition: cell_prism.h:89
libMesh::triple_product
T triple_product(const TypeVector< T > &a, const TypeVector< T > &b, const TypeVector< T > &c)
Definition: type_vector.h:1106
libMesh::SECOND
Definition: enum_order.h:43
libMesh::Prism::_second_order_vertex_child_number
static const unsigned short int _second_order_vertex_child_number[18]
Vector that names a child sharing each second order node.
Definition: cell_prism.h:159
libMesh::Prism18::n_second_order_adjacent_vertices
virtual unsigned int n_second_order_adjacent_vertices(const unsigned int) const override
Definition: cell_prism18.C:485
libMesh::Elem::_nodes
Node ** _nodes
Pointers to the nodes we are connected to.
Definition: elem.h:1743
libMesh::Prism18::_remaining_second_order_adjacent_vertices
static const unsigned short int _remaining_second_order_adjacent_vertices[3][4]
Matrix that tells which vertices define the location of mid-side (or second-order) nodes.
Definition: cell_prism18.h:275
libMesh::Elem::point
const Point & point(const unsigned int i) const
Definition: elem.h:1955
libMesh::Prism18::edge_nodes_map
static const unsigned int edge_nodes_map[num_edges][nodes_per_edge]
This maps the node of the edge to element node numbers.
Definition: cell_prism18.h:231
libMesh::libmesh_assert
libmesh_assert(ctx)
libMesh::Prism18::connectivity
virtual void connectivity(const unsigned int sc, const IOPackage iop, std::vector< dof_id_type > &conn) const override
Definition: cell_prism18.C:318
libMesh::Prism18::n_nodes
virtual unsigned int n_nodes() const override
Definition: cell_prism18.h:95
libMesh::Prism18::default_order
virtual Order default_order() const override
Definition: cell_prism18.C:153
libMesh::Elem::volume
virtual Real volume() const
Definition: elem.C:2617
libMesh::Prism18::num_children
static const int num_children
Definition: cell_prism18.h:217
libMesh::Prism18::side_nodes_map
static const unsigned int side_nodes_map[num_sides][nodes_per_side]
This maps the node of the side to element node numbers.
Definition: cell_prism18.h:225
libMesh::Point
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
libMesh::Prism18::volume
virtual Real volume() const override
A specialization for computing the volume of a Prism18.
Definition: cell_prism18.C:577
libMesh::Elem::key
virtual dof_id_type key() const
Definition: elem.C:410
libMesh::Elem::mapping_type
ElemMappingType mapping_type() const
Definition: elem.h:2524
libMesh::Prism18::is_vertex
virtual bool is_vertex(const unsigned int i) const override
Definition: cell_prism18.C:69
libMesh::Prism18::n_sub_elem
virtual unsigned int n_sub_elem() const override
Definition: cell_prism18.h:100
libMesh::Prism18::_embedding_matrix
static const float _embedding_matrix[num_children][num_nodes][num_nodes]
Matrix that computes new nodal locations/solution values from current nodes/solution.
Definition: cell_prism18.h:261
libMesh::TRI6
Definition: enum_elem_type.h:40
libMesh::LAGRANGE_MAP
Definition: enum_elem_type.h:83
libMesh::Prism18::is_edge
virtual bool is_edge(const unsigned int i) const override
Definition: cell_prism18.C:76
libMesh::Prism::n_sides
virtual unsigned int n_sides() const override final
Definition: cell_prism.h:79
libMesh::Prism::_second_order_adjacent_vertices
static const unsigned short int _second_order_adjacent_vertices[9][2]
Matrix that tells which vertices define the location of mid-side (or second-order) nodes.
Definition: cell_prism.h:154
libMesh::Elem::subdomain_id
subdomain_id_type subdomain_id() const
Definition: elem.h:2069
libMesh::VTK
Definition: enum_io_package.h:42
libMesh::Prism18::nodes_per_side
static const int nodes_per_side
Definition: cell_prism18.h:218
libMesh::Prism18::num_nodes
static const int num_nodes
Geometric constants for Prism18.
Definition: cell_prism18.h:214
libMesh::QUAD9
Definition: enum_elem_type.h:43
libMesh::Prism18::has_affine_map
virtual bool has_affine_map() const override
Definition: cell_prism18.C:120
libMesh::Elem::node_id
dof_id_type node_id(const unsigned int i) const
Definition: elem.h:1977
libMesh::Real
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Definition: libmesh_common.h:121
libMesh::Elem::node_ptr
const Node * node_ptr(const unsigned int i) const
Definition: elem.h:2009
libMesh::TypeVector::relative_fuzzy_equals
bool relative_fuzzy_equals(const TypeVector< T > &rhs, Real tol=TOLERANCE) const
Definition: type_vector.h:1042
libMesh::TECPLOT
Definition: enum_io_package.h:39
libMesh::Prism18::second_order_adjacent_vertex
virtual unsigned short int second_order_adjacent_vertex(const unsigned int n, const unsigned int v) const override
Definition: cell_prism18.C:514
Generated on Sat Jan 25 2020 12:06:59 for libMesh by   doxygen 1.8.16