libMesh
edge_edge3.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 
19 
20 // Local includes
21 #include "libmesh/edge_edge3.h"
22 
23 namespace libMesh
24 {
25 
26 #ifdef LIBMESH_ENABLE_AMR
27 
28 const float Edge3::_embedding_matrix[2][3][3] =
29  {
30  // embedding matrix for child 0
31  {
32  // 0 1 2
33  {1.0, 0.0, 0.0}, // left
34  {0.0, 0.0, 1.0}, // right
35  {0.375,-0.125,0.75} // middle
36  },
37 
38  // embedding matrix for child 1
39  {
40  // 0 1 2
41  {0.0, 0.0, 1.0}, // left
42  {0.0, 1.0, 0.0}, // right
43  {-0.125,0.375,0.75} // middle
44  }
45  };
46 
47 #endif
48 
49 bool Edge3::is_vertex(const unsigned int i) const
50 {
51  if (i < 2)
52  return true;
53  return false;
54 }
55 
56 bool Edge3::is_edge(const unsigned int i) const
57 {
58  if (i < 2)
59  return false;
60  return true;
61 }
62 
63 bool Edge3::is_face(const unsigned int ) const
64 {
65  return false;
66 }
67 
68 bool Edge3::is_node_on_side(const unsigned int n,
69  const unsigned int s) const
70 {
71  libmesh_assert_less (s, 2);
72  libmesh_assert_less (n, 3);
73  return (s == n);
74 }
75 
76 bool Edge3::is_node_on_edge(const unsigned int,
77  const unsigned int libmesh_dbg_var(e)) const
78 {
79  libmesh_assert_equal_to (e, 0);
80  return true;
81 }
82 
83 
84 
86 {
87  return (this->point(2).relative_fuzzy_equals
88  ((this->point(0) + this->point(1))/2));
89 }
90 
91 
92 
93 void Edge3::connectivity(const unsigned int sc,
94  const IOPackage iop,
95  std::vector<dof_id_type> & conn) const
96 {
97  libmesh_assert_less_equal (sc, 1);
98  libmesh_assert_less (sc, this->n_sub_elem());
99  libmesh_assert_not_equal_to (iop, INVALID_IO_PACKAGE);
100 
101  // Create storage
102  conn.resize(2);
103 
104  switch (iop)
105  {
106  case TECPLOT:
107  {
108  switch (sc)
109  {
110  case 0:
111  conn[0] = this->node_id(0)+1;
112  conn[1] = this->node_id(2)+1;
113  return;
114 
115  case 1:
116  conn[0] = this->node_id(2)+1;
117  conn[1] = this->node_id(1)+1;
118  return;
119 
120  default:
121  libmesh_error_msg("Invalid sc = " << sc);
122  }
123  }
124 
125 
126  case VTK:
127  {
128  conn.resize(3);
129  conn[0] = this->node_id(0);
130  conn[1] = this->node_id(1);
131  conn[2] = this->node_id(2);
132  return;
133 
134  /*
135  switch (sc)
136  {
137  case 0:
138  conn[0] = this->node_id(0);
139  conn[1] = this->node_id(2);
140 
141  return;
142 
143  case 1:
144  conn[0] = this->node_id(2);
145  conn[1] = this->node_id(1);
146 
147  return;
148 
149  default:
150  libmesh_error_msg("Invalid sc = " << sc);
151  }
152  */
153  }
154 
155  default:
156  libmesh_error_msg("Unsupported IO package " << iop);
157  }
158 }
159 
160 
161 
162 std::pair<unsigned short int, unsigned short int>
163 Edge3::second_order_child_vertex (const unsigned int) const
164 {
165  return std::pair<unsigned short int, unsigned short int>(0,0);
166 }
167 
168 
169 
171 {
172  // Finding the (exact) length of a general quadratic element
173  // is a surprisingly complicated formula.
174  Point A = this->point(0) + this->point(1) - 2*this->point(2);
175  Point B = (this->point(1) - this->point(0))/2;
176 
177  const Real a = A.norm_sq();
178  const Real c = B.norm_sq();
179 
180  // Degenerate straight line case
181  if (a < TOLERANCE*TOLERANCE)
182  return 2. * std::sqrt(c);
183 
184  const Real b = 2.*(A*B);
185  const Real ba=b/a;
186  const Real ca=c/a;
187 
188  libmesh_assert (1.-ba+ca>0.);
189 
190  const Real s1 = std::sqrt(1. - ba + ca);
191  const Real s2 = std::sqrt(1. + ba + ca);
192 
193  Real log_term = (1. - 0.5*ba + s1) / (-1. - 0.5*ba + s2);
194  libmesh_assert(!libmesh_isnan(log_term) && log_term > 0.);
195 
196  return 0.5*std::sqrt(a)*((1.-0.5*ba)*s1 +
197  (1.+0.5*ba)*s2 +
198  (ca - 0.25*ba*ba)*std::log(log_term)
199  );
200 }
201 
202 
203 
205 {
206  // This might be a curved line through 2-space or 3-space, in which
207  // case the full bounding box can be larger than the bounding box of
208  // just the nodes.
209  Point pmin, pmax;
210 
211  for (unsigned d=0; d<LIBMESH_DIM; ++d)
212  {
213  Real center = this->point(2)(d);
214  Real hd = std::max(std::abs(center - this->point(0)(d)),
215  std::abs(center - this->point(1)(d)));
216 
217  pmin(d) = center - hd;
218  pmax(d) = center + hd;
219  }
220 
221  return BoundingBox(pmin, pmax);
222 }
223 
224 
225 
227 {
228  return this->compute_key(this->node_id(2));
229 }
230 
231 
232 } // namespace libMesh
virtual void connectivity(const unsigned int sc, const IOPackage iop, std::vector< dof_id_type > &conn) const libmesh_override
Definition: edge_edge3.C:93
double abs(double a)
virtual bool is_face(const unsigned int i) const libmesh_override
Definition: edge_edge3.C:63
virtual Real volume() const libmesh_override
An optimized method for computing the length of a 3-node edge.
Definition: edge_edge3.C:170
virtual bool is_node_on_edge(const unsigned int n, const unsigned int e) const libmesh_override
Definition: edge_edge3.C:76
static const Real TOLERANCE
The libMesh namespace provides an interface to certain functionality in the library.
long double max(long double a, double b)
libmesh_assert(j)
Definition: assembly.h:38
virtual std::pair< unsigned short int, unsigned short int > second_order_child_vertex(const unsigned int n) const libmesh_override
Definition: edge_edge3.C:163
virtual bool has_affine_map() const libmesh_override
Definition: edge_edge3.C:85
virtual bool is_edge(const unsigned int i) const libmesh_override
Definition: edge_edge3.C:56
PetscErrorCode Vec Mat libmesh_dbg_var(j)
Real norm_sq() const
Definition: type_vector.h:940
Defines a Cartesian bounding box by the two corner extremum.
Definition: bounding_box.h:40
virtual unsigned int n_sub_elem() const libmesh_override
Definition: edge_edge3.h:74
virtual bool is_vertex(const unsigned int i) const libmesh_override
Definition: edge_edge3.C:49
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const Point & point(const unsigned int i) const
Definition: elem.h:1809
static const float _embedding_matrix[2][3][3]
Matrix that computes new nodal locations/solution values from current nodes/solution.
Definition: edge_edge3.h:204
virtual dof_id_type key() const libmesh_override
Compute a unique key for this element which is suitable for hashing (not necessarily unique...
Definition: edge_edge3.C:226
bool libmesh_isnan(float a)
static PetscErrorCode Mat * A
virtual BoundingBox loose_bounding_box() const libmesh_override
Definition: edge_edge3.C:204
dof_id_type node_id(const unsigned int i) const
Definition: elem.h:1831
static dof_id_type compute_key(dof_id_type n0)
Definition: elem.h:2621
IOPackage
libMesh interfaces with several different software packages for the purposes of creating, reading, and writing mesh files.
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
virtual bool is_node_on_side(const unsigned int n, const unsigned int s) const libmesh_override
Definition: edge_edge3.C:68
uint8_t dof_id_type
Definition: id_types.h:64