libMesh
cell_inf_prism12.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 // Local includes
19 #include "libmesh/libmesh_config.h"
20 
21 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
22 
23 // C++ includes
24 
25 // Local includes cont'd
26 #include "libmesh/cell_inf_prism12.h"
27 #include "libmesh/edge_edge3.h"
28 #include "libmesh/edge_inf_edge2.h"
29 #include "libmesh/face_tri6.h"
30 #include "libmesh/face_inf_quad6.h"
31 #include "libmesh/side.h"
32 
33 namespace libMesh
34 {
35 
36 
37 // ------------------------------------------------------------
38 // InfPrism12 class static member initializations
39 const unsigned int InfPrism12::side_nodes_map[4][6] =
40  {
41  { 0, 1, 2, 6, 7, 8}, // Side 0
42  { 0, 1, 3, 4, 6, 9}, // Side 1
43  { 1, 2, 4, 5, 7, 10}, // Side 2
44  { 2, 0, 5, 3, 8, 11} // Side 3
45  };
46 
47 const unsigned int InfPrism12::edge_nodes_map[6][3] =
48  {
49  {0, 1, 6}, // Side 0
50  {1, 2, 7}, // Side 1
51  {0, 2, 8}, // Side 2
52  {0, 3, 99}, // Side 3
53  {1, 4, 99}, // Side 4
54  {2, 5, 99} // Side 5
55  };
56 
57 
58 // ------------------------------------------------------------
59 // InfPrism12 class member functions
60 
61 bool InfPrism12::is_vertex(const unsigned int i) const
62 {
63  if (i < 3)
64  return true;
65  return false;
66 }
67 
68 bool InfPrism12::is_edge(const unsigned int i) const
69 {
70  if (i < 3)
71  return false;
72  if (i > 8)
73  return false;
74  return true;
75 }
76 
77 bool InfPrism12::is_face(const unsigned int i) const
78 {
79  if (i > 8)
80  return true;
81  return false;
82 }
83 
84 bool InfPrism12::is_node_on_side(const unsigned int n,
85  const unsigned int s) const
86 {
87  libmesh_assert_less (s, n_sides());
88  for (unsigned int i = 0; i != 6; ++i)
89  if (side_nodes_map[s][i] == n)
90  return true;
91  return false;
92 }
93 
94 bool InfPrism12::is_node_on_edge(const unsigned int n,
95  const unsigned int e) const
96 {
97  libmesh_assert_less (e, n_edges());
98  for (unsigned int i = 0; i != 3; ++i)
99  if (edge_nodes_map[e][i] == n)
100  return true;
101  return false;
102 }
103 
104 
105 
106 unsigned int InfPrism12::which_node_am_i(unsigned int side,
107  unsigned int side_node) const
108 {
109  libmesh_assert_less (side, this->n_sides());
110 
111  // Never more than 6 nodes per side.
112  libmesh_assert_less(side_node, 6);
113 
114  return InfPrism12::side_nodes_map[side][side_node];
115 }
116 
117 
118 
119 UniquePtr<Elem> InfPrism12::build_side_ptr (const unsigned int i,
120  bool proxy)
121 {
122  libmesh_assert_less (i, this->n_sides());
123 
124  if (proxy)
125  {
126  switch (i)
127  {
128  // base
129  case 0:
130  return UniquePtr<Elem>(new Side<Tri6,InfPrism12>(this,i));
131 
132  // ifem sides
133  case 1:
134  case 2:
135  case 3:
136  return UniquePtr<Elem>(new Side<InfQuad6,InfPrism12>(this,i));
137 
138  default:
139  libmesh_error_msg("Invalid side i = " << i);
140  }
141  }
142 
143  else
144  {
145  // Create NULL pointer to be initialized, returned later.
146  Elem * face = libmesh_nullptr;
147 
148  switch (i)
149  {
150  case 0: // the triangular face at z=-1, base face
151  {
152  face = new Tri6;
153  break;
154  }
155 
156  case 1: // the quad face at y=0
157  case 2: // the other quad face
158  case 3: // the quad face at x=0
159  {
160  face = new InfQuad6;
161  break;
162  }
163 
164  default:
165  libmesh_error_msg("Invalid side i = " << i);
166  }
167 
168  face->subdomain_id() = this->subdomain_id();
169 
170  // Set the nodes
171  for (unsigned n=0; n<face->n_nodes(); ++n)
172  face->set_node(n) = this->node_ptr(InfPrism12::side_nodes_map[i][n]);
173 
174  return UniquePtr<Elem>(face);
175  }
176 
177  libmesh_error_msg("We'll never get here!");
178  return UniquePtr<Elem>();
179 }
180 
181 
182 UniquePtr<Elem> InfPrism12::build_edge_ptr (const unsigned int i)
183 {
184  libmesh_assert_less (i, this->n_edges());
185 
186  if (i < 3) // base edges
187  return UniquePtr<Elem>(new SideEdge<Edge3,InfPrism12>(this,i));
188  // infinite edges
189  return UniquePtr<Elem>(new SideEdge<InfEdge2,InfPrism12>(this,i));
190 }
191 
192 
193 void InfPrism12::connectivity(const unsigned int sc,
194  const IOPackage iop,
195  std::vector<dof_id_type> & conn) const
196 {
197  libmesh_assert(_nodes);
198  libmesh_assert_less (sc, this->n_sub_elem());
199  libmesh_assert_not_equal_to (iop, INVALID_IO_PACKAGE);
200 
201  switch (iop)
202  {
203  case TECPLOT:
204  {
205  conn.resize(8);
206  switch (sc)
207  {
208  case 0:
209 
210  // guess this is a collapsed hex8
211  conn[0] = this->node_id(0)+1;
212  conn[1] = this->node_id(6)+1;
213  conn[2] = this->node_id(8)+1;
214  conn[3] = this->node_id(8)+1;
215  conn[4] = this->node_id(3)+1;
216  conn[5] = this->node_id(9)+1;
217  conn[6] = this->node_id(11)+1;
218  conn[7] = this->node_id(11)+1;
219 
220  return;
221 
222  case 1:
223 
224  conn[0] = this->node_id(6)+1;
225  conn[1] = this->node_id(7)+1;
226  conn[2] = this->node_id(8)+1;
227  conn[3] = this->node_id(8)+1;
228  conn[4] = this->node_id(9)+1;
229  conn[5] = this->node_id(10)+1;
230  conn[6] = this->node_id(11)+1;
231  conn[7] = this->node_id(11)+1;
232 
233  return;
234 
235  case 2:
236 
237  conn[0] = this->node_id(6)+1;
238  conn[1] = this->node_id(1)+1;
239  conn[2] = this->node_id(7)+1;
240  conn[3] = this->node_id(7)+1;
241  conn[4] = this->node_id(9)+1;
242  conn[5] = this->node_id(4)+1;
243  conn[6] = this->node_id(10)+1;
244  conn[7] = this->node_id(10)+1;
245 
246  return;
247 
248  case 3:
249 
250  conn[0] = this->node_id(8)+1;
251  conn[1] = this->node_id(7)+1;
252  conn[2] = this->node_id(2)+1;
253  conn[3] = this->node_id(2)+1;
254  conn[4] = this->node_id(11)+1;
255  conn[5] = this->node_id(10)+1;
256  conn[6] = this->node_id(5)+1;
257  conn[7] = this->node_id(5)+1;
258 
259  return;
260 
261  default:
262  libmesh_error_msg("Invalid sc = " << sc);
263  }
264  }
265 
266  default:
267  libmesh_error_msg("Unsupported IO package " << iop);
268  }
269 }
270 
271 
272 
273 
274 
275 unsigned short int InfPrism12::second_order_adjacent_vertex (const unsigned int n,
276  const unsigned int v) const
277 {
278  libmesh_assert_greater_equal (n, this->n_vertices());
279  libmesh_assert_less (n, this->n_nodes());
280  libmesh_assert_less (v, 2);
281  return _second_order_adjacent_vertices[n-this->n_vertices()][v];
282 }
283 
284 
285 
286 const unsigned short int InfPrism12::_second_order_adjacent_vertices[6][2] =
287  {
288  { 0, 1}, // vertices adjacent to node 6
289  { 1, 2}, // vertices adjacent to node 7
290  { 0, 2}, // vertices adjacent to node 8
291 
292  { 3, 4}, // vertices adjacent to node 9
293  { 4, 5}, // vertices adjacent to node 10
294  { 3, 5} // vertices adjacent to node 11
295  };
296 
297 
298 
299 std::pair<unsigned short int, unsigned short int>
300 InfPrism12::second_order_child_vertex (const unsigned int n) const
301 {
302  libmesh_assert_greater_equal (n, this->n_vertices());
303  libmesh_assert_less (n, this->n_nodes());
304 
305  return std::pair<unsigned short int, unsigned short int>
306  (_second_order_vertex_child_number[n],
307  _second_order_vertex_child_index[n]);
308 }
309 
310 
311 
312 const unsigned short int InfPrism12::_second_order_vertex_child_number[12] =
313  {
314  99,99,99,99,99,99, // Vertices
315  0,1,0, // Edges
316  0,1,0 // Faces
317  };
318 
319 
320 
321 const unsigned short int InfPrism12::_second_order_vertex_child_index[12] =
322  {
323  99,99,99,99,99,99, // Vertices
324  1,2,2, // Edges
325  4,5,5 // Faces
326  };
327 
328 
329 
330 #ifdef LIBMESH_ENABLE_AMR
331 
332 const float InfPrism12::_embedding_matrix[4][12][12] =
333  {
334  // embedding matrix for child 0
335  {
336  // 0 1 2 3 4 5 6 7 8 9 10 11 th parent Node
337  { 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}, // 0th child N.
338  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0}, // 1
339  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0}, // 2
340  { 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}, // 3
341  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0}, // 4
342  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0}, // 5
343  { 0.375, -0.125, 0.0, 0.0, 0.0, 0.0, 0.75, 0.0, 0.0, 0.0, 0.0, 0.0}, // 6
344  { 0.0, -0.125, -0.125, 0.0, 0.0, 0.0, 0.5, 0.25, 0.5, 0.0, 0.0, 0.0}, // 7
345  { 0.375, 0.0, -0.125, 0.0, 0.0, 0.0, 0.0, 0.0, 0.75, 0.0, 0.0, 0.0}, // 8
346  { 0.0, 0.0, 0.0, 0.375, -0.125, 0.0, 0.0, 0.0, 0.0, 0.75, 0.0, 0.0}, // 9
347  { 0.0, 0.0, 0.0, 0.0, -0.125, -0.125, 0.0, 0.0, 0.0, 0.5, 0.25, 0.5}, // 10
348  { 0.0, 0.0, 0.0, 0.375, 0.0, -0.125, 0.0, 0.0, 0.0, 0.0, 0.0, 0.75} // 11
349  },
350 
351  // embedding matrix for child 1
352  {
353  // 0 1 2 3 4 5 6 7 8 9 10 11 th parent Node
354  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0}, // 0th child N.
355  { 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}, // 1
356  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0}, // 2
357  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0}, // 3
358  { 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}, // 4
359  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0}, // 5
360  { -0.125, 0.375, 0.0, 0.0, 0.0, 0.0, 0.75, 0.0, 0.0, 0.0, 0.0, 0.0}, // 6
361  { 0.0, 0.375, -0.125, 0.0, 0.0, 0.0, 0.0, 0.75, 0.0, 0.0, 0.0, 0.0}, // 7
362  { -0.125, 0.0, -0.125, 0.0, 0.0, 0.0, 0.5, 0.5, 0.25, 0.0, 0.0, 0.0}, // 8
363  { 0.0, 0.0, 0.0, -0.125, 0.375, 0.0, 0.0, 0.0, 0.0, 0.75, 0.0, 0.0}, // 9
364  { 0.0, 0.0, 0.0, 0.0, 0.375, -0.125, 0.0, 0.0, 0.0, 0.0, 0.75, 0.0}, // 10
365  { 0.0, 0.0, 0.0, -0.125, 0.0, -0.125, 0.0, 0.0, 0.0, 0.5, 0.5, 0.25} // 11
366  },
367 
368  // embedding matrix for child 2
369  {
370  // 0 1 2 3 4 5 6 7 8 9 10 11 th parent Node
371  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0}, // 0th child N.
372  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0}, // 1
373  { 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}, // 2
374  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0}, // 3
375  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0}, // 4
376  { 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}, // 5
377  { -0.125, -0.125, 0.0, 0.0, 0.0, 0.0, 0.25, 0.5, 0.5, 0.0, 0.0, 0.0}, // 6
378  { 0.0, -0.125, 0.375, 0.0, 0.0, 0.0, 0.0, 0.75, 0.0, 0.0, 0.0, 0.0}, // 7
379  { -0.125, 0.0, 0.375, 0.0, 0.0, 0.0, 0.0, 0.0, 0.75, 0.0, 0.0, 0.0}, // 8
380  { 0.0, 0.0, 0.0, -0.125, -0.125, 0.0, 0.0, 0.0, 0.0, 0.25, 0.5, 0.5}, // 9
381  { 0.0, 0.0, 0.0, 0.0, -0.125, 0.375, 0.0, 0.0, 0.0, 0.0, 0.75, 0.0}, // 10
382  { 0.0, 0.0, 0.0, -0.125, 0.0, 0.375, 0.0, 0.0, 0.0, 0.0, 0.0, 0.75} // 11
383  },
384 
385  // embedding matrix for child 3
386  {
387  // 0 1 2 3 4 5 6 7 8 9 10 11 th parent Node
388  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0}, // 0th child N.
389  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0}, // 1
390  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0}, // 2
391  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0}, // 3
392  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0}, // 4
393  { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0}, // 5
394  { -0.125, 0.0, -0.125, 0.0, 0.0, 0.0, 0.5, 0.5, 0.25, 0.0, 0.0, 0.0}, // 6
395  { -0.125, -0.125, 0.0, 0.0, 0.0, 0.0, 0.25, 0.5, 0.5, 0.0, 0.0, 0.0}, // 7
396  { 0.0, -0.125, -0.125, 0.0, 0.0, 0.0, 0.5, 0.25, 0.5, 0.0, 0.0, 0.0}, // 8
397  { 0.0, 0.0, 0.0, -0.125, 0.0, -0.125, 0.0, 0.0, 0.0, 0.5, 0.5, 0.25}, // 9
398  { 0.0, 0.0, 0.0, -0.125, -0.125, 0.0, 0.0, 0.0, 0.0, 0.25, 0.5, 0.5}, // 10
399  { 0.0, 0.0, 0.0, 0.0, -0.125, -0.125, 0.0, 0.0, 0.0, 0.5, 0.25, 0.5} // 11
400  }
401 
402  };
403 
404 
405 
406 
407 #endif
408 
409 } // namespace libMesh
410 
411 #endif // ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
unsigned short int side
Definition: xdr_io.C:49
const class libmesh_nullptr_t libmesh_nullptr
The libMesh namespace provides an interface to certain functionality in the library.
libmesh_assert(j)
const dof_id_type n_nodes
Definition: tecplot_io.C:67
IOPackage
libMesh interfaces with several different software packages for the purposes of creating, reading, and writing mesh files.