libMesh
face_tri.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_tri.h"
22 #include "libmesh/edge_edge2.h"
23 #include "libmesh/face_tri3.h"
24 
25 namespace libMesh
26 {
27 
28 
29 // ------------------------------------------------------------
30 // Tri class static member initializations
31 
32 
33 // We need to require C++11...
34 const Real Tri::_master_points[6][3] =
35  {
36  {0, 0},
37  {1, 0},
38  {0, 1},
39  {0.5, 0},
40  {0.5, 0.5},
41  {0, 0.5}
42  };
43 
44 
45 
46 
47 
48 
49 
50 // ------------------------------------------------------------
51 // Tri class member functions
52 dof_id_type Tri::key (const unsigned int s) const
53 {
54  libmesh_assert_less (s, this->n_sides());
55 
56  return this->compute_key(this->node_id(Tri3::side_nodes_map[s][0]),
57  this->node_id(Tri3::side_nodes_map[s][1]));
58 }
59 
60 
61 
62 unsigned int Tri::which_node_am_i(unsigned int side,
63  unsigned int side_node) const
64 {
65  libmesh_assert_less (side, this->n_sides());
66  libmesh_assert_less (side_node, 2);
67 
68  return Tri3::side_nodes_map[side][side_node];
69 }
70 
71 
72 
74 {
75  return this->compute_key(this->node_id(0),
76  this->node_id(1),
77  this->node_id(2));
78 }
79 
80 
81 
82 UniquePtr<Elem> Tri::side_ptr (const unsigned int i)
83 {
84  libmesh_assert_less (i, this->n_sides());
85 
86  Elem * edge = new Edge2;
87 
88  for (unsigned n=0; n<edge->n_nodes(); ++n)
89  edge->set_node(n) = this->node_ptr(Tri3::side_nodes_map[i][n]);
90 
91  return UniquePtr<Elem>(edge);
92 }
93 
94 
95 
96 bool Tri::is_child_on_side(const unsigned int c,
97  const unsigned int s) const
98 {
99  libmesh_assert_less (c, this->n_children());
100  libmesh_assert_less (s, this->n_sides());
101 
102  return (c == s || c == (s+1)%3);
103 }
104 
105 
106 
108 {
109  switch (q)
110  {
111 
115  case DISTORTION:
116  case STRETCH:
117  {
118  const Point & p1 = this->point(0);
119  const Point & p2 = this->point(1);
120  const Point & p3 = this->point(2);
121 
122  Point v1 = p2 - p1;
123  Point v2 = p3 - p1;
124  Point v3 = p3 - p2;
125  const Real l1 = v1.norm();
126  const Real l2 = v2.norm();
127  const Real l3 = v3.norm();
128 
129  // if one length is 0, quality is quite bad!
130  if ((l1 <=0.) || (l2 <= 0.) || (l3 <= 0.))
131  return 0.;
132 
133  const Real s1 = std::sin(std::acos(v1*v2/l1/l2)/2.);
134  v1 *= -1;
135  const Real s2 = std::sin(std::acos(v1*v3/l1/l3)/2.);
136  const Real s3 = std::sin(std::acos(v2*v3/l2/l3)/2.);
137 
138  return 8. * s1 * s2 * s3;
139 
140  }
141  default:
142  return Elem::quality(q);
143  }
144 
150  return Elem::quality(q);
151 
152 }
153 
154 
155 
156 
157 
158 
159 std::pair<Real, Real> Tri::qual_bounds (const ElemQuality q) const
160 {
161  std::pair<Real, Real> bounds;
162 
163  switch (q)
164  {
165 
166  case MAX_ANGLE:
167  bounds.first = 60.;
168  bounds.second = 90.;
169  break;
170 
171  case MIN_ANGLE:
172  bounds.first = 30.;
173  bounds.second = 60.;
174  break;
175 
176  case CONDITION:
177  bounds.first = 1.;
178  bounds.second = 1.3;
179  break;
180 
181  case JACOBIAN:
182  bounds.first = 0.5;
183  bounds.second = 1.155;
184  break;
185 
186  case SIZE:
187  case SHAPE:
188  bounds.first = 0.25;
189  bounds.second = 1.;
190  break;
191 
192  case DISTORTION:
193  bounds.first = 0.6;
194  bounds.second = 1.;
195  break;
196 
197  default:
198  libMesh::out << "Warning: Invalid quality measure chosen." << std::endl;
199  bounds.first = -1;
200  bounds.second = -1;
201  }
202 
203  return bounds;
204 }
205 
206 } // namespace libMesh
static const Real _master_points[6][3]
Master element node locations.
Definition: face_tri.h:164
virtual bool is_child_on_side(const unsigned int c, const unsigned int s) const libmesh_override
Definition: face_tri.C:96
virtual Node *& set_node(const unsigned int i)
Definition: elem.h:1941
UniquePtr< Elem > side(const unsigned int i) const
Definition: elem.h:2105
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
Real norm() const
Definition: type_vector.h:909
The libMesh namespace provides an interface to certain functionality in the library.
static const unsigned int side_nodes_map[3][2]
This maps the node of the side to element node numbers.
Definition: face_tri3.h:133
virtual Real quality(const ElemQuality q) const libmesh_override
Definition: face_tri.C:107
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual unsigned int n_nodes() const =0
virtual unsigned int n_children() const libmesh_override
Definition: face_tri.h:101
const Node * node_ptr(const unsigned int i) const
Definition: elem.h:1874
virtual Real quality(const ElemQuality q) const
Definition: elem.C:1547
virtual UniquePtr< Elem > side_ptr(const unsigned int i) libmesh_override
Definition: face_tri.C:82
virtual unsigned int n_sides() const libmesh_override
Definition: face_tri.h:86
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual dof_id_type key() const libmesh_override
Definition: face_tri.C:73
const Point & point(const unsigned int i) const
Definition: elem.h:1809
virtual std::pair< Real, Real > qual_bounds(const ElemQuality q) const libmesh_override
Definition: face_tri.C:159
OStreamProxy out
virtual unsigned int which_node_am_i(unsigned int side, unsigned int side_node) const libmesh_override
Definition: face_tri.C:62
The Edge2 is an element in 1D composed of 2 nodes.
Definition: edge_edge2.h:43
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
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
uint8_t dof_id_type
Definition: id_types.h:64