Bembel
GradScalarPotential.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_GRADSCALARPOTENTIAL_HPP_
13 #define BEMBEL_SRC_AUGMENTEDEFIE_GRADSCALARPOTENTIAL_HPP_
14 
15 namespace Bembel {
16 // forward declaration of class GradScalarPotential in order to define
17 // traits
18 template <typename LinOp>
19 class GradScalarPotential;
20 
21 template <typename LinOp>
23  typedef Eigen::VectorXcd::Scalar Scalar;
24  static constexpr int OutputSpaceDimension = 3;
25 };
26 
30 template <typename LinOp>
32  : public PotentialBase<GradScalarPotential<LinOp>, LinOp> {
33  // implementation of the kernel evaluation, which may be based on the
34  // information available from the superSpace
35  public:
37  Eigen::Matrix<typename PotentialReturnScalar<
39  std::complex<double>>::Scalar,
40  3, 1>
41  evaluateIntegrand_impl(const FunctionEvaluator<LinOp> &fun_ev,
42  const ElementTreeNode &element,
43  const Eigen::Vector3d &point,
44  const SurfacePoint &p) const {
45  // get evaluation points on unit square
46  auto s = p.segment<2>(0);
47 
48  // get quadrature weights
49  auto ws = p(2);
50 
51  // get points on geometry and tangential derivatives
52  auto x_f = p.segment<3>(3);
53  auto x_f_dx = p.segment<3>(6);
54  auto x_f_dy = p.segment<3>(9);
55 
56  // compute surface measures from tangential derivatives
57  auto x_kappa = x_f_dx.cross(x_f_dy).norm();
58 
59  // evaluate gradient of kernel
60  auto kernel_gradient = evaluateKernelGrad(point, x_f);
61 
62  // assemble Galerkin solution
63  auto cauchy_data = fun_ev.evaluate(element, p);
64 
65  // integrand without basis functions, note that the surface measure
66  // disappears for the divergence
67  auto integrand = kernel_gradient * cauchy_data * ws * x_kappa;
68 
69  return integrand;
70  }
74  Eigen::VectorXcd evaluateKernelGrad(const Eigen::Vector3d &x,
75  const Eigen::Vector3d &y) const {
76  auto c = x - y;
77  auto r = c.norm();
78  auto r3 = r * r * r;
79  auto i = std::complex<double>(0., 1.);
80  return (std::exp(-i * wavenumber_ * r) * (-1. - i * wavenumber_ * r) / 4. /
81  BEMBEL_PI / r3) *
82  c;
83  }
85  // setters
87  void set_wavenumber(std::complex<double> wavenumber) {
88  wavenumber_ = wavenumber;
89  }
91  // getters
93  std::complex<double> get_wavenumber() { return wavenumber_; }
94 
95  private:
96  std::complex<double> wavenumber_;
97 };
98 
99 } // namespace Bembel
100 #endif // BEMBEL_SRC_AUGMENTEDEFIE_GRADSCALARPOTENTIAL_HPP_
The ElementTreeNode corresponds to an element in the element tree.
Eigen::VectorXcd evaluateKernelGrad(const Eigen::Vector3d &x, const Eigen::Vector3d &y) const
Gradient of fundamental solution of Helmholtz problem.
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