Bembel
EFIE.hpp
1 // This file is part of Bembel, the higher order C++ boundary element library.
2 //
3 // Copyright (C) 2024 see <http://www.bembel.eu>
4 //
5 // It was written as part of a cooperation of J. Doelz, H. Harbrecht, S. Kurz,
6 // M. Multerer, S. Schoeps, and F. Wolf at Technische Universitaet Darmstadt,
7 // Universitaet Basel, and Universita della Svizzera italiana, Lugano. This
8 // source code is subject to the GNU General Public License version 3 and
9 // provided WITHOUT ANY WARRANTY, see <http://www.bembel.eu> for further
10 // information.
11 
12 #ifndef BEMBEL_SRC_AUGMENTEDEFIE_EFIE_HPP_
13 #define BEMBEL_SRC_AUGMENTEDEFIE_EFIE_HPP_
14 
15 namespace Bembel {
16 // forward declaration of class EFIE in order to define
17 // traits
18 template <typename LinOp>
19 class EFIE;
20 
21 template <typename LinOp>
22 struct PotentialTraits<EFIE<LinOp>> {
23  typedef Eigen::VectorXcd::Scalar Scalar;
24  static constexpr int OutputSpaceDimension = 3;
25 };
26 
30 template <typename LinOp>
31 class EFIE : public PotentialBase<EFIE<LinOp>, LinOp> {
32  // implementation of the kernel evaluation, which may be based on the
33  // information available from the superSpace
34  public:
35  EFIE() {}
36  Eigen::Matrix<typename PotentialReturnScalar<
38  std::complex<double>>::Scalar,
39  3, 1>
40  evaluateIntegrand_impl(const FunctionEvaluator<LinOp> &fun_ev,
41  const ElementTreeNode &element,
42  const Eigen::Vector3d &point,
43  const SurfacePoint &p) const {
44  // get evaluation points on unit square
45  auto s = p.segment<2>(0);
46 
47  // get quadrature weights
48  auto ws = p(2);
49 
50  // get points on geometry and tangential derivatives
51  auto x_f = p.segment<3>(3);
52 
53  // compute surface measures from tangential derivatives
54  auto h = element.get_h();
55 
56  // evaluate kernel
57  auto kernel = evaluateKernel(point, x_f);
58  auto kernel_gradient = evaluateKernelGrad(point, x_f);
59 
60  // assemble Galerkin solution
61  auto scalar_part = fun_ev.evaluate(element, p);
62  auto divergence_part = fun_ev.evaluateDiv(element, p);
63 
64  // integrand without basis functions, note that the surface measure
65  // disappears for the divergence
66  std::complex<double> omega =
67  wavenumber_ / std::sqrt(Constants::mu0 * Constants::eps0);
68  auto integrand = -(std::complex<double>(0., 1.) * omega * Constants::mu0 *
69  kernel * scalar_part -
70  kernel_gradient * divergence_part / Constants::eps0 /
71  std::complex<double>(0., 1.) / omega) *
72  ws;
73 
74  return integrand;
75  }
76 
80  std::complex<double> evaluateKernel(const Eigen::Vector3d &x,
81  const Eigen::Vector3d &y) const {
82  auto r = (x - y).norm();
83  return std::exp(-std::complex<double>(0., 1.) * wavenumber_ * r) / 4. /
84  BEMBEL_PI / r;
85  }
89  Eigen::VectorXcd evaluateKernelGrad(const Eigen::Vector3d &x,
90  const Eigen::Vector3d &y) const {
91  auto c = x - y;
92  auto r = c.norm();
93  auto r3 = r * r * r;
94  auto i = std::complex<double>(0., 1.);
95  return (std::exp(-i * wavenumber_ * r) * (-1. - i * wavenumber_ * r) / 4. /
96  BEMBEL_PI / r3) *
97  c;
98  }
100  // setters
102  void set_wavenumber(std::complex<double> wavenumber) {
103  wavenumber_ = wavenumber;
104  }
106  // getters
108  std::complex<double> get_wavenumber() { return wavenumber_; }
109 
110  private:
111  std::complex<double> wavenumber_;
112 };
113 
114 } // namespace Bembel
115 #endif // BEMBEL_SRC_AUGMENTEDEFIE_EFIE_HPP_
Eigen::VectorXcd evaluateKernelGrad(const Eigen::Vector3d &x, const Eigen::Vector3d &y) const
Gradient of fundamental solution of Helmholtz problem.
Definition: EFIE.hpp:89
std::complex< double > evaluateKernel(const Eigen::Vector3d &x, const Eigen::Vector3d &y) const
Fundamental solution of Helmholtz problem.
Definition: EFIE.hpp:80
The ElementTreeNode corresponds to an element in the element tree.
double get_h() const
getter
Eigen::Matrix< double, 12, 1 > SurfacePoint
typedef of SurfacePoint
Routines for the evalutation of pointwise errors.
Definition: AnsatzSpace.hpp:14
struct containing specifications on the linear operator has to be specialized or derived for any part...
functional base class. this serves as a common interface for existing functionals.
Definition: Potential.hpp:81
Base case for specifying the return type of the potential.
Definition: Potential.hpp:36
struct containing specifications on the functional has to be specialized or derived for any particula...
Definition: Potential.hpp:28