NiHu  2.0
elastostatics_singular_integrals.hpp
Go to the documentation of this file.
1 // This file is a part of NiHu, a C++ BEM template library.
2 //
3 // Copyright (C) 2012-2014 Peter Fiala <fiala@hit.bme.hu>
4 // Copyright (C) 2012-2014 Peter Rucz <rucz@hit.bme.hu>
5 //
6 // This program is free software: you can redistribute it and/or modify
7 // it under the terms of the GNU General Public License as published by
8 // the Free Software Foundation, either version 3 of the License, or
9 // (at your option) any later version.
10 //
11 // This program is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 // GNU General Public License for more details.
15 //
16 // You should have received a copy of the GNU General Public License
17 // along with this program. If not, see <http://www.gnu.org/licenses/>.
18 
22 #ifndef ELASTOSTATICS_SINGULAR_INTEGRALS_HPP_INCLUDED
23 #define ELASTOSTATICS_SINGULAR_INTEGRALS_HPP_INCLUDED
24 
25 #include <boost/math/constants/constants.hpp>
26 
27 #include "../core/integral_operator.hpp"
28 #include "../core/singular_integral_shortcut.hpp"
29 
30 #include "elastostatics_kernel.hpp"
31 #include "field_type_helpers.hpp"
32 #include "guiggiani_1992.hpp"
33 #include "lib_element.hpp"
34 
35 
36 namespace NiHu
37 {
38 
41 {
42  typedef Eigen::Matrix<double, 2, 2> result_t;
43 public:
50  static result_t eval(line_1_elem const &elem, double nu, double mu)
51  {
52  using boost::math::double_constants::pi;
53 
54  auto const &C = elem.get_coords();
55  auto rvec = (C.col(1) - C.col(0)).eval();
56  auto r = rvec.norm(); // elem length
57  auto gradr = rvec.normalized();
58 
59  return (
60  - (3.-4.*nu) * result_t::Identity() * (r * (std::log(r/2) - 1.))
61  +
62  (gradr * gradr.transpose()) * r
63  ) / (8. * pi * mu * (1. - nu));
64  }
65 };
66 
71 template <class TestField, class TrialField>
73  elastostatics_2d_U_kernel, TestField, TrialField, match::match_1d_type,
74  typename std::enable_if<
75  is_collocational<TestField, TrialField>::value &&
76  is_constant_line<TrialField>::value
77  >::type
78 >
79 {
80 public:
88  template <class result_t>
89  static result_t &eval(
90  result_t &result,
92  field_base<TestField> const &,
93  field_base<TrialField> const &trial_field,
94  element_match const &)
95  {
97  trial_field.get_elem(), kernel.derived().get_poisson_ratio(), kernel.derived().get_shear_modulus());
98  return result;
99  }
100 };
101 
102 
107 template <class TestField, class TrialField>
109  elastostatics_2d_T_kernel, TestField, TrialField, match::match_1d_type,
110  typename std::enable_if<
111  is_collocational<TestField, TrialField>::value &&
112  is_constant_line<TrialField>::value
113  >::type
114 >
115 {
116 public:
124  template <class result_t>
125  static result_t &eval(
126  result_t &result,
128  field_base<TestField> const &,
129  field_base<TrialField> const &,
130  element_match const &)
131  {
132  result.setZero();
133  return result;
134  }
135 };
136 
137 
142 template <class TestField, class TrialField>
144  elastostatics_3d_T_kernel, TestField, TrialField, match::match_2d_type,
145  typename std::enable_if<
146  std::is_same<typename get_formalism<TestField, TrialField>::type, formalism::collocational>::value
147  >::type
148 >
149 {
150 public:
157  template <class result_t>
158  static result_t &eval(
159  result_t &result,
161  field_base<TestField> const &,
162  field_base<TrialField> const &trial_field,
163  element_match const &)
164  {
166  auto const &elem = trial_field.get_elem();
167  guiggiani_t gui(elem, kernel.derived());
168 
169  auto const &xi0 = TestField::nset_t::corner_at(0);
170  gui.integrate(result, xi0, elem.get_normal(xi0));
171 
172  return result;
173  }
174 };
175 
176 
177 
178 
179 
180 
183 {
184  typedef Eigen::Matrix<double, 2, 2> result_t;
185 public:
192  static result_t eval(line_1_elem const &elem, double nu, double mu)
193  {
194  using namespace boost::math::double_constants;
195 
196  auto const &C = elem.get_coords();
197  auto rvec = (C.col(1) - C.col(0)).eval();
198  auto r = rvec.norm(); // elem length
199  auto gradr = rvec.normalized();
200 
201  return r*r * (
202  -(3. - 4.*nu) * result_t::Identity() * (std::log(r) - 1.5)
203  + gradr * gradr.transpose()
204  ) / (8. * pi * mu * (1. - nu));
205  }
206 };
207 
208 
213 template <class TestField, class TrialField>
215  elastostatics_2d_U_kernel, TestField, TrialField, match::match_1d_type,
216  typename std::enable_if<
217  std::is_same<typename get_formalism<TestField, TrialField>::type, formalism::general>::value &&
218  std::is_same<typename TrialField::elem_t::lset_t, line_1_shape_set>::value &&
219  std::is_same<typename TestField::nset_t, line_0_shape_set>::value &&
220  std::is_same<typename TrialField::nset_t, line_0_shape_set>::value
221  >::type
222 >
223 {
224 public:
232  template <class result_t>
233  static result_t &eval(
234  result_t &result,
236  field_base<TestField> const &,
237  field_base<TrialField> const &trial_field,
238  element_match const &)
239  {
241  trial_field.get_elem(), kernel.derived().get_poisson_ratio(), kernel.derived().get_shear_modulus());
242  return result;
243  }
244 };
245 
246 
247 
248 
253 template <class TestField, class TrialField>
255  elastostatics_2d_T_kernel, TestField, TrialField, match::match_1d_type,
256  typename std::enable_if<
257  std::is_same<typename get_formalism<TestField, TrialField>::type, formalism::general>::value &&
258  std::is_same<typename TestField::elem_t::lset_t, line_1_shape_set>::value &&
259  std::is_same<typename TrialField::elem_t::lset_t, line_1_shape_set>::value &&
260  std::is_same<typename TestField::nset_t, line_0_shape_set>::value &&
261  std::is_same<typename TrialField::nset_t, line_0_shape_set>::value
262  >::type
263 >
264 {
265 public:
274  template <class result_t>
275  static result_t &eval(
276  result_t &result,
278  field_base<TestField> const &test_field,
279  field_base<TrialField> const &trial_field,
280  element_match const &)
281  {
282  result.setZero();
283  return result;
284  }
285 };
286 
287 
288 } // end of namespace NiHu
289 
290 
291 #endif // ELASTOSTATICS_SINGULAR_INTEGRALS_HPP_INCLUDED
292 
elastostatics_kernel.hpp
implementation of fundamental solutions of elastostatics
NiHu::elastostatics_2d_U_galerkin_face_constant_line
Galerkin face match integral of 2D Elastostatics U kernel over constant line.
Definition: elastostatics_singular_integrals.hpp:182
NiHu::field_base::get_elem
const elem_t & get_elem() const
return underlying element
Definition: field.hpp:149
NiHu::elastostatics_2d_T_kernel
2d traction kernel for elastostatics
Definition: elastostatics_kernel.hpp:188
NiHu::kernel_base< elastostatics_2d_U_kernel >
NiHu::singular_integral_shortcut< elastostatics_2d_U_kernel, TestField, TrialField, match::match_1d_type, typename std::enable_if< is_collocational< TestField, TrialField >::value &&is_constant_line< TrialField >::value >::type >::eval
static result_t & eval(result_t &result, kernel_base< elastostatics_2d_U_kernel > const &kernel, field_base< TestField > const &, field_base< TrialField > const &trial_field, element_match const &)
evaluate singular integral
Definition: elastostatics_singular_integrals.hpp:89
NiHu::singular_integral_shortcut< elastostatics_2d_T_kernel, TestField, TrialField, match::match_1d_type, typename std::enable_if< std::is_same< typename get_formalism< TestField, TrialField >::type, formalism::general >::value &&std::is_same< typename TestField::elem_t::lset_t, line_1_shape_set >::value &&std::is_same< typename TrialField::elem_t::lset_t, line_1_shape_set >::value &&std::is_same< typename TestField::nset_t, line_0_shape_set >::value &&std::is_same< typename TrialField::nset_t, line_0_shape_set >::value >::type >::eval
static result_t & eval(result_t &result, kernel_base< elastostatics_2d_T_kernel > const &kernel, field_base< TestField > const &test_field, field_base< TrialField > const &trial_field, element_match const &)
evaluate singular integral
Definition: elastostatics_singular_integrals.hpp:275
NiHu::singular_integral_shortcut
a shortcut for the user to define customised singular integral methods
Definition: double_integral.hpp:709
NiHu::singular_integral_shortcut< elastostatics_2d_U_kernel, TestField, TrialField, match::match_1d_type, typename std::enable_if< std::is_same< typename get_formalism< TestField, TrialField >::type, formalism::general >::value &&std::is_same< typename TrialField::elem_t::lset_t, line_1_shape_set >::value &&std::is_same< typename TestField::nset_t, line_0_shape_set >::value &&std::is_same< typename TrialField::nset_t, line_0_shape_set >::value >::type >::eval
static result_t & eval(result_t &result, kernel_base< elastostatics_2d_U_kernel > const &kernel, field_base< TestField > const &, field_base< TrialField > const &trial_field, element_match const &)
evaluate singular integral
Definition: elastostatics_singular_integrals.hpp:233
NiHu::match::match_2d_type
std::integral_constant< int, 2 > match_2d_type
two elements are adjacent with a 2d (surface) match
Definition: match_types.hpp:42
NiHu::field_base
CRTP base class of all fields.
Definition: field.hpp:111
NiHu::surface_element
class describing a surface element that provides a normal vector
Definition: element.hpp:451
lib_element.hpp
guiggiani_1992.hpp
Guiggiani's method for CPV and HPF collocational integrals.
NiHu::match::match_1d_type
std::integral_constant< int, 1 > match_1d_type
two elements are adjacent with a 1d (line) match
Definition: match_types.hpp:40
NiHu::elastostatics_2d_U_collocation_constant_line::eval
static result_t eval(line_1_elem const &elem, double nu, double mu)
Evaluate the integral.
Definition: elastostatics_singular_integrals.hpp:50
NiHu::singular_integral_shortcut< elastostatics_2d_T_kernel, TestField, TrialField, match::match_1d_type, typename std::enable_if< is_collocational< TestField, TrialField >::value &&is_constant_line< TrialField >::value >::type >::eval
static result_t & eval(result_t &result, kernel_base< elastostatics_2d_T_kernel > const &, field_base< TestField > const &, field_base< TrialField > const &, element_match const &)
evaluate singular integral
Definition: elastostatics_singular_integrals.hpp:125
NiHu::element_match
class describing the adjacency (match) state of two elements
Definition: element_match.hpp:39
NiHu::elastostatics_2d_U_collocation_constant_line
Collocational singular integral of 2D Elastostatics U kernel over constant line.
Definition: elastostatics_singular_integrals.hpp:40
NiHu::elastostatics_2d_U_galerkin_face_constant_line::eval
static result_t eval(line_1_elem const &elem, double nu, double mu)
Evaluate the integral.
Definition: elastostatics_singular_integrals.hpp:192
NiHu::elastostatics_2d_U_kernel
The displacement kernel of 2D elastostatics.
Definition: elastostatics_kernel.hpp:115
NiHu::guiggiani
Implementation of Guiggiani's method.
Definition: guiggiani_1992.hpp:84
NiHu::elastostatics_3d_T_kernel
the traction kernel for 3d isotropic elastostatics
Definition: elastostatics_kernel.hpp:329
NiHu::singular_integral_shortcut< elastostatics_3d_T_kernel, TestField, TrialField, match::match_2d_type, typename std::enable_if< std::is_same< typename get_formalism< TestField, TrialField >::type, formalism::collocational >::value >::type >::eval
static result_t & eval(result_t &result, kernel_base< elastostatics_3d_T_kernel > const &kernel, field_base< TestField > const &, field_base< TrialField > const &trial_field, element_match const &)
evaluate singular integral
Definition: elastostatics_singular_integrals.hpp:158
C
Definition: bbfmm_covariance.cpp:47