6#include <ginkgo/extensions/kokkos.hpp>
7#include <ginkgo/ginkgo.hpp>
12using gko::ext::kokkos::map_data;
18#ifdef KOKKOS_ENABLE_SERIAL
19 if (execName ==
"reference-executor") {
20 return gko::ext::kokkos::create_executor(Kokkos::Serial{});
23#ifdef PRECICE_WITH_OPENMP
24 if (execName ==
"omp-executor") {
25 return gko::ext::kokkos::create_executor(Kokkos::OpenMP{});
28#ifdef PRECICE_WITH_CUDA
29 if (execName ==
"cuda-executor") {
30 if (enableUnifiedMemory) {
31 return gko::ext::kokkos::create_executor(Kokkos::Cuda{}, Kokkos::CudaUVMSpace{});
33 return gko::ext::kokkos::create_executor(Kokkos::Cuda{}, Kokkos::CudaSpace{});
37#ifdef PRECICE_WITH_HIP
38 if (execName ==
"hip-executor") {
39 if (enableUnifiedMemory) {
40 return gko::ext::kokkos::create_executor(Kokkos::HIP{}, Kokkos::HIPManagedSpace{});
42 return gko::ext::kokkos::create_executor(Kokkos::HIP{}, Kokkos::HIPSpace{});
51template <
typename MemorySpace,
typename EvalFunctionType>
53 gko::ptr_param<GinkgoMatrix> mtx,
55 gko::ptr_param<GinkgoMatrix> supportPoints,
56 gko::ptr_param<GinkgoMatrix> targetPoints,
60 unsigned int extraDims)
62 auto k_mtx = map_data<MemorySpace>(mtx.get());
63 auto k_supportPoints = map_data<MemorySpace>(supportPoints.get());
64 auto k_targetPoints = map_data<MemorySpace>(targetPoints.get());
71 "create_rbf_system_matrix_row_major",
72 Kokkos::MDRangePolicy<
typename MemorySpace::execution_space, Kokkos::Rank<2>>{{0, 0}, {mtx->get_size()[0], mtx->get_size()[1]}},
73 KOKKOS_LAMBDA(
const int &i,
const int &j) {
81 const bool *deviceActiveAxis = activeAxis.
data();
84 for (
size_t k = 0; k < activeAxis.
size(); ++k) {
85 if (deviceActiveAxis[k]) {
86 double diff = k_supportPoints(j, k) - k_targetPoints(i, k);
90 dist = Kokkos::sqrt(dist);
91 k_mtx(i, j) = f(dist, rbf_params);
96 "create_rbf_system_matrix_col_major",
97 Kokkos::MDRangePolicy<
typename MemorySpace::execution_space, Kokkos::Rank<2>>{{0, 0}, {mtx->get_size()[0], mtx->get_size()[1]}},
98 KOKKOS_LAMBDA(
const int &i,
const int &j) {
103 const bool *deviceActiveAxis = activeAxis.
data();
106 for (
size_t k = 0; k < activeAxis.
size(); ++k) {
107 if (deviceActiveAxis[k]) {
108 double diff = k_supportPoints(k, j) - k_targetPoints(k, i);
112 dist = Kokkos::sqrt(dist);
113 k_mtx(i, j) = f(dist, rbf_params);
118template <
typename EvalFunctionType>
121 gko::ptr_param<GinkgoMatrix> mtx,
123 gko::ptr_param<GinkgoMatrix> supportPoints,
124 gko::ptr_param<GinkgoMatrix> targetPoints,
128 unsigned int extraDims)
130#ifdef KOKKOS_ENABLE_SERIAL
136#ifdef PRECICE_WITH_OPENMP
142#ifdef PRECICE_WITH_CUDA
152#ifdef PRECICE_WITH_HIP
165#define PRECICE_INSTANTIATE_CREATE_RBF_SYSTEM_MATRIX(_function_type) \
166 template void create_rbf_system_matrix<_function_type>(std::shared_ptr<const gko::Executor> exec, \
167 bool unifiedMemory, \
168 gko::ptr_param<GinkgoMatrix> mtx, const std::array<bool, 3> activeAxis, \
169 gko::ptr_param<GinkgoMatrix> supportPoints, gko::ptr_param<GinkgoMatrix> targetPoints, \
170 _function_type f, RadialBasisParameters rbf_params, bool addPolynomial, unsigned int extraDims)
183#undef PRECICE_INSTANTIATE_CREATE_RBF_SYSTEM_MATRIX
185template <
typename MemorySpace>
187 gko::ptr_param<GinkgoMatrix> mtx,
188 gko::ptr_param<const GinkgoMatrix> x,
189 const unsigned int dims)
192 auto k_mtx = map_data<MemorySpace>(mtx.get());
193 auto k_x = map_data<MemorySpace>(x.get());
196 Kokkos::parallel_for(
197 "fill_polynomial_matrix_row_major",
198 Kokkos::MDRangePolicy<
typename MemorySpace::execution_space, Kokkos::Rank<2>>{{0, 0}, {mtx->get_size()[0], mtx->get_size()[1]}},
199 KOKKOS_LAMBDA(
const int &i,
const int &j) {
209 Kokkos::parallel_for(
210 "fill_polynomial_matrix_col_major",
211 Kokkos::MDRangePolicy<
typename MemorySpace::execution_space, Kokkos::Rank<2>>{{0, 0}, {mtx->get_size()[0], mtx->get_size()[1]}},
212 KOKKOS_LAMBDA(
const int &i,
const int &j) {
226 gko::ptr_param<GinkgoMatrix> mtx,
227 gko::ptr_param<const GinkgoMatrix> x,
228 const unsigned int dims)
230#ifdef KOKKOS_ENABLE_SERIAL
236#ifdef PRECICE_WITH_OPENMP
242#ifdef PRECICE_WITH_CUDA
252#ifdef PRECICE_WITH_HIP
#define PRECICE_INSTANTIATE_CREATE_RBF_SYSTEM_MATRIX(_function_type)
#define PRECICE_UNREACHABLE(...)
Wendland radial basis function with compact support.
Wendland radial basis function with compact support.
Wendland radial basis function with compact support.
Wendland radial basis function with compact support.
Wendland radial basis function with compact support.
Radial basis function with compact support.
Radial basis function with global and compact support.
Radial basis function with global support.
Radial basis function with global support.
Radial basis function with global support.
Radial basis function with global support.
void fill_polynomial_matrix(std::shared_ptr< const gko::Executor > exec, bool unifiedMemory, gko::ptr_param< GinkgoMatrix > mtx, gko::ptr_param< const GinkgoMatrix > x, const unsigned int dims)
void fill_polynomial_matrix_impl(std::shared_ptr< const gko::Executor > exec, gko::ptr_param< GinkgoMatrix > mtx, gko::ptr_param< const GinkgoMatrix > x, const unsigned int dims)
void create_rbf_system_matrix(std::shared_ptr< const gko::Executor > exec, bool unifiedMemory, gko::ptr_param< GinkgoMatrix > mtx, const std::array< bool, 3 > activeAxis, gko::ptr_param< GinkgoMatrix > supportPoints, gko::ptr_param< GinkgoMatrix > targetPoints, EvalFunctionType f, ::precice::mapping::RadialBasisParameters rbf_params, bool addPolynomial, unsigned int extraDims)
void create_rbf_system_matrix_impl(std::shared_ptr< const gko::Executor > exec, gko::ptr_param< GinkgoMatrix > mtx, const std::array< bool, 3 > activeAxis, gko::ptr_param< GinkgoMatrix > supportPoints, gko::ptr_param< GinkgoMatrix > targetPoints, EvalFunctionType f, ::precice::mapping::RadialBasisParameters rbf_params, bool addPolynomial, unsigned int extraDims)
contains data mapping from points to meshes.
std::shared_ptr< gko::Executor > create_device_executor(const std::string &execName, bool enableUnifiedMemory)
constexpr T pow_int(const T x)
Computes the power of a given number by an integral exponent given at compile time,...
T dynamic_pointer_cast(T... args)
Wrapper struct that is used to transfer RBF-specific parameters to the GPU.