GCC Code Coverage Report

Directory: Bembel/src/
File: Bembel/src/Maxwell/SingleLayerPotential.hpp
Date: 2024-03-19 14:38:05
Exec Total Coverage
Lines: 30 30 100.0%
Functions: 5 5 100.0%
Branches: 25 50 50.0%
Line Branch Exec Source
1 // This file is part of Bembel, the higher order C++ boundary element library.
2 //
3 // Copyright (C) 2022 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 #ifndef BEMBEL_SRC_MAXWELL_SINGLELAYERPOTENTIAL_HPP_
12 #define BEMBEL_SRC_MAXWELL_SINGLELAYERPOTENTIAL_HPP_
13
14 namespace Bembel {
15 // forward declaration of class MaxwellSingleLayerPotential in order to define
16 // traits
17 template <typename LinOp>
18 class MaxwellSingleLayerPotential;
19 /**
20 * \brief Specification of the PotentialTraits for Maxwell.
21 */
22 template <typename LinOp>
23 struct PotentialTraits<MaxwellSingleLayerPotential<LinOp>> {
24 typedef Eigen::VectorXcd::Scalar Scalar;
25 static constexpr int OutputSpaceDimension = 3;
26 };
27
28 /**
29 * \ingroup Maxwell
30 * \brief This class implements the specification of the integration for the
31 * single layer potential for Maxwell.
32 */
33 template <typename LinOp>
34 class MaxwellSingleLayerPotential
35 : public PotentialBase<MaxwellSingleLayerPotential<LinOp>, LinOp> {
36 // implementation of the kernel evaluation, which may be based on the
37 // information available from the superSpace
38 public:
39 2 MaxwellSingleLayerPotential() {}
40 Eigen::Matrix<typename PotentialReturnScalar<
41 typename LinearOperatorTraits<LinOp>::Scalar,
42 std::complex<double>>::Scalar,
43 3, 1>
44 91800 evaluateIntegrand_impl(const FunctionEvaluator<LinOp> &fun_ev,
45 const ElementTreeNode &element,
46 const Eigen::Vector3d &point,
47 const SurfacePoint &p) const {
48 // get evaluation points on unit square
49
1/2
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
 91800 auto s = p.segment<2>(0);
50
51 // get quadrature weights
52
1/2
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
 91800 auto ws = p(2);
53
54 // get points on geometry and tangential derivatives
55
1/2
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
 91800 auto x_f = p.segment<3>(3);
56
57 // compute surface measures from tangential derivatives
58 91800 auto h = element.get_h();
59
60 // evaluate kernel
61
2/4
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 91800 times.
✗ Branch 5 not taken.
 91800 auto kernel = evaluateKernel(point, x_f);
62
2/4
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 91800 times.
✗ Branch 5 not taken.
 91800 auto kernel_gradient = evaluateKernelGrad(point, x_f);
63
64 // assemble Galerkin solution
65
1/2
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
 91800 auto scalar_part = fun_ev.evaluate(element, p);
66
1/2
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
 91800 auto divergence_part = fun_ev.evaluateDiv(element, p);
67
68 // integrand without basis functions, note that the surface measure
69 // disappears for the divergence
70 // auto integrand = kernel * scalar_part * ws;
71
5/10
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 91800 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 91800 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 91800 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 91800 times.
✗ Branch 14 not taken.
 91800 auto integrand = (kernel * scalar_part +
72 91800 1. / wavenumber2_ * kernel_gradient * divergence_part) *
73 ws;
74
75
1/2
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
 183600 return integrand;
76 91800 }
77
78 /**
79 * \brief Fundamental solution of Helmholtz problem
80 */
81 91800 std::complex<double> evaluateKernel(const Eigen::Vector3d &x,
82 const Eigen::Vector3d &y) const {
83
2/4
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 91800 times.
✗ Branch 5 not taken.
 91800 auto r = (x - y).norm();
84
1/2
✓ Branch 2 taken 91800 times.
✗ Branch 3 not taken.
 91800 return std::exp(-std::complex<double>(0., 1.) * wavenumber_ * r) / 4. /
85 183600 BEMBEL_PI / r;
86 }
87 /**
88 * \brief Gradient of fundamental solution of Helmholtz problem
89 */
90 91800 Eigen::VectorXcd evaluateKernelGrad(const Eigen::Vector3d &x,
91 const Eigen::Vector3d &y) const {
92
1/2
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
 91800 auto c = x - y;
93
1/2
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
 91800 auto r = c.norm();
94 91800 auto r3 = r * r * r;
95 91800 auto i = std::complex<double>(0., 1.);
96
2/4
✓ Branch 4 taken 91800 times.
✗ Branch 5 not taken.
✓ Branch 9 taken 91800 times.
✗ Branch 10 not taken.
 91800 return (std::exp(-i * wavenumber_ * r) * (-1. - i * wavenumber_ * r) / 4. /
97
1/2
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
 183600 BEMBEL_PI / r3) *
98
2/4
✓ Branch 1 taken 91800 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 91800 times.
✗ Branch 5 not taken.
 91800 c;
99 }
100 //////////////////////////////////////////////////////////////////////////////
101 // setters
102 //////////////////////////////////////////////////////////////////////////////
103 2 void set_wavenumber(std::complex<double> wavenumber) {
104 2 wavenumber_ = wavenumber;
105 2 wavenumber2_ = wavenumber_ * wavenumber_;
106 2 }
107 //////////////////////////////////////////////////////////////////////////////
108 // getters
109 //////////////////////////////////////////////////////////////////////////////
110 std::complex<double> get_wavenumber() { return wavenumber_; }
111
112 private:
113 std::complex<double> wavenumber_;
114 std::complex<double> wavenumber2_;
115 };
116
117 } // namespace Bembel
118 #endif // BEMBEL_SRC_MAXWELL_SINGLELAYERPOTENTIAL_HPP_
119