face_tri3.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/side.h"
22 #include "libmesh/edge_edge2.h"
23 #include "libmesh/face_tri3.h"
24 
25 namespace libMesh
26 {
27 
28 
29 
30 // ------------------------------------------------------------
31 // Tri3 class static member initializations
32 const unsigned int Tri3::side_nodes_map[3][2] =
33 {
34  {0, 1}, // Side 0
35  {1, 2}, // Side 1
36  {2, 0} // Side 2
37 };
38 
39 
40 #ifdef LIBMESH_ENABLE_AMR
41 
42 const float Tri3::_embedding_matrix[4][3][3] =
43 {
44  // embedding matrix for child 0
45  {
46  // 0 1 2
47  {1.0, 0.0, 0.0}, // 0
48  {0.5, 0.5, 0.0}, // 1
49  {0.5, 0.0, 0.5} // 2
50  },
51 
52  // embedding matrix for child 1
53  {
54  // 0 1 2
55  {0.5, 0.5, 0.0}, // 0
56  {0.0, 1.0, 0.0}, // 1
57  {0.0, 0.5, 0.5} // 2
58  },
59 
60  // embedding matrix for child 2
61  {
62  // 0 1 2
63  {0.5, 0.0, 0.5}, // 0
64  {0.0, 0.5, 0.5}, // 1
65  {0.0, 0.0, 1.0} // 2
66  },
67 
68  // embedding matrix for child 3
69  {
70  // 0 1 2
71  {0.5, 0.5, 0.0}, // 0
72  {0.0, 0.5, 0.5}, // 1
73  {0.5, 0.0, 0.5} // 2
74  }
75 };
76 
77 #endif
78 
79 
80 
81 // ------------------------------------------------------------
82 // Tri3 class member functions
83 
84 bool Tri3::is_vertex(const unsigned int) const
85 {
86  return true;
87 }
88 
89 bool Tri3::is_edge(const unsigned int) const
90 {
91  return false;
92 }
93 
94 bool Tri3::is_face(const unsigned int) const
95 {
96  return false;
97 }
98 
99 bool Tri3::is_node_on_side(const unsigned int n,
100  const unsigned int s) const
101 {
102  libmesh_assert_less (s, n_sides());
103  for (unsigned int i = 0; i != 2; ++i)
104  if (side_nodes_map[s][i] == n)
105  return true;
106  return false;
107 }
108 
109 AutoPtr<Elem> Tri3::build_side (const unsigned int i,
110  bool proxy) const
111 {
112  libmesh_assert_less (i, this->n_sides());
113 
114  if (proxy)
115  {
116  AutoPtr<Elem> ap(new Side<Edge2,Tri3>(this,i));
117  return ap;
118  }
119 
120  else
121  {
122  AutoPtr<Elem> edge(new Edge2);
123  edge->subdomain_id() = this->subdomain_id();
124 
125  switch (i)
126  {
127  case 0:
128  {
129  edge->set_node(0) = this->get_node(0);
130  edge->set_node(1) = this->get_node(1);
131 
132  return edge;
133  }
134  case 1:
135  {
136  edge->set_node(0) = this->get_node(1);
137  edge->set_node(1) = this->get_node(2);
138 
139  return edge;
140  }
141  case 2:
142  {
143  edge->set_node(0) = this->get_node(2);
144  edge->set_node(1) = this->get_node(0);
145 
146  return edge;
147  }
148  default:
149  {
150  libmesh_error();
151  }
152  }
153  }
154 
155  // We will never get here... Look at the code above.
156  libmesh_error();
157  AutoPtr<Elem> ap(NULL); return ap;
158 }
159 
160 
161 void Tri3::connectivity(const unsigned int libmesh_dbg_var(sf),
162  const IOPackage iop,
163  std::vector<dof_id_type>& conn) const
164 {
165  libmesh_assert_less (sf, this->n_sub_elem());
166  libmesh_assert_not_equal_to (iop, INVALID_IO_PACKAGE);
167 
168  switch (iop)
169  {
170  case TECPLOT:
171  {
172  conn.resize(4);
173  conn[0] = this->node(0)+1;
174  conn[1] = this->node(1)+1;
175  conn[2] = this->node(2)+1;
176  conn[3] = this->node(2)+1;
177  return;
178  }
179 
180  case VTK:
181  {
182  conn.resize(3);
183  conn[0] = this->node(0);
184  conn[1] = this->node(1);
185  conn[2] = this->node(2);
186  return;
187  }
188 
189  default:
190  libmesh_error();
191  }
192 
193  libmesh_error();
194 }
195 
196 
197 
198 
199 
200 
202 {
203  // 3-node triangles have the following formula for computing the area
204  Point v10 ( *(this->get_node(1)) - *(this->get_node(0)) );
205 
206  Point v20 ( *(this->get_node(2)) - *(this->get_node(0)) );
207 
208  return 0.5 * (v10.cross(v20)).size() ;
209 }
210 
211 
212 
213 std::pair<Real, Real> Tri3::min_and_max_angle() const
214 {
215  Point v10 ( this->point(1) - this->point(0) );
216  Point v20 ( this->point(2) - this->point(0) );
217  Point v21 ( this->point(2) - this->point(1) );
218 
219  const Real
220  len_10=v10.size(),
221  len_20=v20.size(),
222  len_21=v21.size()
223  ;
224 
225  const Real
226  theta0=std::acos(( v10*v20)/len_10/len_20),
227  theta1=std::acos((-v10*v21)/len_10/len_21),
228  theta2=libMesh::pi - theta0 - theta1
229  ;
230 
231  libmesh_assert_greater (theta0, 0.);
232  libmesh_assert_greater (theta1, 0.);
233  libmesh_assert_greater (theta2, 0.);
234 
235  return std::make_pair(std::min(theta0, std::min(theta1,theta2)),
236  std::max(theta0, std::max(theta1,theta2)));
237 }
238 
239 } // namespace libMesh

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

Hosted By:
SourceForge.net Logo