Commit b90fd284 authored by ESTERIE Pierre's avatar ESTERIE Pierre
Browse files

Add log matrix kernel in 2d

parent a8cbd2a0
......@@ -14,7 +14,6 @@
#include <utility>
#include <xsimd/math/xsimd_scalar.hpp>
#include <xtensor/xmath.hpp>
#include <xtensor/xfixed.hpp>
#include <xtensor/xmath.hpp>
......@@ -42,7 +41,7 @@ namespace scalfmm::matrix_kernels::laplace
static constexpr std::size_t km{1};
static constexpr std::size_t kn{1};
template<typename ValueType>
using matrix_type = std::array<ValueType, kn*km>;
using matrix_type = std::array<ValueType, kn * km>;
template<typename ValueType>
using vector_type = std::array<ValueType, kn>;
......@@ -98,7 +97,7 @@ namespace scalfmm::matrix_kernels::laplace
static constexpr std::size_t km{2};
static constexpr std::size_t kn{2};
template<typename ValueType>
using matrix_type = std::array<ValueType, kn*km>;
using matrix_type = std::array<ValueType, kn * km>;
template<typename ValueType>
using vector_type = std::array<ValueType, kn>;
......@@ -160,7 +159,7 @@ namespace scalfmm::matrix_kernels::laplace
static constexpr std::size_t km{2};
static constexpr std::size_t kn{2};
template<typename ValueType>
using matrix_type = std::array<ValueType, kn*km>;
using matrix_type = std::array<ValueType, kn * km>;
template<typename ValueType>
using vector_type = std::array<ValueType, kn>;
......@@ -222,7 +221,7 @@ namespace scalfmm::matrix_kernels::laplace
static constexpr std::size_t km{1};
static constexpr std::size_t kn{Dim};
template<typename ValueType>
using matrix_type = std::array<ValueType, kn*km>;
using matrix_type = std::array<ValueType, kn * km>;
template<typename ValueType>
using vector_type = std::array<ValueType, kn>;
......@@ -292,7 +291,7 @@ namespace scalfmm::matrix_kernels::laplace
static constexpr std::size_t km{1};
static constexpr std::size_t kn{1 + Dim};
template<typename ValueType>
using matrix_type = std::array<ValueType, kn*km>;
using matrix_type = std::array<ValueType, kn * km>;
template<typename ValueType>
using vector_type = std::array<ValueType, kn>;
......@@ -357,7 +356,7 @@ namespace scalfmm::matrix_kernels::laplace
static constexpr std::size_t km{1};
static constexpr std::size_t kn{1};
template<typename ValueType>
using matrix_type = std::array<ValueType, kn*km>;
using matrix_type = std::array<ValueType, kn * km>;
template<typename ValueType>
using vector_type = std::array<ValueType, kn>;
......@@ -372,15 +371,15 @@ namespace scalfmm::matrix_kernels::laplace
container::point<ValueType, 2> const& y) const noexcept
{
using std::log;
return matrix_type<ValueType>{
log(std::sqrt((x.at(0) - y.at(0)) * (x.at(0) - y.at(0)) + (x.at(1) - y.at(1)) * (x.at(1) - y.at(1))))};
return matrix_type<ValueType>{xsimd::log(
xsimd::sqrt((x.at(0) - y.at(0)) * (x.at(0) - y.at(0)) + (x.at(1) - y.at(1)) * (x.at(1) - y.at(1))))};
}
// Not scale-factor for non_homogenous kernels
template<typename ValueType>
[[nodiscard]] inline auto scale_factor(ValueType /*cell_width*/) const noexcept
{
// return 1 because non homogeneous kernel functions cannot be scaled!!!
return vector_type<ValueType>({ValueType(1.)});
return vector_type<ValueType>{ValueType(1.)};
}
template<typename ValueType>
......@@ -400,7 +399,6 @@ namespace scalfmm::matrix_kernels::laplace
/// The kernel is homogeneous
/// \f$k(ax,ay) = 1/a k(x,y)\f$
///
#ifdef TO_COMPILE
struct grad_ln_2d
{
static constexpr auto homogeneity_tag{homogeneity::homogenous};
......@@ -408,9 +406,9 @@ namespace scalfmm::matrix_kernels::laplace
static constexpr std::size_t km{1};
static constexpr std::size_t kn{2};
template<typename ValueType>
using matrix_type = xt::xtensor_fixed<ValueType, xt::xshape<kn, km>>;
using matrix_type = std::array<ValueType, kn * km>;
template<typename ValueType>
using vector_type = xt::xtensor_fixed<ValueType, xt::xshape<kn>>;
using vector_type = std::array<ValueType, kn>;
template<typename ValueType>
[[nodiscard]] inline constexpr auto mutual_coefficient() const
......@@ -422,17 +420,16 @@ namespace scalfmm::matrix_kernels::laplace
[[nodiscard]] inline auto evaluate(container::point<ValueType, 2> const& x,
container::point<ValueType, 2> const& y) const noexcept
{
ValueType tmp = ValueType(1.0) / ((((xs.at(Is) - ys.at(Is)) * (xs.at(Is) - ys.at(Is))) + ...));
return matrix_type<ValueType>({{(tmp * (ys.at(Is) - xs.at(Is)))}...}); // r2 * (xs - ys);
auto diff = x - y;
ValueType tmp = ValueType(1.0) / (diff.at(0) * diff.at(0) + diff.at(1) * diff.at(1));
return matrix_type<ValueType>{tmp * diff.at(0), tmp * diff.at(1)};
}
template<typename ValueType>
[[nodiscard]] inline auto scale_factor(ValueType cell_width) const noexcept
{
auto tmp{ValueType(2.) / cell_width};
vector_type<ValueType> sf;
std::fill(std::begin(sf), std::end(sf), tmp);
return sf;
return vector_type<ValueType>{tmp, tmp};
}
template<typename ValueType>
......@@ -442,7 +439,6 @@ namespace scalfmm::matrix_kernels::laplace
}
const std::size_t separation_criterion{1};
};
#endif
} // namespace scalfmm::matrix_kernels::laplace
#endif // SCALFMM_MATRIX_KERNELS_LAPLACE_HPP
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment