libMesh
face_quad.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 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 // C++ includes
19 
20 // Local includes
21 #include "libmesh/face_quad.h"
22 #include "libmesh/edge_edge2.h"
23 #include "libmesh/face_quad4.h"
24 
25 
26 namespace libMesh
27 {
28 
29 
30 
31 // ------------------------------------------------------------
32 // Quad class static member initializations
33 
34 
35 // We need to require C++11...
36 const Real Quad::_master_points[9][3] =
37  {
38  {-1, -1},
39  {1, -1},
40  {1, 1},
41  {-1, 1},
42  {0, -1},
43  {1, 0},
44  {0, 1},
45  {-1, 0},
46  {0, 0}
47  };
48 
49 
50 
51 // ------------------------------------------------------------
52 // Quad class member functions
53 dof_id_type Quad::key (const unsigned int s) const
54 {
55  libmesh_assert_less (s, this->n_sides());
56 
57  return this->compute_key(this->node_id(Quad4::side_nodes_map[s][0]),
58  this->node_id(Quad4::side_nodes_map[s][1]));
59 }
60 
61 
62 
63 unsigned int Quad::which_node_am_i(unsigned int side,
64  unsigned int side_node) const
65 {
66  libmesh_assert_less (side, this->n_sides());
67  libmesh_assert_less (side_node, 2);
68 
69  return Quad4::side_nodes_map[side][side_node];
70 }
71 
72 
73 
75 {
76  return this->compute_key(this->node_id(0),
77  this->node_id(1),
78  this->node_id(2),
79  this->node_id(3));
80 }
81 
82 
83 
84 UniquePtr<Elem> Quad::side_ptr (const unsigned int i)
85 {
86  libmesh_assert_less (i, this->n_sides());
87 
88  Elem * edge = new Edge2;
89 
90  for (unsigned n=0; n<edge->n_nodes(); ++n)
91  edge->set_node(n) = this->node_ptr(Quad4::side_nodes_map[i][n]);
92 
93  return UniquePtr<Elem>(edge);
94 }
95 
96 
97 
98 bool Quad::is_child_on_side(const unsigned int c,
99  const unsigned int s) const
100 {
101  libmesh_assert_less (c, this->n_children());
102  libmesh_assert_less (s, this->n_sides());
103 
104  // A quad's children and nodes don't share the same ordering:
105  // child 2 and 3 are swapped;
106  unsigned int n = (c < 2) ? c : 5-c;
107  return (n == s || n == (s+1)%4);
108 }
109 
110 
111 
112 unsigned int Quad::opposite_side(const unsigned int side_in) const
113 {
114  libmesh_assert_less (side_in, 4);
115 
116  return (side_in + 2) % 4;
117 }
118 
119 
120 
121 unsigned int Quad::opposite_node(const unsigned int node_in,
122  const unsigned int side_in) const
123 {
124  libmesh_assert_less (node_in, 8);
125  libmesh_assert_less (node_in, this->n_nodes());
126  libmesh_assert_less (side_in, this->n_sides());
127  libmesh_assert(this->is_node_on_side(node_in, side_in));
128 
129  static const unsigned char side02_nodes_map[] =
130  {3, 2, 1, 0, 6, 255, 4, 255};
131  static const unsigned char side13_nodes_map[] =
132  {1, 0, 3, 2, 255, 7, 255, 5};
133 
134  switch (side_in)
135  {
136  case 0:
137  case 2:
138  return side02_nodes_map[node_in];
139  case 1:
140  case 3:
141  return side13_nodes_map[node_in];
142  default:
143  libmesh_error_msg("Unsupported side_in = " << side_in);
144  }
145 
146  libmesh_error_msg("We'll never get here!");
147  return 255;
148 }
149 
150 
152 {
153  switch (q)
154  {
155 
160  case DISTORTION:
161  case DIAGONAL:
162  case STRETCH:
163  {
164  // Diagonal between node 0 and node 2
165  const Real d02 = this->length(0,2);
166 
167  // Diagonal between node 1 and node 3
168  const Real d13 = this->length(1,3);
169 
170  // Find the biggest and smallest diagonals
171  if ((d02 > 0.) && (d13 >0.))
172  if (d02 < d13) return d02 / d13;
173  else return d13 / d02;
174  else
175  return 0.;
176  break;
177  }
178 
179  default:
180  return Elem::quality(q);
181  }
182 
188  return Elem::quality(q);
189 }
190 
191 
192 
193 
194 
195 
196 std::pair<Real, Real> Quad::qual_bounds (const ElemQuality q) const
197 {
198  std::pair<Real, Real> bounds;
199 
200  switch (q)
201  {
202 
203  case ASPECT_RATIO:
204  bounds.first = 1.;
205  bounds.second = 4.;
206  break;
207 
208  case SKEW:
209  bounds.first = 0.;
210  bounds.second = 0.5;
211  break;
212 
213  case TAPER:
214  bounds.first = 0.;
215  bounds.second = 0.7;
216  break;
217 
218  case WARP:
219  bounds.first = 0.9;
220  bounds.second = 1.;
221  break;
222 
223  case STRETCH:
224  bounds.first = 0.25;
225  bounds.second = 1.;
226  break;
227 
228  case MIN_ANGLE:
229  bounds.first = 45.;
230  bounds.second = 90.;
231  break;
232 
233  case MAX_ANGLE:
234  bounds.first = 90.;
235  bounds.second = 135.;
236  break;
237 
238  case CONDITION:
239  bounds.first = 1.;
240  bounds.second = 4.;
241  break;
242 
243  case JACOBIAN:
244  bounds.first = 0.5;
245  bounds.second = 1.;
246  break;
247 
248  case SHEAR:
249  case SHAPE:
250  case SIZE:
251  bounds.first = 0.3;
252  bounds.second = 1.;
253  break;
254 
255  case DISTORTION:
256  bounds.first = 0.6;
257  bounds.second = 1.;
258  break;
259 
260  default:
261  libMesh::out << "Warning: Invalid quality measure chosen." << std::endl;
262  bounds.first = -1;
263  bounds.second = -1;
264  }
265 
266  return bounds;
267 }
268 
269 
270 
271 
272 const unsigned short int Quad::_second_order_adjacent_vertices[4][2] =
273  {
274  {0, 1}, // vertices adjacent to node 4
275  {1, 2}, // vertices adjacent to node 5
276  {2, 3}, // vertices adjacent to node 6
277  {0, 3} // vertices adjacent to node 7
278  };
279 
280 
281 
282 const unsigned short int Quad::_second_order_vertex_child_number[9] =
283  {
284  99,99,99,99, // Vertices
285  0,1,2,0, // Edges
286  0 // Interior
287  };
288 
289 
290 
291 const unsigned short int Quad::_second_order_vertex_child_index[9] =
292  {
293  99,99,99,99, // Vertices
294  1,2,3,3, // Edges
295  2 // Interior
296  };
297 
298 
299 
300 #ifdef LIBMESH_ENABLE_AMR
301 
302 // We number 25 "possible node locations" for a 2x2 refinement of
303 // quads with up to 3x3 nodes each
304 const int Quad::_child_node_lookup[4][9] =
305  {
306  // node lookup for child 0 (near node 0)
307  { 0, 2, 12, 10, 1, 7, 11, 5, 6},
308 
309  // node lookup for child 1 (near node 1)
310  { 2, 4, 14, 12, 3, 9, 13, 7, 8},
311 
312  // node lookup for child 2 (near node 3)
313  { 10, 12, 22, 20, 11, 17, 21, 15, 16},
314 
315  // node lookup for child 3 (near node 2)
316  { 12, 14, 24, 22, 13, 19, 23, 17, 18}
317  };
318 
319 
320 #endif
321 
322 } // namespace libMesh
virtual bool is_node_on_side(const unsigned int n, const unsigned int s) const =0
virtual unsigned int opposite_node(const unsigned int n, const unsigned int s) const libmesh_override
Definition: face_quad.C:121
virtual Node *& set_node(const unsigned int i)
Definition: elem.h:1941
virtual bool is_child_on_side(const unsigned int c, const unsigned int s) const libmesh_override
Definition: face_quad.C:98
static const int _child_node_lookup[4][9]
Lookup table from child id, child node id to "possible node location" (a simple dictionary-index in a...
Definition: face_quad.h:201
virtual dof_id_type key() const libmesh_override
Definition: face_quad.C:74
virtual Real quality(const ElemQuality q) const libmesh_override
Definition: face_quad.C:151
UniquePtr< Elem > side(const unsigned int i) const
Definition: elem.h:2105
static const unsigned int side_nodes_map[4][2]
This maps the node of the side to element node numbers.
Definition: face_quad4.h:122
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
static const Real _master_points[9][3]
Master element node locations.
Definition: face_quad.h:195
The libMesh namespace provides an interface to certain functionality in the library.
virtual UniquePtr< Elem > side_ptr(const unsigned int i) libmesh_override
Definition: face_quad.C:84
libmesh_assert(j)
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual unsigned int n_nodes() const =0
const Node * node_ptr(const unsigned int i) const
Definition: elem.h:1874
virtual unsigned int opposite_side(const unsigned int s) const libmesh_override
Definition: face_quad.C:112
virtual Real quality(const ElemQuality q) const
Definition: elem.C:1547
static const unsigned short int _second_order_vertex_child_index[9]
Vector that names the child vertex index for each second order node.
Definition: face_quad.h:190
static const unsigned short int _second_order_vertex_child_number[9]
Vector that names a child sharing each second order node.
Definition: face_quad.h:185
static const unsigned short int _second_order_adjacent_vertices[4][2]
Matrix that tells which vertices define the location of mid-side (or second-order) nodes...
Definition: face_quad.h:180
virtual unsigned int which_node_am_i(unsigned int side, unsigned int side_node) const libmesh_override
Definition: face_quad.C:63
virtual unsigned int n_nodes() const libmesh_override
Definition: face_quad.h:80
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual unsigned int n_children() const libmesh_override
Definition: face_quad.h:100
OStreamProxy out
virtual unsigned int n_sides() const libmesh_override
Definition: face_quad.h:85
The Edge2 is an element in 1D composed of 2 nodes.
Definition: edge_edge2.h:43
Real length(const unsigned int n1, const unsigned int n2) const
Definition: elem.C:492
dof_id_type node_id(const unsigned int i) const
Definition: elem.h:1831
ElemQuality
Defines an enum for element quality metrics.
static dof_id_type compute_key(dof_id_type n0)
Definition: elem.h:2621
virtual std::pair< Real, Real > qual_bounds(const ElemQuality q) const libmesh_override
Definition: face_quad.C:196
uint8_t dof_id_type
Definition: id_types.h:64