From 023773986dedc3540f3a200af991f88ce8ed66d3 Mon Sep 17 00:00:00 2001 From: jdoelz <48793870+jdoelz@users.noreply.github.com> Date: Tue, 13 Feb 2024 14:18:31 +0100 Subject: [PATCH] Add TangentialTrace to Linear Forms with Unittest (#70) --- Bembel/LinearForm | 1 + Bembel/src/LinearForm/TangentialTrace.hpp | 88 +++++++++++++++++++++++ tests/CMakeLists.txt | 1 + tests/test_TangentialTrace.cpp | 51 +++++++++++++ 4 files changed, 141 insertions(+) create mode 100644 Bembel/src/LinearForm/TangentialTrace.hpp create mode 100644 tests/test_TangentialTrace.cpp diff --git a/Bembel/LinearForm b/Bembel/LinearForm index defb38b15..6ff3d2ecb 100644 --- a/Bembel/LinearForm +++ b/Bembel/LinearForm @@ -26,5 +26,6 @@ #include "src/LinearForm/DirichletTrace.hpp" #include "src/LinearForm/DiscreteLinearForm.hpp" #include "src/LinearForm/RotatedTangentialTrace.hpp" +#include "src/LinearForm/TangentialTrace.hpp" #endif // BEMBEL_LINEARFORM_MODULE_ diff --git a/Bembel/src/LinearForm/TangentialTrace.hpp b/Bembel/src/LinearForm/TangentialTrace.hpp new file mode 100644 index 000000000..fc79f9e6f --- /dev/null +++ b/Bembel/src/LinearForm/TangentialTrace.hpp @@ -0,0 +1,88 @@ +// This file is part of Bembel, the higher order C++ boundary element library. +// +// Copyright (C) 2022 see +// +// It was written as part of a cooperation of J. Doelz, H. Harbrecht, S. Kurz, +// M. Multerer, S. Schoeps, and F. Wolf at Technische Universitaet Darmstadt, +// Universitaet Basel, and Universita della Svizzera italiana, Lugano. This +// source code is subject to the GNU General Public License version 3 and +// provided WITHOUT ANY WARRANTY, see for further +// information. +#ifndef BEMBEL_SRC_LINEARFORM_TANGENTIALTRACE_HPP_ +#define BEMBEL_SRC_LINEARFORM_TANGENTIALTRACE_HPP_ + +namespace Bembel { + +template +class TangentialTrace; + +template +struct LinearFormTraits> { + typedef ScalarT Scalar; +}; + +/** + * \ingroup LinearForm + * \brief This class provides a specialization of the linear form required + *for the solution of the electric field integral equation. + **/ +template +class TangentialTrace : public LinearFormBase, Scalar> { + public: + TangentialTrace() {} + void set_function( + const std::function(Eigen::Vector3d)> + &function) { + function_ = function; + } + template + void evaluateIntegrand_impl( + const T &super_space, const SurfacePoint &p, + Eigen::Matrix *intval) const { + int polynomial_degree = super_space.get_polynomial_degree(); + int polynomial_degree_plus_one_squared = + (polynomial_degree + 1) * (polynomial_degree + 1); + + // get evaluation points on unit square + Eigen::Vector2d s = p.segment<2>(0); + + // get quadrature weights + double ws = p(2); + + // get points on geometry and tangential derivatives + Eigen::Vector3d x_f = p.segment<3>(3); + Eigen::Vector3d x_f_dx = p.segment<3>(6); + Eigen::Vector3d x_f_dy = p.segment<3>(9); + + // compute surface measures from tangential derivatives + Eigen::Vector3d x_n = x_f_dx.cross(x_f_dy).normalized(); + + // tangential component + quadrature weights + Eigen::Matrix fun_x_f = function_(x_f); + Eigen::Matrix tangential_component = fun_x_f.cross(x_n) * ws; + + // extract tangential component + Scalar component_x = x_f_dx.dot(tangential_component); + Scalar component_y = x_f_dy.dot(tangential_component); + + // evaluate shape functions + Eigen::Matrix phiPhiVec = + super_space.basis(s); + + // multiply basis functions with integrand + Eigen::Matrix phiPhiMat( + polynomial_degree_plus_one_squared, 2); + phiPhiMat.col(0) = component_x * phiPhiVec; + phiPhiMat.col(1) = component_y * phiPhiVec; + + // compute integrals + (*intval) += phiPhiMat; + return; + } + + private: + std::function(Eigen::Vector3d)> function_; +}; +} // namespace Bembel + +#endif // BEMBEL_SRC_LINEARFORM_TANGENTIALTRACE_HPP_ diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index f64e48356..85ecf5e2c 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -4,6 +4,7 @@ set(UNITTESTS test_SurfacePointUpdate test_Projector test_Glue + test_TangentialTrace test_FMMTransferMatrices test_FMMForwardTransformation test_FMMBackwardTransformation diff --git a/tests/test_TangentialTrace.cpp b/tests/test_TangentialTrace.cpp new file mode 100644 index 000000000..744f9e9a6 --- /dev/null +++ b/tests/test_TangentialTrace.cpp @@ -0,0 +1,51 @@ +// This file is part of Bembel, the higher order C++ boundary element library. +// +// Copyright (C) 2024 see +// +// It was written as part of a cooperation of J. Doelz, H. Harbrecht, S. Kurz, +// M. Multerer, S. Schoeps, and F. Wolf at Technische Universitaet Darmstadt, +// Universitaet Basel, and Universita della Svizzera italiana, Lugano. This +// source code is subject to the GNU General Public License version 3 and +// provided WITHOUT ANY WARRANTY, see for further +// information. + +#include +#include +#include + +#include "tests/TestGeometries.hpp" + +class TestOperatorDivC; +template <> +struct Bembel::LinearOperatorTraits { + typedef Eigen::VectorXd EigenType; + typedef Eigen::VectorXd::Scalar Scalar; + enum { OperatorOrder = 0, Form = DifferentialForm::DivConforming }; +}; + +int main() { + Test::TestGeometryWriter::writeScreen(); + Bembel::Geometry geometry("test_Screen.dat"); + + const int refinement_level = 0; + const int polynomial_degree = 1; + Bembel::AnsatzSpace ansatz_space(geometry, refinement_level, + polynomial_degree); + + std::function fun = [](Eigen::Vector3d in) { + return Eigen::Vector3d(in(0), in(1), in(2)); + }; + + Bembel::DiscreteLinearForm, TestOperatorDivC> + disc_lf(ansatz_space); + disc_lf.get_linear_form().set_function(fun); + disc_lf.compute(); + + Eigen::VectorXd ref_sol(4); + ref_sol << 0.25, 0.25, -0.25, -0.25; + + BEMBEL_TEST_IF((disc_lf.get_discrete_linear_form() - ref_sol).norm() < + Test::Constants::coefficient_accuracy); + + return 0; +}