libMesh
cell_tet4.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/cell_tet4.h"
21 #include "libmesh/edge_edge2.h"
22 #include "libmesh/face_tri3.h"
23 #include "libmesh/tensor_value.h"
24 #include "libmesh/enum_io_package.h"
25 #include "libmesh/enum_order.h"
26 
27 namespace libMesh
28 {
29 
30 
31 
32 // ------------------------------------------------------------
33 // Tet4 class static member initializations
34 const int Tet4::num_nodes;
35 const int Tet4::nodes_per_side;
36 const int Tet4::nodes_per_edge;
37 
39  {
40  {0, 2, 1}, // Side 0
41  {0, 1, 3}, // Side 1
42  {1, 2, 3}, // Side 2
43  {2, 0, 3} // Side 3
44  };
45 
47  {
48  {0, 1}, // Edge 0
49  {1, 2}, // Edge 1
50  {0, 2}, // Edge 2
51  {0, 3}, // Edge 3
52  {1, 3}, // Edge 4
53  {2, 3} // Edge 5
54  };
55 
56 // ------------------------------------------------------------
57 // Tet4 class member functions
58 
59 bool Tet4::is_vertex(const unsigned int libmesh_dbg_var(n)) const
60 {
61  libmesh_assert_not_equal_to (n, invalid_uint);
62  return true;
63 }
64 
65 bool Tet4::is_edge(const unsigned int) const
66 {
67  return false;
68 }
69 
70 bool Tet4::is_face(const unsigned int) const
71 {
72  return false;
73 }
74 
75 bool Tet4::is_node_on_edge(const unsigned int n,
76  const unsigned int e) const
77 {
78  libmesh_assert_less (e, n_edges());
79  return std::find(std::begin(edge_nodes_map[e]),
80  std::end(edge_nodes_map[e]),
81  n) != std::end(edge_nodes_map[e]);
82 }
83 
84 
85 
86 
87 #ifdef LIBMESH_ENABLE_AMR
88 
89 // This function only works if LIBMESH_ENABLE_AMR...
90 bool Tet4::is_child_on_side(const unsigned int c,
91  const unsigned int s) const
92 {
93  // OK, for the Tet4, this is pretty obvious... it is sets of nodes
94  // not equal to the current node. But if we want this algorithm to
95  // be generic and work for Tet10 also it helps to do it this way.
96  const unsigned int nodes_opposite[4][3] =
97  {
98  {1,2,3}, // nodes opposite node 0
99  {0,2,3}, // nodes opposite node 1
100  {0,1,3}, // nodes opposite node 2
101  {0,1,2} // nodes opposite node 3
102  };
103 
104  // Call the base class helper function
105  return Tet::is_child_on_side_helper(c, s, nodes_opposite);
106 }
107 
108 #else
109 
110 bool Tet4::is_child_on_side(const unsigned int /*c*/,
111  const unsigned int /*s*/) const
112 {
113  libmesh_not_implemented();
114  return false;
115 }
116 
117 #endif //LIBMESH_ENABLE_AMR
118 
119 
120 
122 {
123  return this->volume() > tol;
124 }
125 
126 
127 
128 bool Tet4::is_node_on_side(const unsigned int n,
129  const unsigned int s) const
130 {
131  libmesh_assert_less (s, n_sides());
132  return std::find(std::begin(side_nodes_map[s]),
133  std::end(side_nodes_map[s]),
134  n) != std::end(side_nodes_map[s]);
135 }
136 
137 std::vector<unsigned>
138 Tet4::nodes_on_side(const unsigned int s) const
139 {
140  libmesh_assert_less(s, n_sides());
141  return {std::begin(side_nodes_map[s]), std::end(side_nodes_map[s])};
142 }
143 
144 std::vector<unsigned>
145 Tet4::nodes_on_edge(const unsigned int e) const
146 {
147  libmesh_assert_less(e, n_edges());
148  return {std::begin(edge_nodes_map[e]), std::end(edge_nodes_map[e])};
149 }
150 
152 {
153  return FIRST;
154 }
155 
156 std::unique_ptr<Elem> Tet4::build_side_ptr (const unsigned int i)
157 {
158  return this->simple_build_side_ptr<Tri3, Tet4>(i);
159 }
160 
161 
162 
163 void Tet4::build_side_ptr (std::unique_ptr<Elem> & side,
164  const unsigned int i)
165 {
166  this->simple_build_side_ptr<Tet4>(side, i, TRI3);
167 }
168 
169 
170 
171 std::unique_ptr<Elem> Tet4::build_edge_ptr (const unsigned int i)
172 {
173  return this->simple_build_edge_ptr<Edge2,Tet4>(i);
174 }
175 
176 
177 
178 void Tet4::build_edge_ptr (std::unique_ptr<Elem> & edge, const unsigned int i)
179 {
180  this->simple_build_edge_ptr<Tet4>(edge, i, EDGE2);
181 }
182 
183 
184 
185 void Tet4::connectivity(const unsigned int libmesh_dbg_var(sc),
186  const IOPackage iop,
187  std::vector<dof_id_type> & conn) const
188 {
190  libmesh_assert_less (sc, this->n_sub_elem());
191  libmesh_assert_not_equal_to (iop, INVALID_IO_PACKAGE);
192 
193 
194  switch (iop)
195  {
196  case TECPLOT:
197  {
198  conn.resize(8);
199  conn[0] = this->node_id(0)+1;
200  conn[1] = this->node_id(1)+1;
201  conn[2] = this->node_id(2)+1;
202  conn[3] = this->node_id(2)+1;
203  conn[4] = this->node_id(3)+1;
204  conn[5] = this->node_id(3)+1;
205  conn[6] = this->node_id(3)+1;
206  conn[7] = this->node_id(3)+1;
207  return;
208  }
209 
210  case VTK:
211  {
212  conn.resize(4);
213  conn[0] = this->node_id(0);
214  conn[1] = this->node_id(1);
215  conn[2] = this->node_id(2);
216  conn[3] = this->node_id(3);
217  return;
218  }
219 
220  default:
221  libmesh_error_msg("Unsupported IO package " << iop);
222  }
223 }
224 
225 
226 
227 #ifdef LIBMESH_ENABLE_AMR
228 
230  {
231  // embedding matrix for child 0
232  {
233  // 0 1 2 3
234  {1.0, 0.0, 0.0, 0.0}, // 0
235  {0.5, 0.5, 0.0, 0.0}, // 1
236  {0.5, 0.0, 0.5, 0.0}, // 2
237  {0.5, 0.0, 0.0, 0.5} // 3
238  },
239 
240  // embedding matrix for child 1
241  {
242  // 0 1 2 3
243  {0.5, 0.5, 0.0, 0.0}, // 0
244  {0.0, 1.0, 0.0, 0.0}, // 1
245  {0.0, 0.5, 0.5, 0.0}, // 2
246  {0.0, 0.5, 0.0, 0.5} // 3
247  },
248 
249  // embedding matrix for child 2
250  {
251  // 0 1 2 3
252  {0.5, 0.0, 0.5, 0.0}, // 0
253  {0.0, 0.5, 0.5, 0.0}, // 1
254  {0.0, 0.0, 1.0, 0.0}, // 2
255  {0.0, 0.0, 0.5, 0.5} // 3
256  },
257 
258  // embedding matrix for child 3
259  {
260  // 0 1 2 3
261  {0.5, 0.0, 0.0, 0.5}, // 0
262  {0.0, 0.5, 0.0, 0.5}, // 1
263  {0.0, 0.0, 0.5, 0.5}, // 2
264  {0.0, 0.0, 0.0, 1.0} // 3
265  },
266 
267  // embedding matrix for child 4
268  {
269  // 0 1 2 3
270  {0.5, 0.5, 0.0, 0.0}, // 0
271  {0.0, 0.5, 0.0, 0.5}, // 1
272  {0.5, 0.0, 0.5, 0.0}, // 2
273  {0.5, 0.0, 0.0, 0.5} // 3
274  },
275 
276  // embedding matrix for child 5
277  {
278  // 0 1 2 3
279  {0.5, 0.5, 0.0, 0.0}, // 0
280  {0.0, 0.5, 0.5, 0.0}, // 1
281  {0.5, 0.0, 0.5, 0.0}, // 2
282  {0.0, 0.5, 0.0, 0.5} // 3
283  },
284 
285  // embedding matrix for child 6
286  {
287  // 0 1 2 3
288  {0.5, 0.0, 0.5, 0.0}, // 0
289  {0.0, 0.5, 0.5, 0.0}, // 1
290  {0.0, 0.0, 0.5, 0.5}, // 2
291  {0.0, 0.5, 0.0, 0.5} // 3
292  },
293 
294  // embedding matrix for child 7
295  {
296  // 0 1 2 3
297  {0.5, 0.0, 0.5, 0.0}, // 0
298  {0.0, 0.5, 0.0, 0.5}, // 1
299  {0.0, 0.0, 0.5, 0.5}, // 2
300  {0.5, 0.0, 0.0, 0.5} // 3
301  }
302  };
303 
304 #endif // #ifdef LIBMESH_ENABLE_AMR
305 
306 
307 
308 
309 
311 {
312  return Elem::vertex_average();
313 }
314 
315 
316 
318 {
319  // The volume of a tetrahedron is 1/6 the box product formed
320  // by its base and apex vectors
321  Point a = point(3) - point(0);
322 
323  // b is the vector pointing from 0 to 1
324  Point b = point(1) - point(0);
325 
326  // c is the vector pointing from 0 to 2
327  Point c = point(2) - point(0);
328 
329  return triple_product(a, b, c) / 6.;
330 }
331 
332 
333 
334 
335 std::pair<Real, Real> Tet4::min_and_max_angle() const
336 {
337  Point n[4];
338 
339  // Compute the outward normal vectors on each face
340  n[0] = (this->point(2) - this->point(0)).cross(this->point(1) - this->point(0));
341  n[1] = (this->point(1) - this->point(0)).cross(this->point(3) - this->point(0));
342  n[2] = (this->point(2) - this->point(1)).cross(this->point(3) - this->point(1));
343  n[3] = (this->point(0) - this->point(2)).cross(this->point(3) - this->point(2));
344 
345  Real dihedral_angles[6]; // 01, 02, 03, 12, 13, 23
346 
347  // Compute dihedral angles
348  for (unsigned int k=0,i=0; i<4; ++i)
349  for (unsigned int j=i+1; j<4; ++j,k+=1)
350  dihedral_angles[k] = std::acos(n[i]*n[j] / n[i].norm() / n[j].norm()); // return value is between 0 and PI
351 
352  // Return max/min dihedral angles
353  return std::make_pair(*std::min_element(dihedral_angles, dihedral_angles+6),
354  *std::max_element(dihedral_angles, dihedral_angles+6));
355 
356 }
357 
358 
359 
361 {
362  return this->compute_key(this->node_id(0),
363  this->node_id(1),
364  this->node_id(2),
365  this->node_id(3));
366 }
367 
368 
369 
370 bool Tet4::contains_point (const Point & p, Real tol) const
371 {
372  // See the response by Tony Noe on this thread.
373  // http://bbs.dartmouth.edu/~fangq/MATH/download/source/Point_in_Tetrahedron.htm
374  Point
375  col1 = point(1) - point(0),
376  col2 = point(2) - point(0),
377  col3 = point(3) - point(0);
378 
379  Point r;
380 
381  libmesh_try
382  {
383  RealTensorValue(col1(0), col2(0), col3(0),
384  col1(1), col2(1), col3(1),
385  col1(2), col2(2), col3(2)).solve(p - point(0), r);
386  }
387  libmesh_catch (ConvergenceFailure &)
388  {
389  // The Jacobian was singular, therefore the Tet had
390  // zero volume. In this degenerate case, it is not
391  // possible for the Tet to contain any points.
392  return false;
393  }
394 
395  // The point p is inside the tetrahedron if r1 + r2 + r3 < 1 and
396  // 0 <= ri <= 1 for i = 1,2,3.
397  return
398  r(0) > -tol &&
399  r(1) > -tol &&
400  r(2) > -tol &&
401  r(0) + r(1) + r(2) < 1.0 + tol;
402 }
403 
404 
405 
406 #ifdef LIBMESH_ENABLE_AMR
407 Real Tet4::embedding_matrix (const unsigned int i,
408  const unsigned int j,
409  const unsigned int k) const
410 {
411  // Choose an optimal diagonal, if one has not already been selected
412  this->choose_diagonal();
413 
414  // Permuted j and k indices
415  unsigned int
416  jp=j,
417  kp=k;
418 
419  if ((i>3) && (this->_diagonal_selection!=DIAG_02_13))
420  {
421  // Just the enum value cast to an unsigned int...
422  const unsigned ds = static_cast<unsigned int>(this->_diagonal_selection);
423 
424  // Permute j, k:
425  // ds==1 ds==2
426  // 0 -> 1 0 -> 2
427  // 1 -> 2 1 -> 0
428  // 2 -> 0 2 -> 1
429  if (jp != 3)
430  jp = (jp+ds)%3;
431 
432  if (kp != 3)
433  kp = (kp+ds)%3;
434  }
435 
436  // Debugging
437  // libMesh::err << "Selected diagonal " << _diagonal_selection << std::endl;
438  // libMesh::err << "j=" << j << std::endl;
439  // libMesh::err << "k=" << k << std::endl;
440  // libMesh::err << "jp=" << jp << std::endl;
441  // libMesh::err << "kp=" << kp << std::endl;
442 
443  // Call embedding matrix with permuted indices
444  return this->_embedding_matrix[i][jp][kp];
445 }
446 
447 
448 
449 // void Tet4::reselect_diagonal (const Diagonal diag)
450 // {
451 // /* Make sure that the element has just been refined. */
452 // libmesh_assert(_children);
453 // libmesh_assert_equal_to (n_children(), 8);
454 // libmesh_assert_equal_to (_children[0]->refinement_flag(), JUST_REFINED);
455 // libmesh_assert_equal_to (_children[1]->refinement_flag(), JUST_REFINED);
456 // libmesh_assert_equal_to (_children[2]->refinement_flag(), JUST_REFINED);
457 // libmesh_assert_equal_to (_children[3]->refinement_flag(), JUST_REFINED);
458 // libmesh_assert_equal_to (_children[4]->refinement_flag(), JUST_REFINED);
459 // libmesh_assert_equal_to (_children[5]->refinement_flag(), JUST_REFINED);
460 // libmesh_assert_equal_to (_children[6]->refinement_flag(), JUST_REFINED);
461 // libmesh_assert_equal_to (_children[7]->refinement_flag(), JUST_REFINED);
462 //
463 // /* Check whether anything has to be changed. */
464 // if (_diagonal_selection!=diag)
465 // {
466 // /* Set new diagonal selection. */
467 // _diagonal_selection = diag;
468 //
469 // /* The first four children do not have to be changed. For the
470 // others, only the nodes have to be changed. Note that we have
471 // to keep track of the nodes ourselves since there is no \p
472 // MeshRefinement object with a valid \p _new_nodes_map
473 // available. */
474 // for (unsigned int c=4; c<this->n_children(); c++)
475 // {
476 // Elem * child = this->child_ptr(c);
477 // for (unsigned int nc=0; nc<child->n_nodes(); nc++)
478 // {
479 // /* Unassign the current node. */
480 // child->set_node(nc, nullptr);
481 //
482 // /* We have to find the correct new node now. We know
483 // that it exists somewhere. We make use of the fact
484 // that the embedding matrix for these children consists
485 // of entries 0.0 and 0.5 only. Also, we make use of
486 // the properties of the embedding matrix for the first
487 // (unchanged) four children, which allow us to use a
488 // simple mechanism to find the required node. */
489 //
490 //
491 // unsigned int first_05_in_embedding_matrix = libMesh::invalid_uint;
492 // for (unsigned int n=0; n<this->n_nodes(); n++)
493 // {
494 // if (this->embedding_matrix(c,nc,n) != 0.0)
495 // {
496 // /* It must be 0.5 then. Check whether it's the
497 // first or second time that we get a 0.5
498 // value. */
499 // if (first_05_in_embedding_matrix==libMesh::invalid_uint)
500 // {
501 // /* First time, so just remember this position. */
502 // first_05_in_embedding_matrix = n;
503 // }
504 // else
505 // {
506 // /* Second time, so we know now which node to
507 // use. */
508 // child->set_node(nc, this->child_ptr(n)->node_ptr(first_05_in_embedding_matrix));
509 // }
510 //
511 // }
512 // }
513 //
514 // /* Make sure that a node has been found. */
515 // libmesh_assert(child->node_ptr(nc));
516 // }
517 // }
518 // }
519 // }
520 
521 
522 
523 // void Tet4::reselect_optimal_diagonal (const Diagonal exclude_this)
524 // {
525 // Real diag_01_23 = (this->point(0)+this->point(1)-this->point(2)-this->point(3)).norm_sq();
526 // Real diag_02_13 = (this->point(0)-this->point(1)+this->point(2)-this->point(3)).norm_sq();
527 // Real diag_03_12 = (this->point(0)-this->point(1)-this->point(2)+this->point(3)).norm_sq();
528 //
529 // Diagonal use_this = INVALID_DIAG;
530 // switch (exclude_this)
531 // {
532 // case DIAG_01_23:
533 // use_this = DIAG_02_13;
534 // if (diag_03_12 < diag_02_13)
535 // {
536 // use_this = DIAG_03_12;
537 // }
538 // break;
539 //
540 // case DIAG_02_13:
541 // use_this = DIAG_03_12;
542 // if (diag_01_23 < diag_03_12)
543 // {
544 // use_this = DIAG_01_23;
545 // }
546 // break;
547 //
548 // case DIAG_03_12:
549 // use_this = DIAG_02_13;
550 // if (diag_01_23 < diag_02_13)
551 // {
552 // use_this = DIAG_01_23;
553 // }
554 // break;
555 //
556 // default:
557 // use_this = DIAG_02_13;
558 // if (diag_01_23 < diag_02_13 || diag_03_12 < diag_02_13)
559 // {
560 // if (diag_01_23 < diag_03_12)
561 // {
562 // use_this = DIAG_01_23;
563 // }
564 // else
565 // {
566 // use_this = DIAG_03_12;
567 // }
568 // }
569 // break;
570 // }
571 //
572 // reselect_diagonal (use_this);
573 // }
574 #endif // #ifdef LIBMESH_ENABLE_AMR
575 
576 void Tet4::permute(unsigned int perm_num)
577 {
578  libmesh_assert_less (perm_num, 12);
579 
580  const unsigned int side = perm_num % 4;
581  const unsigned int rotate = perm_num / 4;
582 
583  for (unsigned int i = 0; i != rotate; ++i)
584  {
585  swap3nodes(0,1,2);
586  swap3neighbors(1,2,3);
587  }
588 
589  switch (side) {
590  case 0:
591  break;
592  case 1:
593  swap3nodes(0,2,3);
594  swap3neighbors(0,2,1);
595  break;
596  case 2:
597  swap3nodes(2,0,3);
598  swap3neighbors(0,1,2);
599  break;
600  case 3:
601  swap3nodes(2,1,3);
602  swap3neighbors(0,1,3);
603  break;
604  default:
605  libmesh_error();
606  }
607 }
608 
609 
610 void Tet4::flip(BoundaryInfo * boundary_info)
611 {
612  libmesh_assert(boundary_info);
613 
614  swap2nodes(0,2);
615  swap2neighbors(1,2);
616  swap2boundarysides(1,2,boundary_info);
617  swap2boundaryedges(0,1,boundary_info);
618  swap2boundaryedges(3,5,boundary_info);
619 }
620 
621 
622 ElemType Tet4::side_type (const unsigned int libmesh_dbg_var(s)) const
623 {
624  libmesh_assert_less (s, 4);
625  return TRI3;
626 }
627 
628 
629 } // namespace libMesh
virtual void permute(unsigned int perm_num) override final
Permutes the element (by swapping node and neighbor pointers) according to the specified index...
Definition: cell_tet4.C:576
ElemType
Defines an enum for geometric element types.
void swap2boundaryedges(unsigned short e1, unsigned short e2, BoundaryInfo *boundary_info) const
Swaps two edges in boundary_info, if it is non-null.
Definition: elem.C:3550
Order
defines an enum for polynomial orders.
Definition: enum_order.h:40
virtual std::vector< unsigned int > nodes_on_edge(const unsigned int e) const override
Definition: cell_tet4.C:145
Node ** _nodes
Pointers to the nodes we are connected to.
Definition: elem.h:2245
virtual Point true_centroid() const override
The centroid of a 4-node tetrahedron is simply given by the average of its vertex positions...
Definition: cell_tet4.C:310
std::pair< Real, Real > min_and_max_angle() const
Definition: cell_tet4.C:335
auto norm() const -> decltype(std::norm(T()))
Definition: type_vector.h:907
const unsigned int invalid_uint
A number which is used quite often to represent an invalid or uninitialized value for an unsigned int...
Definition: libmesh.h:310
virtual std::vector< unsigned int > nodes_on_side(const unsigned int s) const override
Definition: cell_tet4.C:138
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_tet4.h:195
virtual unsigned int n_sides() const override final
Definition: cell_tet.h:81
static const int num_sides
Geometric constants for all Tets.
Definition: cell_tet.h:74
IOPackage
libMesh interfaces with several different software packages for the purposes of creating, reading, and writing mesh files.
void swap2boundarysides(unsigned short s1, unsigned short s2, BoundaryInfo *boundary_info) const
Swaps two sides in boundary_info, if it is non-null.
Definition: elem.C:3534
static const int num_edges
Definition: cell_tet.h:75
virtual bool is_edge(const unsigned int i) const override
Definition: cell_tet4.C:65
Diagonal _diagonal_selection
The currently-selected diagonal used during refinement.
Definition: cell_tet.h:249
virtual unsigned int n_sub_elem() const override
Definition: cell_tet4.h:90
The libMesh namespace provides an interface to certain functionality in the library.
virtual void connectivity(const unsigned int sc, const IOPackage iop, std::vector< dof_id_type > &conn) const override
Definition: cell_tet4.C:185
virtual bool is_vertex(const unsigned int i) const override
Definition: cell_tet4.C:59
virtual dof_id_type key() const override
Definition: cell_tet4.C:360
void swap3neighbors(unsigned int n1, unsigned int n2, unsigned int n3)
Swaps three neighbor_ptrs, "rotating" them.
Definition: elem.h:2133
TensorValue< Real > RealTensorValue
Useful typedefs to allow transparent switching between Real and Complex data types.
void swap3nodes(unsigned int n1, unsigned int n2, unsigned int n3)
Swaps three node_ptrs, "rotating" them.
Definition: elem.h:2124
T triple_product(const TypeVector< T > &a, const TypeVector< T > &b, const TypeVector< T > &c)
Definition: type_vector.h:1029
static const int nodes_per_side
Definition: cell_tet4.h:188
static const int num_nodes
Geometric constants for Tet4.
Definition: cell_tet4.h:187
void swap2nodes(unsigned int n1, unsigned int n2)
Swaps two node_ptrs.
Definition: elem.h:2092
virtual Order default_order() const override
Definition: cell_tet4.C:151
static const Real _embedding_matrix[num_children][num_nodes][num_nodes]
Matrix that computes new nodal locations/solution values from current nodes/solution.
Definition: cell_tet4.h:272
The BoundaryInfo class contains information relevant to boundary conditions including storing faces...
Definition: boundary_info.h:57
libmesh_assert(ctx)
virtual unsigned int n_edges() const override final
Definition: cell_tet.h:91
ElemType side_type(const unsigned int s) const override final
Definition: cell_tet4.C:622
virtual void flip(BoundaryInfo *) override final
Flips the element (by swapping node and neighbor pointers) to have a mapping Jacobian of opposite sig...
Definition: cell_tet4.C:610
virtual bool is_node_on_edge(const unsigned int n, const unsigned int e) const override
Definition: cell_tet4.C:75
virtual std::unique_ptr< Elem > build_side_ptr(const unsigned int i) override
Builds a TRI3 built coincident with face i.
Definition: cell_tet4.C:156
auto norm(const T &a) -> decltype(std::abs(a))
Definition: tensor_tools.h:74
RealTensorValue rotate(MeshBase &mesh, const Real phi, const Real theta=0., const Real psi=0.)
Rotates the mesh in the xy plane.
virtual bool is_face(const unsigned int i) const override
Definition: cell_tet4.C:70
virtual bool has_invertible_map(Real tol) const override
Definition: cell_tet4.C:121
void swap2neighbors(unsigned int n1, unsigned int n2)
Swaps two neighbor_ptrs.
Definition: elem.h:2102
virtual bool is_node_on_side(const unsigned int n, const unsigned int s) const override
Definition: cell_tet4.C:128
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual bool is_child_on_side(const unsigned int c, const unsigned int s) const override
Definition: cell_tet4.C:90
virtual Real volume() const override
An optimized method for computing the area of a 4-node tetrahedron.
Definition: cell_tet4.C:317
A class representing a solver&#39;s failure to converge, to be thrown by "libmesh_convergence_failure();"...
virtual bool contains_point(const Point &p, Real tol) const override
Uses simple geometric tests to determine if the point p is inside the tetrahedron.
Definition: cell_tet4.C:370
void choose_diagonal() const
Derived classes use this function to select an initial diagonal during refinement.
Definition: cell_tet.C:204
static const int num_children
Definition: cell_tet.h:76
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_tet4.h:201
static dof_id_type compute_key(dof_id_type n0)
Definition: elem.h:3294
virtual std::unique_ptr< Elem > build_edge_ptr(const unsigned int i) override
Builds a EDGE2 built coincident with face i.
Definition: cell_tet4.C:171
bool is_child_on_side_helper(const unsigned int c, const unsigned int s, const unsigned int checked_nodes[][3]) const
Called by descendant classes with appropriate data to determine if child c is on side s...
Definition: cell_tet.C:152
virtual Real embedding_matrix(const unsigned int i, const unsigned int j, const unsigned int k) const override
Matrix used to create the elements children.
Definition: cell_tet4.C:407
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39
dof_id_type node_id(const unsigned int i) const
Definition: elem.h:2475
const Point & point(const unsigned int i) const
Definition: elem.h:2453
static const int nodes_per_edge
Definition: cell_tet4.h:189
Point vertex_average() const
Definition: elem.C:688
uint8_t dof_id_type
Definition: id_types.h:67