Line data Source code
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 : // C++ includes
20 : #include <algorithm> // for std::min, std::max
21 :
22 : // Local includes
23 : #include "libmesh/cell_hex.h"
24 : #include "libmesh/cell_hex8.h"
25 : #include "libmesh/face_quad4.h"
26 : #include "libmesh/enum_elem_quality.h"
27 : #include "libmesh/tensor_value.h"
28 :
29 : #include <array>
30 :
31 : namespace libMesh
32 : {
33 :
34 :
35 :
36 : // ------------------------------------------------------------
37 : // Hex class static member initializations
38 : const int Hex::num_sides;
39 : const int Hex::num_edges;
40 : const int Hex::num_children;
41 :
42 : const Real Hex::_master_points[27][3] =
43 : {
44 : {-1, -1, -1},
45 : {1, -1, -1},
46 : {1, 1, -1},
47 : {-1, 1, -1},
48 : {-1, -1, 1},
49 : {1, -1, 1},
50 : {1, 1, 1},
51 : {-1, 1, 1},
52 : {0, -1, -1},
53 : {1, 0, -1},
54 : {0, 1, -1},
55 : {-1, 0, -1},
56 : {-1, -1, 0},
57 : {1, -1, 0},
58 : {1, 1, 0},
59 : {-1, 1, 0},
60 : {0, -1, 1},
61 : {1, 0, 1},
62 : {0, 1, 1},
63 : {-1, 0, 1},
64 : {0, 0, -1},
65 : {0, -1, 0},
66 : {1, 0, 0},
67 : {0, 1, 0},
68 : {-1, 0, 0},
69 : {0, 0, 1}
70 : };
71 :
72 : const unsigned int Hex::edge_sides_map[12][2] =
73 : {
74 : {0, 1}, // Edge 0
75 : {0, 2}, // Edge 1
76 : {0, 3}, // Edge 2
77 : {0, 4}, // Edge 3
78 : {1, 4}, // Edge 4
79 : {1, 2}, // Edge 5
80 : {2, 3}, // Edge 6
81 : {3, 4}, // Edge 7
82 : {1, 5}, // Edge 8
83 : {2, 5}, // Edge 9
84 : {3, 5}, // Edge 10
85 : {4, 5} // Edge 11
86 : };
87 :
88 : const unsigned int Hex::adjacent_edges_map[/*num_vertices*/8][/*n_adjacent_edges*/3] =
89 : {
90 : {0, 3, 4}, // Edges adjacent to node 0
91 : {0, 1, 5}, // Edges adjacent to node 1
92 : {1, 2, 6}, // Edges adjacent to node 2
93 : {2, 3, 7}, // Edges adjacent to node 3
94 : {4, 8, 11}, // Edges adjacent to node 4
95 : {5, 8, 9}, // Edges adjacent to node 5
96 : {6, 9, 10}, // Edges adjacent to node 6
97 : {7, 10, 11} // Edges adjacent to node 7
98 : };
99 :
100 : // ------------------------------------------------------------
101 : // Hex class member functions
102 144576 : dof_id_type Hex::key (const unsigned int s) const
103 : {
104 12048 : libmesh_assert_less (s, this->n_sides());
105 :
106 156624 : return this->compute_key(this->node_id(Hex8::side_nodes_map[s][0]),
107 144576 : this->node_id(Hex8::side_nodes_map[s][1]),
108 144576 : this->node_id(Hex8::side_nodes_map[s][2]),
109 156624 : this->node_id(Hex8::side_nodes_map[s][3]));
110 : }
111 :
112 :
113 :
114 32354539 : dof_id_type Hex::low_order_key (const unsigned int s) const
115 : {
116 899692 : libmesh_assert_less (s, this->n_sides());
117 :
118 33254219 : return this->compute_key(this->node_id(Hex8::side_nodes_map[s][0]),
119 32354539 : this->node_id(Hex8::side_nodes_map[s][1]),
120 32354539 : this->node_id(Hex8::side_nodes_map[s][2]),
121 33254231 : this->node_id(Hex8::side_nodes_map[s][3]));
122 : }
123 :
124 :
125 :
126 3970039 : unsigned int Hex::local_side_node(unsigned int side,
127 : unsigned int side_node) const
128 : {
129 330322 : libmesh_assert_less (side, this->n_sides());
130 330322 : libmesh_assert_less (side_node, Hex8::nodes_per_side);
131 :
132 3970039 : return Hex8::side_nodes_map[side][side_node];
133 : }
134 :
135 :
136 :
137 157743520 : unsigned int Hex::local_edge_node(unsigned int edge,
138 : unsigned int edge_node) const
139 : {
140 11805302 : libmesh_assert_less (edge, this->n_edges());
141 11805302 : libmesh_assert_less (edge_node, Hex8::nodes_per_edge);
142 :
143 157743520 : return Hex8::edge_nodes_map[edge][edge_node];
144 : }
145 :
146 :
147 :
148 546721277 : std::unique_ptr<Elem> Hex::side_ptr (const unsigned int i)
149 : {
150 44766268 : libmesh_assert_less (i, this->n_sides());
151 :
152 546721277 : std::unique_ptr<Elem> face = std::make_unique<Quad4>();
153 :
154 2733606385 : for (auto n : face->node_index_range())
155 2365950180 : face->set_node(n, this->node_ptr(Hex8::side_nodes_map[i][n]));
156 :
157 546721277 : return face;
158 0 : }
159 :
160 :
161 :
162 28156828 : void Hex::side_ptr (std::unique_ptr<Elem> & side,
163 : const unsigned int i)
164 : {
165 28156828 : this->simple_side_ptr<Hex,Hex8>(side, i, QUAD4);
166 28156828 : }
167 :
168 :
169 :
170 3509190 : bool Hex::is_child_on_side(const unsigned int c,
171 : const unsigned int s) const
172 : {
173 1285960 : libmesh_assert_less (c, this->n_children());
174 1285960 : libmesh_assert_less (s, this->n_sides());
175 :
176 : // This array maps the Hex8 node numbering to the Hex8 child
177 : // numbering. I.e.
178 : // node 6 touches child 7, and
179 : // node 7 touches child 6, etc.
180 3509190 : const unsigned int node_child_map[8] = { 0, 1, 3, 2, 4, 5, 7, 6 };
181 :
182 13524906 : for (unsigned int i = 0; i != 4; ++i)
183 11624470 : if (node_child_map[Hex8::side_nodes_map[s][i]] == c)
184 585656 : return true;
185 :
186 700304 : return false;
187 : }
188 :
189 :
190 :
191 18727108 : bool Hex::is_edge_on_side(const unsigned int e,
192 : const unsigned int s) const
193 : {
194 11132724 : libmesh_assert_less (e, this->n_edges());
195 11132724 : libmesh_assert_less (s, this->n_sides());
196 :
197 18727108 : return (edge_sides_map[e][0] == s || edge_sides_map[e][1] == s);
198 : }
199 :
200 :
201 :
202 13824 : std::vector<unsigned int> Hex::sides_on_edge(const unsigned int e) const
203 : {
204 1152 : libmesh_assert_less(e, this->n_edges());
205 :
206 13824 : return {edge_sides_map[e][0], edge_sides_map[e][1]};
207 : }
208 :
209 :
210 :
211 0 : unsigned int Hex::opposite_side(const unsigned int side_in) const
212 : {
213 0 : libmesh_assert_less (side_in, 6);
214 : static const unsigned char hex_opposites[6] = {5, 3, 4, 1, 2, 0};
215 0 : return hex_opposites[side_in];
216 : }
217 :
218 :
219 :
220 0 : unsigned int Hex::opposite_node(const unsigned int node_in,
221 : const unsigned int side_in) const
222 : {
223 0 : libmesh_assert_less (node_in, 26);
224 0 : libmesh_assert_less (node_in, this->n_nodes());
225 0 : libmesh_assert_less (side_in, this->n_sides());
226 0 : libmesh_assert(this->is_node_on_side(node_in, side_in));
227 :
228 : static const unsigned char side05_nodes_map[] =
229 : {4, 5, 6, 7, 0, 1, 2, 3, 16, 17, 18, 19, 255, 255, 255, 255, 8, 9, 10, 11, 25, 255, 255, 255, 255, 20};
230 : static const unsigned char side13_nodes_map[] =
231 : {3, 2, 1, 0, 7, 6, 5, 4, 10, 255, 8, 255, 15, 14, 13, 12, 18, 255, 16, 255, 255, 23, 255, 21, 255, 255};
232 : static const unsigned char side24_nodes_map[] =
233 : {1, 0, 3, 2, 5, 4, 7, 6, 255, 11, 255, 9, 13, 12, 15, 14, 255, 19, 255, 17, 255, 255, 24, 255, 22, 255};
234 :
235 0 : switch (side_in)
236 : {
237 0 : case 0:
238 : case 5:
239 0 : return side05_nodes_map[node_in];
240 0 : case 1:
241 : case 3:
242 0 : return side13_nodes_map[node_in];
243 0 : case 2:
244 : case 4:
245 0 : return side24_nodes_map[node_in];
246 0 : default:
247 0 : libmesh_error_msg("Unsupported side_in = " << side_in);
248 : }
249 : }
250 :
251 :
252 :
253 : bool
254 2760 : Hex::is_flipped() const
255 : {
256 2592 : return (triple_product(this->point(1)-this->point(0),
257 2592 : this->point(3)-this->point(0),
258 3096 : this->point(4)-this->point(0)) < 0);
259 : }
260 :
261 :
262 : std::vector<unsigned int>
263 26400 : Hex::edges_adjacent_to_node(const unsigned int n) const
264 : {
265 2200 : libmesh_assert_less(n, this->n_nodes());
266 :
267 : // For vertices, we use the Hex::adjacent_edges_map, otherwise each
268 : // of the mid-edge nodes is adjacent only to the edge it is on, and
269 : // face/internal nodes are not adjacent to any edge.
270 26400 : if (this->is_vertex(n))
271 12480 : return {std::begin(adjacent_edges_map[n]), std::end(adjacent_edges_map[n])};
272 14880 : else if (this->is_edge(n))
273 11520 : return {n - this->n_vertices()};
274 :
275 280 : libmesh_assert(this->is_face(n) || this->is_internal(n));
276 3080 : return {};
277 : }
278 :
279 :
280 2016 : Real Hex::quality (const ElemQuality q) const
281 : {
282 2016 : switch (q)
283 : {
284 :
285 : #if LIBMESH_DIM >= 3
286 : /**
287 : * Compute the min/max diagonal ratio.
288 : * Source: CUBIT User's Manual.
289 : */
290 0 : case DIAGONAL:
291 : {
292 : // Diagonal between node 0 and node 6
293 0 : const Real d06 = this->length(0,6);
294 :
295 : // Diagonal between node 3 and node 5
296 0 : const Real d35 = this->length(3,5);
297 :
298 : // Diagonal between node 1 and node 7
299 0 : const Real d17 = this->length(1,7);
300 :
301 : // Diagonal between node 2 and node 4
302 0 : const Real d24 = this->length(2,4);
303 :
304 : // Find the biggest and smallest diagonals
305 0 : const Real min = std::min(d06, std::min(d35, std::min(d17, d24)));
306 0 : const Real max = std::max(d06, std::max(d35, std::max(d17, d24)));
307 :
308 0 : libmesh_assert_not_equal_to (max, 0.0);
309 :
310 0 : return min / max;
311 :
312 : break;
313 : }
314 :
315 : /**
316 : * Minimum ratio of lengths derived from opposite edges.
317 : * Source: CUBIT User's Manual.
318 : */
319 0 : case TAPER:
320 : {
321 :
322 : /**
323 : * Compute the side lengths.
324 : */
325 0 : const Real d01 = this->length(0,1);
326 0 : const Real d12 = this->length(1,2);
327 0 : const Real d23 = this->length(2,3);
328 0 : const Real d03 = this->length(0,3);
329 0 : const Real d45 = this->length(4,5);
330 0 : const Real d56 = this->length(5,6);
331 0 : const Real d67 = this->length(6,7);
332 0 : const Real d47 = this->length(4,7);
333 0 : const Real d04 = this->length(0,4);
334 0 : const Real d15 = this->length(1,5);
335 0 : const Real d37 = this->length(3,7);
336 0 : const Real d26 = this->length(2,6);
337 :
338 0 : std::vector<Real> edge_ratios(12);
339 : // Front
340 0 : edge_ratios[0] = std::min(d01, d45) / std::max(d01, d45);
341 0 : edge_ratios[1] = std::min(d04, d15) / std::max(d04, d15);
342 :
343 : // Right
344 0 : edge_ratios[2] = std::min(d15, d26) / std::max(d15, d26);
345 0 : edge_ratios[3] = std::min(d12, d56) / std::max(d12, d56);
346 :
347 : // Back
348 0 : edge_ratios[4] = std::min(d67, d23) / std::max(d67, d23);
349 0 : edge_ratios[5] = std::min(d26, d37) / std::max(d26, d37);
350 :
351 : // Left
352 0 : edge_ratios[6] = std::min(d04, d37) / std::max(d04, d37);
353 0 : edge_ratios[7] = std::min(d03, d47) / std::max(d03, d47);
354 :
355 : // Bottom
356 0 : edge_ratios[8] = std::min(d01, d23) / std::max(d01, d23);
357 0 : edge_ratios[9] = std::min(d03, d12) / std::max(d03, d12);
358 :
359 : // Top
360 0 : edge_ratios[10] = std::min(d45, d67) / std::max(d45, d67);
361 0 : edge_ratios[11] = std::min(d56, d47) / std::max(d56, d47);
362 :
363 0 : return *(std::min_element(edge_ratios.begin(), edge_ratios.end())) ;
364 :
365 : break;
366 : }
367 :
368 :
369 : /**
370 : * Minimum edge length divided by max diagonal length.
371 : * Source: CUBIT User's Manual.
372 : */
373 0 : case STRETCH:
374 : {
375 0 : const Real sqrt3 = 1.73205080756888;
376 :
377 : /**
378 : * Compute the maximum diagonal.
379 : */
380 0 : const Real d06 = this->length(0,6);
381 0 : const Real d17 = this->length(1,7);
382 0 : const Real d35 = this->length(3,5);
383 0 : const Real d24 = this->length(2,4);
384 0 : const Real max_diag = std::max(d06, std::max(d17, std::max(d35, d24)));
385 :
386 0 : libmesh_assert_not_equal_to ( max_diag, 0.0 );
387 :
388 : /**
389 : * Compute the minimum edge length.
390 : */
391 0 : std::vector<Real> edges(12);
392 0 : edges[0] = this->length(0,1);
393 0 : edges[1] = this->length(1,2);
394 0 : edges[2] = this->length(2,3);
395 0 : edges[3] = this->length(0,3);
396 0 : edges[4] = this->length(4,5);
397 0 : edges[5] = this->length(5,6);
398 0 : edges[6] = this->length(6,7);
399 0 : edges[7] = this->length(4,7);
400 0 : edges[8] = this->length(0,4);
401 0 : edges[9] = this->length(1,5);
402 0 : edges[10] = this->length(2,6);
403 0 : edges[11] = this->length(3,7);
404 :
405 0 : const Real min_edge = *(std::min_element(edges.begin(), edges.end()));
406 0 : return sqrt3 * min_edge / max_diag ;
407 : }
408 :
409 :
410 0 : case SHAPE:
411 : case SKEW:
412 : {
413 : // From: P. Knupp, "Algebraic mesh quality metrics for
414 : // unstructured initial meshes," Finite Elements in Analysis
415 : // and Design 39, 2003, p. 217-241, Sections 6.2 and 6.3.
416 :
417 : // Make local copies of points, we will access these several
418 : // times below.
419 : const Point
420 0 : x0 = point(0), x1 = point(1), x2 = point(2), x3 = point(3),
421 0 : x4 = point(4), x5 = point(5), x6 = point(6), x7 = point(7);
422 :
423 : // The columns of the Jacobian matrices are:
424 : // \vec{x}_{\xi} = \vec{a1}*eta*zeta + \vec{b1}*eta + \vec{c1}*zeta + \vec{d1}
425 : // \vec{x}_{\eta} = \vec{a2}*xi*zeta + \vec{b2}*xi + \vec{c2}*zeta + \vec{d2}
426 : // \vec{x}_{\zeta} = \vec{a3}*xi*eta + \vec{b3}*xi + \vec{c3}*eta + \vec{d3}
427 : // where the ai, bi, ci, and di are constants defined below.
428 0 : const Point a1 = -x0 + x1 - x2 + x3 + x4 - x5 + x6 - x7;
429 0 : const Point b1 = x0 - x1 + x2 - x3 + x4 - x5 + x6 - x7;
430 0 : const Point c1 = x0 - x1 - x2 + x3 - x4 + x5 + x6 - x7;
431 0 : const Point d1 = -x0 + x1 + x2 - x3 - x4 + x5 + x6 - x7;
432 :
433 0 : const Point a2 = a1;
434 0 : const Point b2 = b1;
435 0 : const Point c2 = x0 + x1 - x2 - x3 - x4 - x5 + x6 + x7;
436 0 : const Point d2 = -x0 - x1 + x2 + x3 - x4 - x5 + x6 + x7;
437 :
438 0 : const Point a3 = a1;
439 0 : const Point b3 = c1;
440 0 : const Point c3 = c2;
441 0 : const Point d3 = -x0 - x1 - x2 - x3 + x4 + x5 + x6 + x7;
442 :
443 : // Form the nodal Jacobians. These were computed using a
444 : // Python script and the formulas above. Note that we are
445 : // actually computing the Jacobian _columns_ and passing them
446 : // to the RealTensor constructor which expects _rows_, but
447 : // it's OK because we are only interested in determinants and
448 : // products which are not affected by taking the transpose.
449 : std::array<RealTensor, 8> A =
450 : {{
451 : RealTensor(d1, d2, d3),
452 0 : RealTensor(d1, b2 + d2, b3 + d3),
453 0 : RealTensor(b1 + d1, b2 + d2, a3 + b3 + c3 + d3),
454 0 : RealTensor(b1 + d1, d2, c3 + d3),
455 0 : RealTensor(c1 + d1, c2 + d2, d3),
456 0 : RealTensor(c1 + d1, a2 + b2 + c2 + d2, b3 + d3),
457 0 : RealTensor(a1 + b1 + c1 + d1, a2 + b2 + c2 + d2, a3 + b3 + c3 + d3),
458 0 : RealTensor(a1 + b1 + c1 + d1, c2 + d2, c3 + d3)
459 0 : }};
460 :
461 : // Compute Nodal areas, alpha_k = det(A_k).
462 : // If any of these are zero or negative, we return zero
463 : // (lowest possible value) for the quality, since the formulas
464 : // below require positive nodal areas.
465 : std::array<Real, 8> alpha;
466 0 : for (unsigned int k=0; k<alpha.size(); ++k)
467 : {
468 0 : alpha[k] = A[k].det();
469 0 : if (alpha[k] <= 0.)
470 0 : return 0.;
471 : }
472 :
473 : // Compute metric tensors, T_k = A_k^T * A_k.
474 0 : std::array<RealTensor, 8> T;
475 0 : for (unsigned int k=0; k<T.size(); ++k)
476 0 : T[k] = A[k] * A[k].transpose();
477 :
478 : // Compute and return the shape metric. These only use the
479 : // diagonal entries of the T_k.
480 0 : Real den = 0.;
481 0 : if (q == SHAPE)
482 : {
483 0 : for (unsigned int k=0; k<T.size(); ++k)
484 0 : den += T[k].tr() / std::pow(alpha[k], 2./3.);
485 0 : return (den == 0.) ? 0 : (24. / den);
486 : }
487 : else
488 : {
489 0 : for (unsigned int k=0; k<T.size(); ++k)
490 0 : den += std::pow(std::sqrt(T[k](0,0) * T[k](1,1) * T[k](2,2)) / alpha[k], 2./3.);
491 0 : return (den == 0.) ? 0 : (8. / den);
492 : }
493 : }
494 : #endif // LIBMESH_DIM >= 3
495 :
496 : /**
497 : * I don't know what to do for this metric.
498 : * Maybe the base class knows...
499 : */
500 2016 : default:
501 2016 : return Elem::quality(q);
502 : }
503 : }
504 :
505 :
506 :
507 0 : std::pair<Real, Real> Hex::qual_bounds (const ElemQuality q) const
508 : {
509 0 : std::pair<Real, Real> bounds;
510 :
511 0 : switch (q)
512 : {
513 :
514 0 : case ASPECT_RATIO:
515 0 : bounds.first = 1.;
516 0 : bounds.second = 4.;
517 0 : break;
518 :
519 0 : case SKEW:
520 0 : bounds.first = 0.;
521 0 : bounds.second = 0.5;
522 0 : break;
523 :
524 0 : case SHEAR:
525 : case SHAPE:
526 0 : bounds.first = 0.3;
527 0 : bounds.second = 1.;
528 0 : break;
529 :
530 0 : case CONDITION:
531 0 : bounds.first = 1.;
532 0 : bounds.second = 8.;
533 0 : break;
534 :
535 0 : case SCALED_JACOBIAN:
536 : case JACOBIAN:
537 0 : bounds.first = 0.5;
538 0 : bounds.second = 1.;
539 0 : break;
540 :
541 0 : case DISTORTION:
542 0 : bounds.first = 0.6;
543 0 : bounds.second = 1.;
544 0 : break;
545 :
546 0 : case TAPER:
547 0 : bounds.first = 0.;
548 0 : bounds.second = 0.4;
549 0 : break;
550 :
551 0 : case STRETCH:
552 0 : bounds.first = 0.25;
553 0 : bounds.second = 1.;
554 0 : break;
555 :
556 0 : case DIAGONAL:
557 0 : bounds.first = 0.65;
558 0 : bounds.second = 1.;
559 0 : break;
560 :
561 0 : case SIZE:
562 0 : bounds.first = 0.5;
563 0 : bounds.second = 1.;
564 0 : break;
565 :
566 0 : default:
567 0 : libMesh::out << "Warning: Invalid quality measure chosen." << std::endl;
568 0 : bounds.first = -1;
569 0 : bounds.second = -1;
570 : }
571 :
572 0 : return bounds;
573 : }
574 :
575 :
576 :
577 13777552 : bool Hex::on_reference_element(const Point & p,
578 : const Real eps) const
579 : {
580 1473852 : const Real & xi = p(0);
581 1473852 : const Real & eta = p(1);
582 1473852 : const Real & zeta = p(2);
583 :
584 : // The reference hexahedral element is [-1,1]^3.
585 25671034 : return ((xi >= -1.-eps) &&
586 11893482 : (xi <= 1.+eps) &&
587 9157600 : (eta >= -1.-eps) &&
588 7932478 : (eta <= 1.+eps) &&
589 23068562 : (zeta >= -1.-eps) &&
590 14970844 : (zeta <= 1.+eps));
591 : }
592 :
593 :
594 :
595 : const unsigned short int Hex::_second_order_vertex_child_number[27] =
596 : {
597 : 99,99,99,99,99,99,99,99, // Vertices
598 : 0,1,2,0,0,1,2,3,4,5,6,5, // Edges
599 : 0,0,1,2,0,4, // Faces
600 : 0 // Interior
601 : };
602 :
603 :
604 :
605 : const unsigned short int Hex::_second_order_vertex_child_index[27] =
606 : {
607 : 99,99,99,99,99,99,99,99, // Vertices
608 : 1,2,3,3,4,5,6,7,5,6,7,7, // Edges
609 : 2,5,6,7,7,6, // Faces
610 : 6 // Interior
611 : };
612 :
613 :
614 : const unsigned short int Hex::_second_order_adjacent_vertices[12][2] =
615 : {
616 : { 0, 1}, // vertices adjacent to node 8
617 : { 1, 2}, // vertices adjacent to node 9
618 : { 2, 3}, // vertices adjacent to node 10
619 : { 0, 3}, // vertices adjacent to node 11
620 :
621 : { 0, 4}, // vertices adjacent to node 12
622 : { 1, 5}, // vertices adjacent to node 13
623 : { 2, 6}, // vertices adjacent to node 14
624 : { 3, 7}, // vertices adjacent to node 15
625 :
626 : { 4, 5}, // vertices adjacent to node 16
627 : { 5, 6}, // vertices adjacent to node 17
628 : { 6, 7}, // vertices adjacent to node 18
629 : { 4, 7} // vertices adjacent to node 19
630 : };
631 :
632 :
633 : #ifdef LIBMESH_ENABLE_AMR
634 :
635 : // We number 125 "possible node locations" for a 2x2x2 refinement of
636 : // hexes with up to 3x3x3 nodes each
637 : const int Hex::_child_node_lookup[8][27] =
638 : {
639 : // node lookup for child 0 (near node 0)
640 : { 0, 2, 12, 10, 50, 52, 62, 60, 1, 7, 11, 5, 25, 27, 37, 35,
641 : 51, 57, 61, 55, 6, 26, 32, 36, 30, 56, 31},
642 :
643 : // node lookup for child 1 (near node 1)
644 : { 2, 4, 14, 12, 52, 54, 64, 62, 3, 9, 13, 7, 27, 29, 39, 37,
645 : 53, 59, 63, 57, 8, 28, 34, 38, 32, 58, 33},
646 :
647 : // node lookup for child 2 (near node 3)
648 : { 10, 12, 22, 20, 60, 62, 72, 70, 11, 17, 21, 15, 35, 37, 47, 45,
649 : 61, 67, 71, 65, 16, 36, 42, 46, 40, 66, 41},
650 :
651 : // node lookup for child 3 (near node 2)
652 : { 12, 14, 24, 22, 62, 64, 74, 72, 13, 19, 23, 17, 37, 39, 49, 47,
653 : 63, 69, 73, 67, 18, 38, 44, 48, 42, 68, 43},
654 :
655 : // node lookup for child 4 (near node 4)
656 : { 50, 52, 62, 60, 100, 102, 112, 110, 51, 57, 61, 55, 75, 77, 87, 85,
657 : 101, 107, 111, 105, 56, 76, 82, 86, 80, 106, 81},
658 :
659 : // node lookup for child 5 (near node 5)
660 : { 52, 54, 64, 62, 102, 104, 114, 112, 53, 59, 63, 57, 77, 79, 89, 87,
661 : 103, 109, 113, 107, 58, 78, 84, 88, 82, 108, 93},
662 :
663 : // node lookup for child 6 (near node 7)
664 : { 60, 62, 72, 70, 110, 112, 122, 120, 61, 67, 71, 65, 85, 87, 97, 95,
665 : 111, 117, 121, 115, 66, 86, 92, 96, 90, 116, 91},
666 :
667 : // node lookup for child 7 (near node 6)
668 : { 62, 64, 74, 72, 112, 114, 124, 122, 63, 69, 73, 67, 87, 89, 99, 97,
669 : 113, 119, 123, 117, 68, 88, 94, 98, 92, 118, 103}
670 : };
671 :
672 : #endif // LIBMESH_ENABLE_AMR
673 :
674 :
675 : } // namespace libMesh
|
