libMesh
cell_inf_prism6.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_prism6.h"
27 #include "libmesh/edge_edge2.h"
28 #include "libmesh/edge_inf_edge2.h"
29 #include "libmesh/fe_interface.h"
30 #include "libmesh/side.h"
31 #include "libmesh/face_inf_quad4.h"
32 #include "libmesh/face_tri3.h"
33 
34 namespace libMesh
35 {
36 
37 
38 // ------------------------------------------------------------
39 // InfPrism6 class static member initializations
40 const unsigned int InfPrism6::side_nodes_map[4][4] =
41  {
42  { 0, 1, 2, 99}, // Side 0
43  { 0, 1, 3, 4}, // Side 1
44  { 1, 2, 4, 5}, // Side 2
45  { 2, 0, 5, 3} // Side 3
46  };
47 
48 const unsigned int InfPrism6::edge_nodes_map[6][2] =
49  {
50  {0, 1}, // Edge 0
51  {1, 2}, // Edge 1
52  {0, 2}, // Edge 2
53  {0, 3}, // Edge 3
54  {1, 4}, // Edge 4
55  {2, 5} // Edge 5
56  };
57 
58 
59 // ------------------------------------------------------------
60 // InfPrism6 class member functions
61 
62 bool InfPrism6::is_vertex(const unsigned int i) const
63 {
64  if (i < 3)
65  return true;
66  return false;
67 }
68 
69 bool InfPrism6::is_edge(const unsigned int i) const
70 {
71  if (i < 3)
72  return false;
73  return true;
74 }
75 
76 bool InfPrism6::is_face(const unsigned int) const
77 {
78  return false;
79 }
80 
81 bool InfPrism6::is_node_on_side(const unsigned int n,
82  const unsigned int s) const
83 {
84  libmesh_assert_less (s, n_sides());
85  for (unsigned int i = 0; i != 4; ++i)
86  if (side_nodes_map[s][i] == n)
87  return true;
88  return false;
89 }
90 
91 bool InfPrism6::is_node_on_edge(const unsigned int n,
92  const unsigned int e) const
93 {
94  libmesh_assert_less (e, n_edges());
95  for (unsigned int i = 0; i != 2; ++i)
96  if (edge_nodes_map[e][i] == n)
97  return true;
98  return false;
99 }
100 
101 
102 UniquePtr<Elem> InfPrism6::build_side_ptr (const unsigned int i,
103  bool proxy)
104 {
105  libmesh_assert_less (i, this->n_sides());
106 
107  if (proxy)
108  {
109  switch (i)
110  {
111  // base
112  case 0:
113  return UniquePtr<Elem>(new Side<Tri3,InfPrism6>(this,i));
114 
115  // ifem sides
116  case 1:
117  case 2:
118  case 3:
119  return UniquePtr<Elem>(new Side<InfQuad4,InfPrism6>(this,i));
120 
121  default:
122  libmesh_error_msg("Invalid side i = " << i);
123  }
124  }
125 
126  else
127  {
128  // Create NULL pointer to be initialized, returned later.
129  Elem * face = libmesh_nullptr;
130 
131  switch (i)
132  {
133  case 0: // the triangular face at z=-1, base face
134  {
135  face = new Tri3;
136  break;
137  }
138 
139  case 1: // the quad face at y=0
140  case 2: // the other quad face
141  case 3: // the quad face at x=0
142  {
143  face = new InfQuad4;
144  break;
145  }
146 
147  default:
148  libmesh_error_msg("Invalid side i = " << i);
149  }
150 
151  face->subdomain_id() = this->subdomain_id();
152 
153  // Set the nodes
154  for (unsigned n=0; n<face->n_nodes(); ++n)
155  face->set_node(n) = this->node_ptr(InfPrism6::side_nodes_map[i][n]);
156 
157  return UniquePtr<Elem>(face);
158  }
159 
160  libmesh_error_msg("We'll never get here!");
161  return UniquePtr<Elem>();
162 }
163 
164 
165 UniquePtr<Elem> InfPrism6::build_edge_ptr (const unsigned int i)
166 {
167  libmesh_assert_less (i, n_edges());
168 
169  if (i < 3)
170  return UniquePtr<Elem>(new SideEdge<Edge2,InfPrism6>(this,i));
171  return UniquePtr<Elem>(new SideEdge<InfEdge2,InfPrism6>(this,i));
172 }
173 
174 void InfPrism6::connectivity(const unsigned int libmesh_dbg_var(sc),
175  const IOPackage iop,
176  std::vector<dof_id_type> & conn) const
177 {
178  libmesh_assert(_nodes);
179  libmesh_assert_less (sc, this->n_sub_elem());
180  libmesh_assert_not_equal_to (iop, INVALID_IO_PACKAGE);
181 
182  switch (iop)
183  {
184  case TECPLOT:
185  {
186  conn.resize(8);
187  conn[0] = this->node_id(0)+1;
188  conn[1] = this->node_id(1)+1;
189  conn[2] = this->node_id(2)+1;
190  conn[3] = this->node_id(2)+1;
191  conn[4] = this->node_id(3)+1;
192  conn[5] = this->node_id(4)+1;
193  conn[6] = this->node_id(5)+1;
194  conn[7] = this->node_id(5)+1;
195  return;
196  }
197 
198  default:
199  libmesh_error_msg("Unsupported IO package " << iop);
200  }
201 }
202 
203 
204 
205 
206 
207 #ifdef LIBMESH_ENABLE_AMR
208 
209 const float InfPrism6::_embedding_matrix[4][6][6] =
210  {
211  // embedding matrix for child 0
212  {
213  // 0 1 2 3 4 5 th parent Node
214  { 1.0, 0.0, 0.0, 0.0, 0.0, 0.0}, // 0th child N.
215  { 0.5, 0.5, 0.0, 0.0, 0.0, 0.0}, // 1
216  { 0.5, 0.0, 0.5, 0.0, 0.0, 0.0}, // 2
217  { 0.0, 0.0, 0.0, 1.0, 0.0, 0.0}, // 3
218  { 0.0, 0.0, 0.0, 0.5, 0.5, 0.0}, // 4
219  { 0.0, 0.0, 0.0, 0.5, 0.0, 0.5} // 5
220  },
221 
222  // embedding matrix for child 1
223  {
224  // 0 1 2 3 4 5 th parent Node
225  { 0.5, 0.5, 0.0, 0.0, 0.0, 0.0}, // 0th child N.
226  { 0.0, 1.0, 0.0, 0.0, 0.0, 0.0}, // 1
227  { 0.0, 0.5, 0.5, 0.0, 0.0, 0.0}, // 2
228  { 0.0, 0.0, 0.0, 0.5, 0.5, 0.0}, // 3
229  { 0.0, 0.0, 0.0, 0.0, 1.0, 0.0}, // 4
230  { 0.0, 0.0, 0.0, 0.0, 0.5, 0.5} // 5
231  },
232 
233  // embedding matrix for child 2
234  {
235  // 0 1 2 3 4 5 th parent Node
236  { 0.5, 0.0, 0.5, 0.0, 0.0, 0.0}, // 0th child N.
237  { 0.0, 0.5, 0.5, 0.0, 0.0, 0.0}, // 1
238  { 0.0, 0.0, 1.0, 0.0, 0.0, 0.0}, // 2
239  { 0.0, 0.0, 0.0, 0.5, 0.0, 0.5}, // 3
240  { 0.0, 0.0, 0.0, 0.0, 0.5, 0.5}, // 4
241  { 0.0, 0.0, 0.0, 0.0, 0.0, 1.0} // 5
242  },
243 
244  // embedding matrix for child 3
245  {
246  // 0 1 2 3 4 5 th parent Node
247  { 0.5, 0.5, 0.0, 0.0, 0.0, 0.0}, // 0th child N.
248  { 0.0, 0.5, 0.5, 0.0, 0.0, 0.0}, // 1
249  { 0.5, 0.0, 0.5, 0.0, 0.0, 0.0}, // 2
250  { 0.0, 0.0, 0.0, 0.5, 0.5, 0.0}, // 3
251  { 0.0, 0.0, 0.0, 0.0, 0.5, 0.5}, // 4
252  { 0.0, 0.0, 0.0, 0.5, 0.0, 0.5} // 5
253  }
254  };
255 
256 
257 
258 #endif
259 
260 } // namespace libMesh
261 
262 #endif // ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
const class libmesh_nullptr_t libmesh_nullptr
The libMesh namespace provides an interface to certain functionality in the library.
libmesh_assert(j)
PetscErrorCode Vec Mat libmesh_dbg_var(j)
IOPackage
libMesh interfaces with several different software packages for the purposes of creating, reading, and writing mesh files.