face_quad.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 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 
24 namespace libMesh
25 {
26 
27 
28 // ------------------------------------------------------------
29 // Quad class member functions
30 dof_id_type Quad::key (const unsigned int s) const
31 {
32  libmesh_assert_less (s, this->n_sides());
33 
34  switch (s)
35  {
36  case 0:
37  return
38  this->compute_key (this->node(0),
39  this->node(1));
40 
41  case 1:
42  return
43  this->compute_key (this->node(1),
44  this->node(2));
45 
46  case 2:
47  return
48  this->compute_key (this->node(2),
49  this->node(3));
50 
51  case 3:
52  return
53  this->compute_key (this->node(3),
54  this->node(0));
55  }
56 
57 
58  // We will never get here... Look at the code above.
59  libmesh_error();
60  return 0;
61 }
62 
63 
64 
65 AutoPtr<Elem> Quad::side (const unsigned int i) const
66 {
67  libmesh_assert_less (i, this->n_sides());
68 
69  Elem* edge = new Edge2;
70 
71  switch (i)
72  {
73  case 0:
74  {
75  edge->set_node(0) = this->get_node(0);
76  edge->set_node(1) = this->get_node(1);
77 
78  AutoPtr<Elem> ap_edge(edge);
79  return ap_edge;
80  }
81  case 1:
82  {
83  edge->set_node(0) = this->get_node(1);
84  edge->set_node(1) = this->get_node(2);
85 
86  AutoPtr<Elem> ap_edge(edge);
87  return ap_edge;
88  }
89  case 2:
90  {
91  edge->set_node(0) = this->get_node(2);
92  edge->set_node(1) = this->get_node(3);
93 
94  AutoPtr<Elem> ap_edge(edge);
95  return ap_edge;
96  }
97  case 3:
98  {
99  edge->set_node(0) = this->get_node(3);
100  edge->set_node(1) = this->get_node(0);
101 
102  AutoPtr<Elem> ap_edge(edge);
103  return ap_edge;
104  }
105  default:
106  {
107  libmesh_error();
108  }
109  }
110 
111 
112  // We will never get here... Look at the code above.
113  libmesh_error();
114  AutoPtr<Elem> ap_edge(edge);
115  return ap_edge;
116 }
117 
118 
119 
120 bool Quad::is_child_on_side(const unsigned int c,
121  const unsigned int s) const
122 {
123  libmesh_assert_less (c, this->n_children());
124  libmesh_assert_less (s, this->n_sides());
125 
126  // A quad's children and nodes don't share the same ordering:
127  // child 2 and 3 are swapped;
128  unsigned int n = (c < 2) ? c : 5-c;
129  return (n == s || n == (s+1)%4);
130 }
131 
132 
133 
134 unsigned int Quad::opposite_side(const unsigned int side_in) const
135 {
136  libmesh_assert_less (side_in, 4);
137 
138  return (side_in + 2) % 4;
139 }
140 
141 
142 
143 unsigned int Quad::opposite_node(const unsigned int node_in,
144  const unsigned int side_in) const
145 {
146  libmesh_assert_less (node_in, 8);
147  libmesh_assert_less (node_in, this->n_nodes());
148  libmesh_assert_less (side_in, this->n_sides());
149  libmesh_assert(this->is_node_on_side(node_in, side_in));
150 
151  //unsigned int opposite;
152 
153  static const unsigned char side02_nodes_map[] =
154  {3, 2, 1, 0, 6, 255, 4, 255};
155  static const unsigned char side13_nodes_map[] =
156  {1, 0, 3, 2, 255, 7, 255, 5};
157 
158  switch (side_in)
159  {
160  case 0:
161  case 2:
162  return side02_nodes_map[node_in];
163  break;
164  case 1:
165  case 3:
166  return side13_nodes_map[node_in];
167  break;
168  }
169 
170  libmesh_error();
171  return 255;
172 }
173 
174 
176 {
177  switch (q)
178  {
179 
184  case DISTORTION:
185  case DIAGONAL:
186  case STRETCH:
187  {
188  // Diagonal between node 0 and node 2
189  const Real d02 = this->length(0,2);
190 
191  // Diagonal between node 1 and node 3
192  const Real d13 = this->length(1,3);
193 
194  // Find the biggest and smallest diagonals
195  if ( (d02 > 0.) && (d13 >0.) )
196  if (d02 < d13) return d02 / d13;
197  else return d13 / d02;
198  else
199  return 0.;
200  break;
201  }
202 
203  default:
204  return Elem::quality(q);
205  }
206 
212  return Elem::quality(q);
213 }
214 
215 
216 
217 
218 
219 
220 std::pair<Real, Real> Quad::qual_bounds (const ElemQuality q) const
221 {
222  std::pair<Real, Real> bounds;
223 
224  switch (q)
225  {
226 
227  case ASPECT_RATIO:
228  bounds.first = 1.;
229  bounds.second = 4.;
230  break;
231 
232  case SKEW:
233  bounds.first = 0.;
234  bounds.second = 0.5;
235  break;
236 
237  case TAPER:
238  bounds.first = 0.;
239  bounds.second = 0.7;
240  break;
241 
242  case WARP:
243  bounds.first = 0.9;
244  bounds.second = 1.;
245  break;
246 
247  case STRETCH:
248  bounds.first = 0.25;
249  bounds.second = 1.;
250  break;
251 
252  case MIN_ANGLE:
253  bounds.first = 45.;
254  bounds.second = 90.;
255  break;
256 
257  case MAX_ANGLE:
258  bounds.first = 90.;
259  bounds.second = 135.;
260  break;
261 
262  case CONDITION:
263  bounds.first = 1.;
264  bounds.second = 4.;
265  break;
266 
267  case JACOBIAN:
268  bounds.first = 0.5;
269  bounds.second = 1.;
270  break;
271 
272  case SHEAR:
273  case SHAPE:
274  case SIZE:
275  bounds.first = 0.3;
276  bounds.second = 1.;
277  break;
278 
279  case DISTORTION:
280  bounds.first = 0.6;
281  bounds.second = 1.;
282  break;
283 
284  default:
285  libMesh::out << "Warning: Invalid quality measure chosen." << std::endl;
286  bounds.first = -1;
287  bounds.second = -1;
288  }
289 
290  return bounds;
291 }
292 
293 
294 
295 
296 const unsigned short int Quad::_second_order_adjacent_vertices[4][2] =
297 {
298  {0, 1}, // vertices adjacent to node 4
299  {1, 2}, // vertices adjacent to node 5
300  {2, 3}, // vertices adjacent to node 6
301  {0, 3} // vertices adjacent to node 7
302 };
303 
304 
305 
306 const unsigned short int Quad::_second_order_vertex_child_number[9] =
307 {
308  99,99,99,99, // Vertices
309  0,1,2,0, // Edges
310  0 // Interior
311 };
312 
313 
314 
315 const unsigned short int Quad::_second_order_vertex_child_index[9] =
316 {
317  99,99,99,99, // Vertices
318  1,2,3,3, // Edges
319  2 // Interior
320 };
321 
322 } // namespace libMesh

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:04 UTC

Hosted By:
SourceForge.net Logo