12 #ifndef BEMBEL_SRC_AUGMENTEDEFIE_VECTORPOTENTIAL_HPP_
13 #define BEMBEL_SRC_AUGMENTEDEFIE_VECTORPOTENTIAL_HPP_
18 template <
typename LinOp>
19 class VectorPotential;
21 template <
typename LinOp>
23 typedef Eigen::VectorXcd::Scalar Scalar;
24 static constexpr
int OutputSpaceDimension = 3;
30 template <
typename LinOp>
38 std::complex<double>>::Scalar,
42 const Eigen::Vector3d &point,
45 auto s = p.segment<2>(0);
51 auto x_f = p.segment<3>(3);
57 auto scalar_part = fun_ev.evaluate(element, p);
61 auto integrand = kernel * scalar_part * ws;
70 const Eigen::Vector3d &y)
const {
71 auto r = (x - y).norm();
72 return std::exp(-std::complex<double>(0., 1.) * wavenumber_ * r) / 4. /
78 void set_wavenumber(std::complex<double> wavenumber) {
79 wavenumber_ = wavenumber;
84 std::complex<double> get_wavenumber() {
return wavenumber_; }
87 std::complex<double> wavenumber_;
The ElementTreeNode corresponds to an element in the element tree.
std::complex< double > evaluateKernel(const Eigen::Vector3d &x, const Eigen::Vector3d &y) const
Fundamental solution of Helmholtz problem.
Eigen::Matrix< double, 12, 1 > SurfacePoint
typedef of SurfacePoint
Routines for the evalutation of pointwise errors.
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.
Base case for specifying the return type of the potential.
struct containing specifications on the functional has to be specialized or derived for any particula...