From a06fdcd443c7a067460e54797fc3ea9736245733 Mon Sep 17 00:00:00 2001 From: Baptiste Legouix Date: Thu, 11 Jul 2024 13:51:32 +0200 Subject: [PATCH 1/6] Left/Right, Xmin/Xmax -> Lower/Upper (#534) --- .../splines/greville_interpolation_points.hpp | 16 +-- .../splines/knots_as_interpolation_points.hpp | 10 +- .../ddc/kernels/splines/spline_builder.hpp | 100 ++++++++--------- .../ddc/kernels/splines/spline_builder_2d.hpp | 40 +++---- .../ddc/kernels/splines/spline_evaluator.hpp | 52 ++++----- .../kernels/splines/spline_evaluator_2d.hpp | 104 +++++++++--------- 6 files changed, 161 insertions(+), 161 deletions(-) diff --git a/include/ddc/kernels/splines/greville_interpolation_points.hpp b/include/ddc/kernels/splines/greville_interpolation_points.hpp index 869e34b50..19f2e9f82 100644 --- a/include/ddc/kernels/splines/greville_interpolation_points.hpp +++ b/include/ddc/kernels/splines/greville_interpolation_points.hpp @@ -20,10 +20,10 @@ namespace ddc { * A class which provides helper functions to initialise the Greville points from a B-Spline definition. * * @tparam BSplines The bspline class relative to which the Greville points will be calculated. - * @tparam BcXmin The (left) boundary condition that will be used to build the splines. - * @tparam BcXmax The (right) boundary condition that will be used to build the splines. + * @tparam BcLower The lower boundary condition that will be used to build the splines. + * @tparam BcUpper The upper boundary condition that will be used to build the splines. */ -template +template class GrevilleInterpolationPoints { using continuous_dimension_type = typename BSplines::continuous_dimension_type; @@ -114,8 +114,8 @@ class GrevilleInterpolationPoints return SamplingImpl(greville_points); } - static constexpr int N_BE_MIN = n_boundary_equations(BcXmin, BSplines::degree()); - static constexpr int N_BE_MAX = n_boundary_equations(BcXmax, BSplines::degree()); + static constexpr int N_BE_MIN = n_boundary_equations(BcLower, BSplines::degree()); + static constexpr int N_BE_MAX = n_boundary_equations(BcUpper, BSplines::degree()); template static constexpr bool is_uniform_mesh_v = U::is_uniform() && ((N_BE_MIN != 0 && N_BE_MAX != 0) || U::is_periodic()); @@ -167,7 +167,7 @@ class GrevilleInterpolationPoints std::vector points_with_bcs(npoints); // Construct Greville-like points at the edge - if constexpr (BcXmin == ddc::BoundCond::GREVILLE) { + if constexpr (BcLower == ddc::BoundCond::GREVILLE) { for (std::size_t i(0); i < BSplines::degree() / 2 + 1; ++i) { points_with_bcs[i] = (BSplines::degree() - i) * ddc::discrete_space().rmin(); @@ -190,7 +190,7 @@ class GrevilleInterpolationPoints } int const n_start - = (BcXmin == ddc::BoundCond::GREVILLE) ? BSplines::degree() / 2 + 1 : 1; + = (BcLower == ddc::BoundCond::GREVILLE) ? BSplines::degree() / 2 + 1 : 1; int const domain_size = n_break_points - 2; ddc::DiscreteElement domain_start(1); ddc::DiscreteDomain const @@ -202,7 +202,7 @@ class GrevilleInterpolationPoints }); // Construct Greville-like points at the edge - if constexpr (BcXmax == ddc::BoundCond::GREVILLE) { + if constexpr (BcUpper == ddc::BoundCond::GREVILLE) { for (std::size_t i(0); i < BSplines::degree() / 2 + 1; ++i) { points_with_bcs[npoints - 1 - i] = (BSplines::degree() - i) * ddc::discrete_space().rmax(); diff --git a/include/ddc/kernels/splines/knots_as_interpolation_points.hpp b/include/ddc/kernels/splines/knots_as_interpolation_points.hpp index fa0e7f8d7..9bd2b2b81 100644 --- a/include/ddc/kernels/splines/knots_as_interpolation_points.hpp +++ b/include/ddc/kernels/splines/knots_as_interpolation_points.hpp @@ -26,14 +26,14 @@ namespace ddc { * with this choice of interpolation points. * * @tparam BSplines The type of the uniform or non-uniform spline basis whose knots are used as interpolation points. - * @tparam BcXmin The lower boundary condition. - * @tparam BcXmin The upper boundary condition. + * @tparam BcLower The lower boundary condition. + * @tparam BcLower The upper boundary condition. */ -template +template class KnotsAsInterpolationPoints { - static_assert(BcXmin != ddc::BoundCond::GREVILLE); - static_assert(BcXmax != ddc::BoundCond::GREVILLE); + static_assert(BcLower != ddc::BoundCond::GREVILLE); + static_assert(BcUpper != ddc::BoundCond::GREVILLE); using continuous_dimension_type = typename BSplines::continuous_dimension_type; diff --git a/include/ddc/kernels/splines/spline_builder.hpp b/include/ddc/kernels/splines/spline_builder.hpp index 7cf19fed3..d5bad9617 100644 --- a/include/ddc/kernels/splines/spline_builder.hpp +++ b/include/ddc/kernels/splines/spline_builder.hpp @@ -32,14 +32,14 @@ enum class SplineSolver { * A class which contains an operator () which can be used to build a spline approximation * of a function. A spline approximation is represented by coefficients stored in a Chunk * of B-splines. The spline is constructed such that it respects the boundary conditions - * BcXmin and BcXmax, and it interpolates the function at the points on the interpolation_mesh + * BcLower and BcUpper, and it interpolates the function at the points on the interpolation_mesh * associated with interpolation_mesh_type. * @tparam ExecSpace The Kokkos execution space on which the spline approximation is performed. * @tparam MemorySpace The Kokkos memory space on which the data (interpolation function and splines coefficients) is stored. * @tparam BSplines The discrete dimension representing the B-splines. * @tparam InterpolationMesh The discrete dimension on which interpolation points are defined. - * @tparam BcXmin The lower boundary condition. - * @tparam BcXmax The upper boundary condition. + * @tparam BcLower The lower boundary condition. + * @tparam BcUpper The upper boundary condition. * @tparam Solver The SplineSolver giving the backend used to perform the spline approximation. * @tparam IDimX A variadic template of all the discrete dimensions forming the full space (InterpolationMesh + batched dimensions). */ @@ -48,17 +48,17 @@ template < class MemorySpace, class BSplines, class InterpolationMesh, - ddc::BoundCond BcXmin, - ddc::BoundCond BcXmax, + ddc::BoundCond BcLower, + ddc::BoundCond BcUpper, SplineSolver Solver, class... IDimX> class SplineBuilder { static_assert( - (BSplines::is_periodic() && (BcXmin == ddc::BoundCond::PERIODIC) - && (BcXmax == ddc::BoundCond::PERIODIC)) - || (!BSplines::is_periodic() && (BcXmin != ddc::BoundCond::PERIODIC) - && (BcXmax != ddc::BoundCond::PERIODIC))); + (BSplines::is_periodic() && (BcLower == ddc::BoundCond::PERIODIC) + && (BcUpper == ddc::BoundCond::PERIODIC)) + || (!BSplines::is_periodic() && (BcLower != ddc::BoundCond::PERIODIC) + && (BcUpper != ddc::BoundCond::PERIODIC))); private: using continuous_dimension_type = typename InterpolationMesh::continuous_dimension_type; @@ -143,16 +143,16 @@ class SplineBuilder static constexpr bool s_odd = BSplines::degree() % 2; /// @brief The number of equations defining the boundary condition at the lower bound. - static constexpr int s_nbc_xmin = n_boundary_equations(BcXmin, BSplines::degree()); + static constexpr int s_nbc_xmin = n_boundary_equations(BcLower, BSplines::degree()); /// @brief The number of equations defining the boundary condition at the upper bound. - static constexpr int s_nbc_xmax = n_boundary_equations(BcXmax, BSplines::degree()); + static constexpr int s_nbc_xmax = n_boundary_equations(BcUpper, BSplines::degree()); /// @brief The boundary condition implemented at the lower bound. - static constexpr ddc::BoundCond s_bc_xmin = BcXmin; + static constexpr ddc::BoundCond s_bc_xmin = BcLower; /// @brief The boundary condition implemented at the upper bound. - static constexpr ddc::BoundCond s_bc_xmax = BcXmax; + static constexpr ddc::BoundCond s_bc_xmax = BcUpper; private: batched_interpolation_domain_type m_batched_interpolation_domain; @@ -194,7 +194,7 @@ class SplineBuilder / ddc::discrete_space().ncells()) { static_assert( - ((BcXmin == BoundCond::PERIODIC) == (BcXmax == BoundCond::PERIODIC)), + ((BcLower == BoundCond::PERIODIC) == (BcUpper == BoundCond::PERIODIC)), "Incompatible boundary conditions"); // Calculate block sizes @@ -415,8 +415,8 @@ template < class MemorySpace, class BSplines, class InterpolationMesh, - ddc::BoundCond BcXmin, - ddc::BoundCond BcXmax, + ddc::BoundCond BcLower, + ddc::BoundCond BcUpper, SplineSolver Solver, class... IDimX> int SplineBuilder< @@ -424,8 +424,8 @@ int SplineBuilder< MemorySpace, BSplines, InterpolationMesh, - BcXmin, - BcXmax, + BcLower, + BcUpper, Solver, IDimX...>::compute_offset(interpolation_domain_type const& interpolation_domain) { @@ -458,8 +458,8 @@ template < class MemorySpace, class BSplines, class InterpolationMesh, - ddc::BoundCond BcXmin, - ddc::BoundCond BcXmax, + ddc::BoundCond BcLower, + ddc::BoundCond BcUpper, SplineSolver Solver, class... IDimX> void SplineBuilder< @@ -467,12 +467,12 @@ void SplineBuilder< MemorySpace, BSplines, InterpolationMesh, - BcXmin, - BcXmax, + BcLower, + BcUpper, Solver, IDimX...>::compute_block_sizes_uniform(int& lower_block_size, int& upper_block_size) const { - switch (BcXmin) { + switch (BcLower) { case ddc::BoundCond::PERIODIC: upper_block_size = (bsplines_type::degree()) / 2; break; @@ -485,7 +485,7 @@ void SplineBuilder< default: throw std::runtime_error("ddc::BoundCond not handled"); } - switch (BcXmax) { + switch (BcUpper) { case ddc::BoundCond::PERIODIC: lower_block_size = (bsplines_type::degree()) / 2; break; @@ -505,8 +505,8 @@ template < class MemorySpace, class BSplines, class InterpolationMesh, - ddc::BoundCond BcXmin, - ddc::BoundCond BcXmax, + ddc::BoundCond BcLower, + ddc::BoundCond BcUpper, SplineSolver Solver, class... IDimX> void SplineBuilder< @@ -514,13 +514,13 @@ void SplineBuilder< MemorySpace, BSplines, InterpolationMesh, - BcXmin, - BcXmax, + BcLower, + BcUpper, Solver, IDimX...>::compute_block_sizes_non_uniform(int& lower_block_size, int& upper_block_size) const { - switch (BcXmin) { + switch (BcLower) { case ddc::BoundCond::PERIODIC: upper_block_size = bsplines_type::degree() - 1; break; @@ -533,7 +533,7 @@ void SplineBuilder< default: throw std::runtime_error("ddc::BoundCond not handled"); } - switch (BcXmax) { + switch (BcUpper) { case ddc::BoundCond::PERIODIC: lower_block_size = bsplines_type::degree() - 1; break; @@ -553,8 +553,8 @@ template < class MemorySpace, class BSplines, class InterpolationMesh, - ddc::BoundCond BcXmin, - ddc::BoundCond BcXmax, + ddc::BoundCond BcLower, + ddc::BoundCond BcUpper, SplineSolver Solver, class... IDimX> void SplineBuilder< @@ -562,8 +562,8 @@ void SplineBuilder< MemorySpace, BSplines, InterpolationMesh, - BcXmin, - BcXmax, + BcLower, + BcUpper, Solver, IDimX...>:: allocate_matrix( @@ -621,8 +621,8 @@ template < class MemorySpace, class BSplines, class InterpolationMesh, - ddc::BoundCond BcXmin, - ddc::BoundCond BcXmax, + ddc::BoundCond BcLower, + ddc::BoundCond BcUpper, SplineSolver Solver, class... IDimX> void SplineBuilder< @@ -630,13 +630,13 @@ void SplineBuilder< MemorySpace, BSplines, InterpolationMesh, - BcXmin, - BcXmax, + BcLower, + BcUpper, Solver, IDimX...>::build_matrix_system() { // Hermite boundary conditions at xmin, if any - if constexpr (BcXmin == ddc::BoundCond::HERMITE) { + if constexpr (BcLower == ddc::BoundCond::HERMITE) { std::array derivs_ptr; ddc::DSpan2D @@ -685,7 +685,7 @@ void SplineBuilder< }); // Hermite boundary conditions at xmax, if any - if constexpr (BcXmax == ddc::BoundCond::HERMITE) { + if constexpr (BcUpper == ddc::BoundCond::HERMITE) { std::array derivs_ptr; std::experimental::mdspan< @@ -723,8 +723,8 @@ template < class MemorySpace, class BSplines, class InterpolationMesh, - ddc::BoundCond BcXmin, - ddc::BoundCond BcXmax, + ddc::BoundCond BcLower, + ddc::BoundCond BcUpper, SplineSolver Solver, class... IDimX> template @@ -733,8 +733,8 @@ void SplineBuilder< MemorySpace, BSplines, InterpolationMesh, - BcXmin, - BcXmax, + BcLower, + BcUpper, Solver, IDimX...>:: operator()( @@ -754,21 +754,21 @@ operator()( assert(vals.template extent() == ddc::discrete_space().nbasis() - s_nbc_xmin - s_nbc_xmax); - assert((BcXmin == ddc::BoundCond::HERMITE) + assert((BcLower == ddc::BoundCond::HERMITE) != (!derivs_xmin.has_value() || derivs_xmin->template extent() == 0)); - assert((BcXmax == ddc::BoundCond::HERMITE) + assert((BcUpper == ddc::BoundCond::HERMITE) != (!derivs_xmax.has_value() || derivs_xmax->template extent() == 0)); - if constexpr (BcXmin == BoundCond::HERMITE) { + if constexpr (BcLower == BoundCond::HERMITE) { assert(ddc::DiscreteElement(derivs_xmin->domain().front()).uid() == 1); } - if constexpr (BcXmax == BoundCond::HERMITE) { + if constexpr (BcUpper == BoundCond::HERMITE) { assert(ddc::DiscreteElement(derivs_xmax->domain().front()).uid() == 1); } // Hermite boundary conditions at xmin, if any // NOTE: For consistency with the linear system, the i-th derivative // provided by the user must be multiplied by dx^i - if constexpr (BcXmin == BoundCond::HERMITE) { + if constexpr (BcLower == BoundCond::HERMITE) { assert(derivs_xmin->template extent() == s_nbc_xmin); auto derivs_xmin_values = *derivs_xmin; auto const dx_proxy = m_dx; @@ -809,7 +809,7 @@ operator()( // NOTE: For consistency with the linear system, the i-th derivative // provided by the user must be multiplied by dx^i auto const& nbasis_proxy = ddc::discrete_space().nbasis(); - if constexpr (BcXmax == BoundCond::HERMITE) { + if constexpr (BcUpper == BoundCond::HERMITE) { assert(derivs_xmax->template extent() == s_nbc_xmax); auto derivs_xmax_values = *derivs_xmax; auto const dx_proxy = m_dx; diff --git a/include/ddc/kernels/splines/spline_builder_2d.hpp b/include/ddc/kernels/splines/spline_builder_2d.hpp index 3ffb2c034..6d2b64fdc 100644 --- a/include/ddc/kernels/splines/spline_builder_2d.hpp +++ b/include/ddc/kernels/splines/spline_builder_2d.hpp @@ -25,10 +25,10 @@ template < class BSpline2, class IDimI1, class IDimI2, - ddc::BoundCond BcXmin1, - ddc::BoundCond BcXmax1, - ddc::BoundCond BcXmin2, - ddc::BoundCond BcXmax2, + ddc::BoundCond BcLower1, + ddc::BoundCond BcUpper1, + ddc::BoundCond BcLower2, + ddc::BoundCond BcUpper2, ddc::SplineSolver Solver, class... IDimX> class SplineBuilder2D @@ -46,8 +46,8 @@ class SplineBuilder2D MemorySpace, BSpline1, IDimI1, - BcXmin1, - BcXmax1, + BcLower1, + BcUpper1, Solver, IDimX...>; @@ -57,8 +57,8 @@ class SplineBuilder2D MemorySpace, BSpline2, IDimI2, - BcXmin2, - BcXmax2, + BcLower2, + BcUpper2, Solver, std::conditional_t, BSpline1, IDimX>...>; @@ -68,8 +68,8 @@ class SplineBuilder2D MemorySpace, BSpline1, IDimI1, - BcXmin1, - BcXmax1, + BcLower1, + BcUpper1, Solver, std::conditional_t< std::is_same_v, @@ -394,10 +394,10 @@ template < class BSpline2, class IDimI1, class IDimI2, - ddc::BoundCond BcXmin1, - ddc::BoundCond BcXmax1, - ddc::BoundCond BcXmin2, - ddc::BoundCond BcXmax2, + ddc::BoundCond BcLower1, + ddc::BoundCond BcUpper1, + ddc::BoundCond BcLower2, + ddc::BoundCond BcUpper2, ddc::SplineSolver Solver, class... IDimX> template @@ -408,10 +408,10 @@ void SplineBuilder2D< BSpline2, IDimI1, IDimI2, - BcXmin1, - BcXmax1, - BcXmin2, - BcXmax2, + BcLower1, + BcUpper1, + BcLower2, + BcUpper2, Solver, IDimX...>:: operator()( @@ -465,7 +465,7 @@ operator()( ddc::KokkosAllocator()); auto spline1_deriv_min = spline1_deriv_min_alloc.span_view(); auto spline1_deriv_min_opt = std::optional(spline1_deriv_min.span_cview()); - if constexpr (BcXmin2 == ddc::BoundCond::HERMITE) { + if constexpr (BcLower2 == ddc::BoundCond::HERMITE) { m_spline_builder_deriv1( spline1_deriv_min, *derivs_min2, @@ -489,7 +489,7 @@ operator()( ddc::KokkosAllocator()); auto spline1_deriv_max = spline1_deriv_max_alloc.span_view(); auto spline1_deriv_max_opt = std::optional(spline1_deriv_max.span_cview()); - if constexpr (BcXmax2 == ddc::BoundCond::HERMITE) { + if constexpr (BcUpper2 == ddc::BoundCond::HERMITE) { m_spline_builder_deriv1( spline1_deriv_max, *derivs_max2, diff --git a/include/ddc/kernels/splines/spline_evaluator.hpp b/include/ddc/kernels/splines/spline_evaluator.hpp index efe57eccf..82d051464 100644 --- a/include/ddc/kernels/splines/spline_evaluator.hpp +++ b/include/ddc/kernels/splines/spline_evaluator.hpp @@ -24,8 +24,8 @@ namespace ddc { * @tparam MemorySpace The Kokkos memory space on which the data (spline coefficients and evaluation) is stored. * @tparam BSplines The discrete dimension representing the B-splines. * @tparam EvaluationMesh The discrete dimension on which evaluation points are defined. - * @tparam LeftExtrapolationRule The lower extrapolation rule type. - * @tparam RightExtrapolationRule The upper extrapolation rule type. + * @tparam LowerExtrapolationRule The lower extrapolation rule type. + * @tparam UpperExtrapolationRule The upper extrapolation rule type. * @tparam IDimX A variadic template of all the discrete dimensions forming the full space (EvaluationMesh + batched dimensions). */ template < @@ -33,8 +33,8 @@ template < class MemorySpace, class BSplines, class EvaluationMesh, - class LeftExtrapolationRule, - class RightExtrapolationRule, + class LowerExtrapolationRule, + class UpperExtrapolationRule, class... IDimX> class SplineEvaluator { @@ -97,63 +97,63 @@ class SplineEvaluator ddc::detail::TypeSeq>>; /// @brief The type of the extrapolation rule at the lower boundary. - using left_extrapolation_rule_type = LeftExtrapolationRule; + using lower_extrapolation_rule_type = LowerExtrapolationRule; /// @brief The type of the extrapolation rule at the upper boundary. - using right_extrapolation_rule_type = RightExtrapolationRule; + using upper_extrapolation_rule_type = UpperExtrapolationRule; private: - LeftExtrapolationRule m_left_extrap_rule; + LowerExtrapolationRule m_lower_extrap_rule; - RightExtrapolationRule m_right_extrap_rule; + UpperExtrapolationRule m_upper_extrap_rule; public: static_assert( - std::is_same_v> == bsplines_type::is_periodic() && std::is_same_v< - RightExtrapolationRule, + UpperExtrapolationRule, typename ddc::PeriodicExtrapolationRule< continuous_dimension_type>> == bsplines_type::is_periodic(), "PeriodicExtrapolationRule has to be used if and only if dimension is periodic"); static_assert( std::is_invocable_r_v< double, - LeftExtrapolationRule, + LowerExtrapolationRule, ddc::Coordinate, ddc::ChunkSpan< double const, spline_domain_type, std::experimental::layout_right, memory_space>>, - "LeftExtrapolationRule::operator() has to be callable with usual arguments."); + "LowerExtrapolationRule::operator() has to be callable with usual arguments."); static_assert( std::is_invocable_r_v< double, - RightExtrapolationRule, + UpperExtrapolationRule, ddc::Coordinate, ddc::ChunkSpan< double const, spline_domain_type, std::experimental::layout_right, memory_space>>, - "RightExtrapolationRule::operator() has to be callable with usual arguments."); + "UpperExtrapolationRule::operator() has to be callable with usual arguments."); /** * @brief Build a SplineEvaluator acting on batched_spline_domain. * - * @param left_extrap_rule The extrapolation rule at the lower boundary. - * @param right_extrap_rule The extrapolation rule at the upper boundary. + * @param lower_extrap_rule The extrapolation rule at the lower boundary. + * @param upper_extrap_rule The extrapolation rule at the upper boundary. * * @see NullExtrapolationRule ConstantExtrapolationRule PeriodicExtrapolationRule */ explicit SplineEvaluator( - LeftExtrapolationRule const& left_extrap_rule, - RightExtrapolationRule const& right_extrap_rule) - : m_left_extrap_rule(left_extrap_rule) - , m_right_extrap_rule(right_extrap_rule) + LowerExtrapolationRule const& lower_extrap_rule, + UpperExtrapolationRule const& upper_extrap_rule) + : m_lower_extrap_rule(lower_extrap_rule) + , m_upper_extrap_rule(upper_extrap_rule) { } @@ -199,9 +199,9 @@ class SplineEvaluator * * @see NullExtrapolationRule ConstantExtrapolationRule PeriodicExtrapolationRule */ - left_extrapolation_rule_type left_extrapolation_rule() const + lower_extrapolation_rule_type lower_extrapolation_rule() const { - return m_left_extrap_rule; + return m_lower_extrap_rule; } /** @@ -213,9 +213,9 @@ class SplineEvaluator * * @see NullExtrapolationRule ConstantExtrapolationRule PeriodicExtrapolationRule */ - right_extrapolation_rule_type right_extrapolation_rule() const + upper_extrapolation_rule_type upper_extrapolation_rule() const { - return m_right_extrap_rule; + return m_upper_extrap_rule; } /** @@ -419,10 +419,10 @@ class SplineEvaluator } } else { if (coord_eval_interest < ddc::discrete_space().rmin()) { - return m_left_extrap_rule(coord_eval_interest, spline_coef); + return m_lower_extrap_rule(coord_eval_interest, spline_coef); } if (coord_eval_interest > ddc::discrete_space().rmax()) { - return m_right_extrap_rule(coord_eval_interest, spline_coef); + return m_upper_extrap_rule(coord_eval_interest, spline_coef); } } return eval_no_bc(coord_eval_interest, spline_coef); diff --git a/include/ddc/kernels/splines/spline_evaluator_2d.hpp b/include/ddc/kernels/splines/spline_evaluator_2d.hpp index 8263cda3b..fff76e40e 100644 --- a/include/ddc/kernels/splines/spline_evaluator_2d.hpp +++ b/include/ddc/kernels/splines/spline_evaluator_2d.hpp @@ -26,10 +26,10 @@ namespace ddc { * @tparam BSplines2 The discrete dimension representing the B-splines along the second dimension of interest. * @tparam EvaluationMesh1 The first discrete dimension on which evaluation points are defined. * @tparam EvaluationMesh2 The second discrete dimension on which evaluation points are defined. - * @tparam LeftExtrapolationRule1 The lower extrapolation rule type along first dimension of interest. - * @tparam RightExtrapolationRule1 The upper extrapolation rule type along first dimension of interest. - * @tparam LeftExtrapolationRule2 The lower extrapolation rule type along second dimension of interest. - * @tparam RightExtrapolationRule2 The upper extrapolation rule type along second dimension of interest. + * @tparam LowerExtrapolationRule1 The lower extrapolation rule type along first dimension of interest. + * @tparam UpperExtrapolationRule1 The upper extrapolation rule type along first dimension of interest. + * @tparam LowerExtrapolationRule2 The lower extrapolation rule type along second dimension of interest. + * @tparam UpperExtrapolationRule2 The upper extrapolation rule type along second dimension of interest. * @tparam IDimX A variadic template of all the discrete dimensions forming the full space (EvaluationMesh1 + EvaluationMesh2 + batched dimensions). */ template < @@ -39,10 +39,10 @@ template < class BSplines2, class EvaluationMesh1, class EvaluationMesh2, - class LeftExtrapolationRule1, - class RightExtrapolationRule1, - class LeftExtrapolationRule2, - class RightExtrapolationRule2, + class LowerExtrapolationRule1, + class UpperExtrapolationRule1, + class LowerExtrapolationRule2, + class UpperExtrapolationRule2, class... IDimX> class SplineEvaluator2D { @@ -125,112 +125,112 @@ class SplineEvaluator2D ddc::detail::TypeSeq>>; /// @brief The type of the extrapolation rule at the lower boundary along the first dimension. - using left_extrapolation_rule_1_type = LeftExtrapolationRule1; + using lower_extrapolation_rule_1_type = LowerExtrapolationRule1; /// @brief The type of the extrapolation rule at the upper boundary along the first dimension. - using right_extrapolation_rule_1_type = RightExtrapolationRule1; + using upper_extrapolation_rule_1_type = UpperExtrapolationRule1; /// @brief The type of the extrapolation rule at the lower boundary along the second dimension. - using left_extrapolation_rule_2_type = LeftExtrapolationRule2; + using lower_extrapolation_rule_2_type = LowerExtrapolationRule2; /// @brief The type of the extrapolation rule at the upper boundary along the second dimension. - using right_extrapolation_rule_2_type = RightExtrapolationRule2; + using upper_extrapolation_rule_2_type = UpperExtrapolationRule2; private: - LeftExtrapolationRule1 m_left_extrap_rule_1; + LowerExtrapolationRule1 m_lower_extrap_rule_1; - RightExtrapolationRule1 m_right_extrap_rule_1; + UpperExtrapolationRule1 m_upper_extrap_rule_1; - LeftExtrapolationRule2 m_left_extrap_rule_2; + LowerExtrapolationRule2 m_lower_extrap_rule_2; - RightExtrapolationRule2 m_right_extrap_rule_2; + UpperExtrapolationRule2 m_upper_extrap_rule_2; public: static_assert( - std::is_same_v> == bsplines_type1::is_periodic() && std::is_same_v< - RightExtrapolationRule1, + UpperExtrapolationRule1, typename ddc::PeriodicExtrapolationRule< continuous_dimension_type1>> == bsplines_type1::is_periodic() && std::is_same_v< - LeftExtrapolationRule2, + LowerExtrapolationRule2, typename ddc::PeriodicExtrapolationRule< continuous_dimension_type2>> == bsplines_type2::is_periodic() && std::is_same_v< - RightExtrapolationRule2, + UpperExtrapolationRule2, typename ddc::PeriodicExtrapolationRule< continuous_dimension_type2>> == bsplines_type2::is_periodic(), "PeriodicExtrapolationRule has to be used if and only if dimension is periodic"); static_assert( std::is_invocable_r_v< double, - LeftExtrapolationRule1, + LowerExtrapolationRule1, ddc::Coordinate, ddc::ChunkSpan< double const, spline_domain_type, std::experimental::layout_right, memory_space>>, - "LeftExtrapolationRule1::operator() has to be callable " + "LowerExtrapolationRule1::operator() has to be callable " "with usual arguments."); static_assert( std::is_invocable_r_v< double, - RightExtrapolationRule1, + UpperExtrapolationRule1, ddc::Coordinate, ddc::ChunkSpan< double const, spline_domain_type, std::experimental::layout_right, memory_space>>, - "RightExtrapolationRule1::operator() has to be callable " + "UpperExtrapolationRule1::operator() has to be callable " "with usual arguments."); static_assert( std::is_invocable_r_v< double, - LeftExtrapolationRule2, + LowerExtrapolationRule2, ddc::Coordinate, ddc::ChunkSpan< double const, spline_domain_type, std::experimental::layout_right, memory_space>>, - "LeftExtrapolationRule2::operator() has to be callable " + "LowerExtrapolationRule2::operator() has to be callable " "with usual arguments."); static_assert( std::is_invocable_r_v< double, - RightExtrapolationRule2, + UpperExtrapolationRule2, ddc::Coordinate, ddc::ChunkSpan< double const, spline_domain_type, std::experimental::layout_right, memory_space>>, - "RightExtrapolationRule2::operator() has to be callable " + "UpperExtrapolationRule2::operator() has to be callable " "with usual arguments."); /** * @brief Build a SplineEvaluator2D acting on batched_spline_domain. * - * @param left_extrap_rule1 The extrapolation rule at the lower boundary along the first dimension. - * @param right_extrap_rule1 The extrapolation rule at the upper boundary along the first dimension. - * @param left_extrap_rule2 The extrapolation rule at the lower boundary along the second dimension. - * @param right_extrap_rule2 The extrapolation rule at the upper boundary along the second dimension. + * @param lower_extrap_rule1 The extrapolation rule at the lower boundary along the first dimension. + * @param upper_extrap_rule1 The extrapolation rule at the upper boundary along the first dimension. + * @param lower_extrap_rule2 The extrapolation rule at the lower boundary along the second dimension. + * @param upper_extrap_rule2 The extrapolation rule at the upper boundary along the second dimension. * * @see NullExtrapolationRule ConstantExtrapolationRule PeriodicExtrapolationRule */ explicit SplineEvaluator2D( - LeftExtrapolationRule1 const& left_extrap_rule1, - RightExtrapolationRule1 const& right_extrap_rule1, - LeftExtrapolationRule2 const& left_extrap_rule2, - RightExtrapolationRule2 const& right_extrap_rule2) - : m_left_extrap_rule_1(left_extrap_rule1) - , m_right_extrap_rule_1(right_extrap_rule1) - , m_left_extrap_rule_2(left_extrap_rule2) - , m_right_extrap_rule_2(right_extrap_rule2) + LowerExtrapolationRule1 const& lower_extrap_rule1, + UpperExtrapolationRule1 const& upper_extrap_rule1, + LowerExtrapolationRule2 const& lower_extrap_rule2, + UpperExtrapolationRule2 const& upper_extrap_rule2) + : m_lower_extrap_rule_1(lower_extrap_rule1) + , m_upper_extrap_rule_1(upper_extrap_rule1) + , m_lower_extrap_rule_2(lower_extrap_rule2) + , m_upper_extrap_rule_2(upper_extrap_rule2) { } @@ -276,9 +276,9 @@ class SplineEvaluator2D * * @see NullExtrapolationRule ConstantExtrapolationRule PeriodicExtrapolationRule */ - left_extrapolation_rule_1_type left_extrapolation_rule_dim_1() const + lower_extrapolation_rule_1_type lower_extrapolation_rule_dim_1() const { - return m_left_extrap_rule_1; + return m_lower_extrap_rule_1; } /** @@ -290,9 +290,9 @@ class SplineEvaluator2D * * @see NullExtrapolationRule ConstantExtrapolationRule PeriodicExtrapolationRule */ - right_extrapolation_rule_1_type right_extrapolation_rule_dim_1() const + upper_extrapolation_rule_1_type upper_extrapolation_rule_dim_1() const { - return m_right_extrap_rule_1; + return m_upper_extrap_rule_1; } /** @@ -304,9 +304,9 @@ class SplineEvaluator2D * * @see NullExtrapolationRule ConstantExtrapolationRule PeriodicExtrapolationRule */ - left_extrapolation_rule_2_type left_extrapolation_rule_dim_2() const + lower_extrapolation_rule_2_type lower_extrapolation_rule_dim_2() const { - return m_left_extrap_rule_2; + return m_lower_extrap_rule_2; } /** @@ -318,9 +318,9 @@ class SplineEvaluator2D * * @see NullExtrapolationRule ConstantExtrapolationRule PeriodicExtrapolationRule */ - right_extrapolation_rule_2_type right_extrapolation_rule_dim_2() const + upper_extrapolation_rule_2_type upper_extrapolation_rule_dim_2() const { - return m_right_extrap_rule_2; + return m_upper_extrap_rule_2; } /** @@ -883,18 +883,18 @@ class SplineEvaluator2D } if constexpr (!bsplines_type1::is_periodic()) { if (ddc::get(coord_eval) < ddc::discrete_space().rmin()) { - return m_left_extrap_rule_1(coord_eval, spline_coef); + return m_lower_extrap_rule_1(coord_eval, spline_coef); } if (ddc::get(coord_eval) > ddc::discrete_space().rmax()) { - return m_right_extrap_rule_1(coord_eval, spline_coef); + return m_upper_extrap_rule_1(coord_eval, spline_coef); } } if constexpr (!bsplines_type2::is_periodic()) { if (ddc::get(coord_eval) < ddc::discrete_space().rmin()) { - return m_left_extrap_rule_2(coord_eval, spline_coef); + return m_lower_extrap_rule_2(coord_eval, spline_coef); } if (ddc::get(coord_eval) > ddc::discrete_space().rmax()) { - return m_right_extrap_rule_2(coord_eval, spline_coef); + return m_upper_extrap_rule_2(coord_eval, spline_coef); } } return eval_no_bc( From 30215e4e2272e657141f13180a0fd93a6152c5f9 Mon Sep 17 00:00:00 2001 From: Baptiste Legouix Date: Fri, 12 Jul 2024 07:56:09 +0200 Subject: [PATCH 2/6] Rename Matrix tests (#538) --- tests/splines/CMakeLists.txt | 10 +- ...{matrix.cpp => splines_linear_problem.cpp} | 211 +++++++++--------- 2 files changed, 114 insertions(+), 107 deletions(-) rename tests/splines/{matrix.cpp => splines_linear_problem.cpp} (65%) diff --git a/tests/splines/CMakeLists.txt b/tests/splines/CMakeLists.txt index de6476f24..f8959ffc2 100644 --- a/tests/splines/CMakeLists.txt +++ b/tests/splines/CMakeLists.txt @@ -33,17 +33,17 @@ target_link_libraries(bsplines_tests DDC::DDC ) -add_executable(matrix_tests +add_executable(splines_linear_problem_tests ../main.cpp - matrix.cpp + splines_linear_problem.cpp ) -target_compile_features(matrix_tests PUBLIC cxx_std_17) -target_link_libraries(matrix_tests +target_compile_features(splines_linear_problem_tests PUBLIC cxx_std_17) +target_link_libraries(splines_linear_problem_tests PUBLIC GTest::gtest DDC::DDC ) -gtest_discover_tests(matrix_tests DISCOVERY_MODE PRE_TEST) +gtest_discover_tests(splines_linear_problem_tests DISCOVERY_MODE PRE_TEST) foreach(DEGREE_X RANGE "${SPLINES_TEST_DEGREE_MIN}" "${SPLINES_TEST_DEGREE_MAX}") foreach(BSPLINES_TYPE "BSPLINES_TYPE_UNIFORM" "BSPLINES_TYPE_NON_UNIFORM") diff --git a/tests/splines/matrix.cpp b/tests/splines/splines_linear_problem.cpp similarity index 65% rename from tests/splines/matrix.cpp rename to tests/splines/splines_linear_problem.cpp index 19f5d6ae0..5b1f139a1 100644 --- a/tests/splines/matrix.cpp +++ b/tests/splines/splines_linear_problem.cpp @@ -80,17 +80,18 @@ namespace DDC_HIP_5_7_ANONYMOUS_NAMESPACE_WORKAROUND(MATRIX_CPP) } void solve_and_validate( - ddc::detail::SplinesLinearProblem & matrix) + ddc::detail::SplinesLinearProblem< + Kokkos::DefaultExecutionSpace> & splines_linear_problem) { - const std::size_t N = matrix.size(); + const std::size_t N = splines_linear_problem.size(); std::vector val_ptr(N * N); ddc::detail::SplinesLinearProblem::MultiRHS val(val_ptr.data(), N, N); - copy_matrix(val, matrix); + copy_matrix(val, splines_linear_problem); - matrix.setup_solver(); + splines_linear_problem.setup_solver(); Kokkos::DualView inv_ptr("inv_ptr", N * N); ddc::detail::SplinesLinearProblem::MultiRHS @@ -98,8 +99,9 @@ namespace DDC_HIP_5_7_ANONYMOUS_NAMESPACE_WORKAROUND(MATRIX_CPP) fill_identity(inv); inv_ptr.modify_host(); inv_ptr.sync_device(); - matrix.solve(ddc::detail::SplinesLinearProblem< - Kokkos::DefaultExecutionSpace>::MultiRHS(inv_ptr.d_view.data(), N, N)); + splines_linear_problem.solve( + ddc::detail::SplinesLinearProblem< + Kokkos::DefaultExecutionSpace>::MultiRHS(inv_ptr.d_view.data(), N, N)); inv_ptr.modify_device(); inv_ptr.sync_host(); @@ -109,7 +111,7 @@ namespace DDC_HIP_5_7_ANONYMOUS_NAMESPACE_WORKAROUND(MATRIX_CPP) fill_identity(inv_tr); inv_tr_ptr.modify_host(); inv_tr_ptr.sync_device(); - matrix + splines_linear_problem .solve(ddc::detail::SplinesLinearProblem:: MultiRHS(inv_tr_ptr.d_view.data(), N, N), true); @@ -122,276 +124,281 @@ namespace DDC_HIP_5_7_ANONYMOUS_NAMESPACE_WORKAROUND(MATRIX_CPP) } // namespace ) -TEST(MatrixSparse, Formatting) +TEST(SplinesLinearProblemSparse, Formatting) { - ddc::detail::SplinesLinearProblemSparse matrix(2); - matrix.set_element(0, 0, 1); - matrix.set_element(0, 1, 2); - matrix.set_element(1, 0, 3); - matrix.set_element(1, 1, 4); + ddc::detail::SplinesLinearProblemSparse splines_linear_problem( + 2); + splines_linear_problem.set_element(0, 0, 1); + splines_linear_problem.set_element(0, 1, 2); + splines_linear_problem.set_element(1, 0, 3); + splines_linear_problem.set_element(1, 1, 4); std::stringstream ss; - ss << matrix; + ss << splines_linear_problem; EXPECT_EQ(ss.str(), " 1.000 2.000\n 3.000 4.000\n"); } -TEST(Matrix, Dense) +TEST(SplinesLinearProblem, Dense) { std::size_t const N = 10; std::size_t const k = 10; - std::unique_ptr> matrix - = std::make_unique< + std::unique_ptr> + splines_linear_problem = std::make_unique< ddc::detail::SplinesLinearProblemDense>(N); // Build a non-symmetric full-rank matrix (without zero) for (std::size_t i(0); i < N; ++i) { - matrix->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); + splines_linear_problem->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); for (std::size_t j(std::max(0, int(i) - int(k))); j < i; ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } for (std::size_t j(i + 1); j < std::min(N, i + k + 1); ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } } - solve_and_validate(*matrix); + solve_and_validate(*splines_linear_problem); } -TEST(Matrix, Band) +TEST(SplinesLinearProblem, Band) { std::size_t const N = 10; std::size_t const k = 3; - std::unique_ptr> matrix - = std::make_unique< + std::unique_ptr> + splines_linear_problem = std::make_unique< ddc::detail::SplinesLinearProblemBand>(N, k, k); // Build a non-symmetric full-rank band matrix for (std::size_t i(0); i < N; ++i) { - matrix->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); + splines_linear_problem->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); for (std::size_t j(std::max(0, int(i) - int(k))); j < i; ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } for (std::size_t j(i + 1); j < std::min(N, i + k + 1); ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } } - solve_and_validate(*matrix); + solve_and_validate(*splines_linear_problem); } -TEST(Matrix, PDSBand) +TEST(SplinesLinearProblem, PDSBand) { std::size_t const N = 10; std::size_t const k = 3; - std::unique_ptr> matrix - = std::make_unique< + std::unique_ptr> + splines_linear_problem = std::make_unique< ddc::detail::SplinesLinearProblemPDSBand>(N, k); // Build a positive-definite symmetric full-rank band matrix for (std::size_t i(0); i < N; ++i) { - matrix->set_element(i, i, 2.0 * k + 1); + splines_linear_problem->set_element(i, i, 2.0 * k + 1); for (std::size_t j(std::max(0, int(i) - int(k))); j < i; ++j) { - matrix->set_element(i, j, -1.0); + splines_linear_problem->set_element(i, j, -1.0); } for (std::size_t j(i + 1); j < std::min(N, i + k + 1); ++j) { - matrix->set_element(i, j, -1.0); + splines_linear_problem->set_element(i, j, -1.0); } } - solve_and_validate(*matrix); + solve_and_validate(*splines_linear_problem); } -TEST(Matrix, PDSTridiag) +TEST(SplinesLinearProblem, PDSTridiag) { std::size_t const N = 10; std::size_t const k = 1; - std::unique_ptr> matrix - = std::make_unique< + std::unique_ptr> + splines_linear_problem = std::make_unique< ddc::detail::SplinesLinearProblemPDSTridiag>(N); // Build a positive-definite symmetric full-rank tridiagonal matrix for (std::size_t i(0); i < N; ++i) { - matrix->set_element(i, i, 2.0 * k + 1); + splines_linear_problem->set_element(i, i, 2.0 * k + 1); for (std::size_t j(std::max(0, int(i) - int(k))); j < i; ++j) { - matrix->set_element(i, j, -1.0); + splines_linear_problem->set_element(i, j, -1.0); } for (std::size_t j(i + 1); j < std::min(N, i + k + 1); ++j) { - matrix->set_element(i, j, -1.0); + splines_linear_problem->set_element(i, j, -1.0); } } - solve_and_validate(*matrix); + solve_and_validate(*splines_linear_problem); } -TEST(Matrix, 2x2Blocks) +TEST(SplinesLinearProblem, 2x2Blocks) { std::size_t const N = 10; std::size_t const k = 10; std::unique_ptr> top_left_block = std::make_unique< ddc::detail::SplinesLinearProblemDense>(7); - std::unique_ptr> matrix - = std::make_unique> + splines_linear_problem = std::make_unique>(N, std::move(top_left_block)); // Build a non-symmetric full-rank matrix (without zero) for (std::size_t i(0); i < N; ++i) { - matrix->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); + splines_linear_problem->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); for (std::size_t j(std::max(0, int(i) - int(k))); j < i; ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } for (std::size_t j(i + 1); j < std::min(N, i + k + 1); ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } } - solve_and_validate(*matrix); + solve_and_validate(*splines_linear_problem); } -TEST(Matrix, 3x3Blocks) +TEST(SplinesLinearProblem, 3x3Blocks) { std::size_t const N = 10; std::size_t const k = 10; std::unique_ptr> center_block = std::make_unique< ddc::detail::SplinesLinearProblemDense>(N - 5); - std::unique_ptr> matrix - = std::make_unique> + splines_linear_problem = std::make_unique>(N, 2, std::move(center_block)); // Build a non-symmetric full-rank matrix (without zero) for (std::size_t i(0); i < N; ++i) { std::cout << i; - matrix->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); + splines_linear_problem->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); for (std::size_t j(std::max(0, int(i) - int(k))); j < i; ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } for (std::size_t j(i + 1); j < std::min(N, i + k + 1); ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } } - solve_and_validate(*matrix); + solve_and_validate(*splines_linear_problem); } -class MatrixSizesFixture : public testing::TestWithParam> +class SplinesLinearProblemSizesFixture + : public testing::TestWithParam> { }; -TEST_P(MatrixSizesFixture, NonSymmetric) +TEST_P(SplinesLinearProblemSizesFixture, NonSymmetric) { auto const [N, k] = GetParam(); - std::unique_ptr> matrix - = ddc::detail::SplinesLinearProblemMaker::make_new_band< + std::unique_ptr> + splines_linear_problem = ddc::detail::SplinesLinearProblemMaker::make_new_band< Kokkos::DefaultExecutionSpace>(N, k, k, false); // Build a non-symmetric full-rank band matrix for (std::size_t i(0); i < N; ++i) { - matrix->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); + splines_linear_problem->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); for (std::size_t j(std::max(0, int(i) - int(k))); j < i; ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } for (std::size_t j(i + 1); j < std::min(N, i + k + 1); ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } } - solve_and_validate(*matrix); + solve_and_validate(*splines_linear_problem); } -TEST_P(MatrixSizesFixture, PositiveDefiniteSymmetric) +TEST_P(SplinesLinearProblemSizesFixture, PositiveDefiniteSymmetric) { auto const [N, k] = GetParam(); - std::unique_ptr> matrix - = ddc::detail::SplinesLinearProblemMaker::make_new_band< + std::unique_ptr> + splines_linear_problem = ddc::detail::SplinesLinearProblemMaker::make_new_band< Kokkos::DefaultExecutionSpace>(N, k, k, true); // Build a positive-definite symmetric full-rank band matrix for (std::size_t i(0); i < N; ++i) { - matrix->set_element(i, i, 2.0 * k + 1); + splines_linear_problem->set_element(i, i, 2.0 * k + 1); for (std::size_t j(std::max(0, int(i) - int(k))); j < i; ++j) { - matrix->set_element(i, j, -1.0); + splines_linear_problem->set_element(i, j, -1.0); } for (std::size_t j(i + 1); j < std::min(N, i + k + 1); ++j) { - matrix->set_element(i, j, -1.0); + splines_linear_problem->set_element(i, j, -1.0); } } - solve_and_validate(*matrix); + solve_and_validate(*splines_linear_problem); } -TEST_P(MatrixSizesFixture, OffsetBanded) +TEST_P(SplinesLinearProblemSizesFixture, OffsetBanded) { auto const [N, k] = GetParam(); - std::unique_ptr> matrix - = ddc::detail::SplinesLinearProblemMaker::make_new_band< + std::unique_ptr> + splines_linear_problem = ddc::detail::SplinesLinearProblemMaker::make_new_band< Kokkos::DefaultExecutionSpace>(N, 0, 2 * k, true); // Build a positive-definite symmetric full-rank band matrix permuted in such a way the band is shifted for (std::size_t i(0); i < N; ++i) { for (std::size_t j(i); j < std::min(N, i + k); ++j) { - matrix->set_element(i, i, -1.0); + splines_linear_problem->set_element(i, i, -1.0); } if (i + k < N) { - matrix->set_element(i, i + k, 2.0 * k + 1); + splines_linear_problem->set_element(i, i + k, 2.0 * k + 1); } for (std::size_t j(i + k + 1); j < std::min(N, i + k + 1); ++j) { - matrix->set_element(i, j, -1.0); + splines_linear_problem->set_element(i, j, -1.0); } } - solve_and_validate(*matrix); + solve_and_validate(*splines_linear_problem); } -TEST_P(MatrixSizesFixture, 2x2Blocks) +TEST_P(SplinesLinearProblemSizesFixture, 2x2Blocks) { auto const [N, k] = GetParam(); - std::unique_ptr> matrix + std::unique_ptr> + splines_linear_problem = ddc::detail::SplinesLinearProblemMaker::make_new_block_matrix_with_band_main_block< Kokkos::DefaultExecutionSpace>(N, k, k, false, 3); // Build a non-symmetric full-rank band matrix for (std::size_t i(0); i < N; ++i) { - matrix->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); + splines_linear_problem->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); for (std::size_t j(std::max(0, int(i) - int(k))); j < i; ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } for (std::size_t j(i + 1); j < std::min(N, i + k + 1); ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } } - solve_and_validate(*matrix); + solve_and_validate(*splines_linear_problem); } -TEST_P(MatrixSizesFixture, 3x3Blocks) +TEST_P(SplinesLinearProblemSizesFixture, 3x3Blocks) { auto const [N, k] = GetParam(); - std::unique_ptr> matrix + std::unique_ptr> + splines_linear_problem = ddc::detail::SplinesLinearProblemMaker::make_new_block_matrix_with_band_main_block< Kokkos::DefaultExecutionSpace>(N, k, k, false, 3, 2); // Build a non-symmetric full-rank band matrix for (std::size_t i(0); i < N; ++i) { - matrix->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); + splines_linear_problem->set_element(i, i, 3. / 4 * ((N + 1) * i + 1)); for (std::size_t j(std::max(0, int(i) - int(k))); j < i; ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } for (std::size_t j(i + 1); j < std::min(N, i + k + 1); ++j) { - matrix->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); + splines_linear_problem->set_element(i, j, -(1. / 4) / k * (N * i + j + 1)); } } - solve_and_validate(*matrix); + solve_and_validate(*splines_linear_problem); } -TEST_P(MatrixSizesFixture, PeriodicBand) +TEST_P(SplinesLinearProblemSizesFixture, PeriodicBand) { auto const [N, k] = GetParam(); // Build a full-rank periodic band matrix permuted in such a way the band is shifted for (std::ptrdiff_t s(-k + k / 2 + 1); s < static_cast(k - k / 2); ++s) { - std::unique_ptr> matrix + std::unique_ptr> + splines_linear_problem = ddc::detail::SplinesLinearProblemMaker::make_new_periodic_band_matrix< Kokkos::DefaultExecutionSpace>( N, @@ -403,41 +410,41 @@ TEST_P(MatrixSizesFixture, PeriodicBand) std::ptrdiff_t diag = ddc::detail:: modulo(static_cast(j - i), static_cast(N)); if (diag == s || diag == N + s) { - matrix->set_element(i, j, 2.0 * k + 1); + splines_linear_problem->set_element(i, j, 2.0 * k + 1); } else if (diag <= s + k || diag >= N + s - k) { - matrix->set_element(i, j, -1.); + splines_linear_problem->set_element(i, j, -1.); } } } - solve_and_validate(*matrix); + solve_and_validate(*splines_linear_problem); } } -TEST_P(MatrixSizesFixture, Sparse) +TEST_P(SplinesLinearProblemSizesFixture, Sparse) { auto const [N, k] = GetParam(); - std::unique_ptr> matrix - = ddc::detail::SplinesLinearProblemMaker::make_new_sparse< + std::unique_ptr> + splines_linear_problem = ddc::detail::SplinesLinearProblemMaker::make_new_sparse< Kokkos::DefaultExecutionSpace>(N); // Build a positive-definite symmetric diagonal-dominant band matrix (stored as a sparse matrix) for (std::size_t i(0); i < N; ++i) { for (std::size_t j(0); j < N; ++j) { if (i == j) { - matrix->set_element(i, j, 3. / 4); + splines_linear_problem->set_element(i, j, 3. / 4); } else if ( std::abs(static_cast(j - i)) <= static_cast(k)) { - matrix->set_element(i, j, -(1. / 4) / k); + splines_linear_problem->set_element(i, j, -(1. / 4) / k); } } } - solve_and_validate(*matrix); + solve_and_validate(*splines_linear_problem); } INSTANTIATE_TEST_SUITE_P( MyGroup, - MatrixSizesFixture, + SplinesLinearProblemSizesFixture, testing::Combine(testing::Values(10, 20), testing::Range(1, 7))); From 0da9a03f43e761848c57b7536f7217730a54a323 Mon Sep 17 00:00:00 2001 From: Baptiste Legouix Date: Mon, 15 Jul 2024 15:25:56 +0200 Subject: [PATCH 3/6] Optimize 3x3 blocks (#518) --- .../ddc/kernels/splines/spline_builder.hpp | 11 ++- .../splines/splines_linear_problem.hpp | 20 ++++ .../splines_linear_problem_3x3_blocks.hpp | 99 ++++++++++++------- tests/splines/splines_linear_problem.cpp | 43 +++++--- 4 files changed, 121 insertions(+), 52 deletions(-) diff --git a/include/ddc/kernels/splines/spline_builder.hpp b/include/ddc/kernels/splines/spline_builder.hpp index d5bad9617..97bd41b01 100644 --- a/include/ddc/kernels/splines/spline_builder.hpp +++ b/include/ddc/kernels/splines/spline_builder.hpp @@ -300,17 +300,18 @@ class SplineBuilder /** * @brief Get the whole domain on which spline coefficients are defined, with the dimension of interest being the leading dimension. * - * This is used internally due to solver limitation and because it may be beneficial to computation performance. + * This is used internally due to solver limitation and because it may be beneficial to computation performance. For LAPACK backend and non-periodic boundary condition, we are using SplinesLinearSolver3x3Blocks which requires upper_block_size additional rows for internal operations. * * @return The (transposed) domain for the spline coefficients. */ batched_spline_tr_domain_type batched_spline_tr_domain() const noexcept { - return batched_spline_tr_domain_type( - batched_spline_domain().restrict(ddc::DiscreteDomain( + return batched_spline_tr_domain_type(ddc::replace_dim_of( + batched_spline_domain(), + ddc::DiscreteDomain( ddc::DiscreteElement(0), ddc::DiscreteVector( - ddc::discrete_space().nbasis())))); + matrix->required_number_of_rhs_rows())))); } public: @@ -846,7 +847,7 @@ operator()( // Create a 2D Kokkos::View to manage spline_tr as a matrix Kokkos::View bcoef_section( spline_tr.data_handle(), - ddc::discrete_space().nbasis(), + static_cast(spline_tr.template extent()), batch_domain().size()); // Compute spline coef matrix->solve(bcoef_section); diff --git a/include/ddc/kernels/splines/splines_linear_problem.hpp b/include/ddc/kernels/splines/splines_linear_problem.hpp index 11bb921d4..6434db556 100644 --- a/include/ddc/kernels/splines/splines_linear_problem.hpp +++ b/include/ddc/kernels/splines/splines_linear_problem.hpp @@ -78,6 +78,26 @@ class SplinesLinearProblem { return m_size; } + + /** + * @brief Get the required number of rows of the multi-rhs view passed to solve(). + * + * Implementations may require a number of rows larger than what `size` returns for optimization purposes. + * + * @return The required number of rows of the multi-rhs view. It is guaranteed to be greater or equal to `size`. + */ + std::size_t required_number_of_rhs_rows() const + { + std::size_t const nrows = impl_required_number_of_rhs_rows(); + assert(nrows >= size()); + return nrows; + } + +private: + virtual std::size_t impl_required_number_of_rhs_rows() const + { + return m_size; + } }; /** diff --git a/include/ddc/kernels/splines/splines_linear_problem_3x3_blocks.hpp b/include/ddc/kernels/splines/splines_linear_problem_3x3_blocks.hpp index 1dc93e99b..cdcbcfa5a 100644 --- a/include/ddc/kernels/splines/splines_linear_problem_3x3_blocks.hpp +++ b/include/ddc/kernels/splines/splines_linear_problem_3x3_blocks.hpp @@ -96,9 +96,10 @@ class SplinesLinearProblem3x3Blocks : public SplinesLinearProblem2x2Blocks | b_top | -- Considered as a - * | b_bottom | | b_bottom | -- single bottom block + * | b_top | | - | + * | b_center | -> | b_center | + * | b_bottom | | b_top | -- Considered as a + * | - | | b_bottom | -- single bottom block * * @param b The multiple right-hand sides. */ @@ -106,34 +107,43 @@ class SplinesLinearProblem3x3Blocks : public SplinesLinearProblem2x2Blockssize(); // size of the center block - // prevent Kokkos::deep_copy(b_top_dst, b_top) to be a deep_copy between overlapping allocations - assert(nq >= m_top_size); - MultiRHS const b_top = Kokkos:: subview(b, std::pair {0, m_top_size}, Kokkos::ALL); - MultiRHS const b_center = Kokkos:: + MultiRHS const b_bottom = Kokkos:: subview(b, - std::pair {m_top_size, m_top_size + nq}, + std::pair {m_top_size + nq, size()}, Kokkos::ALL); - MultiRHS const b_center_dst - = Kokkos::subview(b, std::pair {0, nq}, Kokkos::ALL); MultiRHS const b_top_dst = Kokkos:: - subview(b, std::pair {nq, nq + m_top_size}, Kokkos::ALL); + subview(b, + std::pair {m_top_size + nq, 2 * m_top_size + nq}, + Kokkos::ALL); + MultiRHS const b_bottom_dst = Kokkos:: + subview(b, + std::pair< + std::size_t, + std::size_t> {2 * m_top_size + nq, m_top_size + size()}, + Kokkos::ALL); - MultiRHS const buffer = Kokkos::create_mirror(ExecSpace(), b_center); + if (b_bottom.extent(0) > b_top.extent(0)) { + // Need a buffer to prevent overlapping + MultiRHS const buffer = Kokkos::create_mirror(ExecSpace(), b_bottom); - Kokkos::deep_copy(buffer, b_center); + Kokkos::deep_copy(buffer, b_bottom); + Kokkos::deep_copy(b_bottom_dst, buffer); + } else { + Kokkos::deep_copy(b_bottom_dst, b_bottom); + } Kokkos::deep_copy(b_top_dst, b_top); - Kokkos::deep_copy(b_center_dst, buffer); } /** * @brief Perform row interchanges on multiple right-hand sides to restore its 3-blocks structure. * - * | b_center | | b_top | - * | b_top | -> | b_center | - * | b_bottom | | b_bottom | + * | - | | b_top | + * | b_center | -> | b_center | + * | b_top | | b_bottom | + * | b_bottom | | - | * * @param b The multiple right-hand sides. */ @@ -141,26 +151,34 @@ class SplinesLinearProblem3x3Blocks : public SplinesLinearProblem2x2Blockssize(); // size of the center block - // prevent Kokkos::deep_copy(b_top, b_top_src) to be a deep_copy between overlapping allocations - assert(nq >= m_top_size); - - MultiRHS const b_center_src - = Kokkos::subview(b, std::pair {0, nq}, Kokkos::ALL); - MultiRHS const b_top_src = Kokkos:: - subview(b, std::pair {nq, nq + m_top_size}, Kokkos::ALL); - MultiRHS const b_top = Kokkos:: subview(b, std::pair {0, m_top_size}, Kokkos::ALL); - MultiRHS const b_center = Kokkos:: + MultiRHS const b_bottom = Kokkos:: subview(b, - std::pair {m_top_size, m_top_size + nq}, + std::pair {m_top_size + nq, size()}, Kokkos::ALL); - MultiRHS const buffer = Kokkos::create_mirror(ExecSpace(), b_center); + MultiRHS const b_top_src = Kokkos:: + subview(b, + std::pair {m_top_size + nq, 2 * m_top_size + nq}, + Kokkos::ALL); + MultiRHS const b_bottom_src = Kokkos:: + subview(b, + std::pair< + std::size_t, + std::size_t> {2 * m_top_size + nq, m_top_size + size()}, + Kokkos::ALL); - Kokkos::deep_copy(buffer, b_center_src); Kokkos::deep_copy(b_top, b_top_src); - Kokkos::deep_copy(b_center, buffer); + if (b_bottom.extent(0) > b_top.extent(0)) { + // Need a buffer to prevent overlapping + MultiRHS const buffer = Kokkos::create_mirror(ExecSpace(), b_bottom); + + Kokkos::deep_copy(buffer, b_bottom_src); + Kokkos::deep_copy(b_bottom, buffer); + } else { + Kokkos::deep_copy(b_bottom, b_bottom_src); + } } public: @@ -169,17 +187,32 @@ class SplinesLinearProblem3x3Blocks : public SplinesLinearProblem2x2Blocks::solve(b, transpose); + SplinesLinearProblem2x2Blocks:: + solve(Kokkos:: + subview(b, + std::pair< + std::size_t, + std::size_t> {m_top_size, m_top_size + size()}, + Kokkos::ALL), + transpose); interchange_rows_from_2_to_3_blocks_rhs(b); } + +private: + std::size_t impl_required_number_of_rhs_rows() const override + { + return size() + m_top_size; + } }; } // namespace ddc::detail diff --git a/tests/splines/splines_linear_problem.cpp b/tests/splines/splines_linear_problem.cpp index 5b1f139a1..ee79962de 100644 --- a/tests/splines/splines_linear_problem.cpp +++ b/tests/splines/splines_linear_problem.cpp @@ -21,7 +21,6 @@ namespace DDC_HIP_5_7_ANONYMOUS_NAMESPACE_WORKAROUND(MATRIX_CPP) void fill_identity( ddc::detail::SplinesLinearProblem::MultiRHS mat) { - assert(mat.extent(0) == mat.extent(1)); for (std::size_t i(0); i < mat.extent(0); ++i) { for (std::size_t j(0); j < mat.extent(1); ++j) { mat(i, j) = int(i == j); @@ -93,33 +92,45 @@ namespace DDC_HIP_5_7_ANONYMOUS_NAMESPACE_WORKAROUND(MATRIX_CPP) splines_linear_problem.setup_solver(); - Kokkos::DualView inv_ptr("inv_ptr", N * N); + Kokkos::DualView + inv_ptr("inv_ptr", splines_linear_problem.required_number_of_rhs_rows() * N); ddc::detail::SplinesLinearProblem::MultiRHS - inv(inv_ptr.h_view.data(), N, N); + inv(inv_ptr.h_view.data(), splines_linear_problem.required_number_of_rhs_rows(), N); fill_identity(inv); inv_ptr.modify_host(); inv_ptr.sync_device(); splines_linear_problem.solve( - ddc::detail::SplinesLinearProblem< - Kokkos::DefaultExecutionSpace>::MultiRHS(inv_ptr.d_view.data(), N, N)); + ddc::detail::SplinesLinearProblem::MultiRHS( + inv_ptr.d_view.data(), + splines_linear_problem.required_number_of_rhs_rows(), + N)); inv_ptr.modify_device(); inv_ptr.sync_host(); - Kokkos::DualView inv_tr_ptr("inv_tr_ptr", N * N); + Kokkos::DualView + inv_tr_ptr("inv_tr_ptr", splines_linear_problem.required_number_of_rhs_rows() * N); ddc::detail::SplinesLinearProblem::MultiRHS - inv_tr(inv_tr_ptr.h_view.data(), N, N); + inv_tr(inv_tr_ptr.h_view.data(), + splines_linear_problem.required_number_of_rhs_rows(), + N); fill_identity(inv_tr); inv_tr_ptr.modify_host(); inv_tr_ptr.sync_device(); splines_linear_problem - .solve(ddc::detail::SplinesLinearProblem:: - MultiRHS(inv_tr_ptr.d_view.data(), N, N), + .solve(ddc::detail::SplinesLinearProblem::MultiRHS( + inv_tr_ptr.d_view.data(), + splines_linear_problem.required_number_of_rhs_rows(), + N), true); inv_tr_ptr.modify_device(); inv_tr_ptr.sync_host(); - check_inverse(val, inv); - check_inverse_transpose(val, inv_tr); + check_inverse( + val, + Kokkos::subview(inv, std::pair {0, N}, Kokkos::ALL)); + check_inverse_transpose( + val, + Kokkos::subview(inv_tr, std::pair {0, N}, Kokkos::ALL)); } } // namespace ) @@ -255,12 +266,13 @@ TEST(SplinesLinearProblem, 3x3Blocks) { std::size_t const N = 10; std::size_t const k = 10; + std::size_t const top_size = 2; std::unique_ptr> center_block = std::make_unique< ddc::detail::SplinesLinearProblemDense>(N - 5); std::unique_ptr> splines_linear_problem = std::make_unique>(N, 2, std::move(center_block)); + Kokkos::DefaultExecutionSpace>>(N, top_size, std::move(center_block)); // Build a non-symmetric full-rank matrix (without zero) for (std::size_t i(0); i < N; ++i) { @@ -350,10 +362,11 @@ TEST_P(SplinesLinearProblemSizesFixture, OffsetBanded) TEST_P(SplinesLinearProblemSizesFixture, 2x2Blocks) { auto const [N, k] = GetParam(); + std::size_t const bottom_size = 3; std::unique_ptr> splines_linear_problem = ddc::detail::SplinesLinearProblemMaker::make_new_block_matrix_with_band_main_block< - Kokkos::DefaultExecutionSpace>(N, k, k, false, 3); + Kokkos::DefaultExecutionSpace>(N, k, k, false, bottom_size); // Build a non-symmetric full-rank band matrix for (std::size_t i(0); i < N; ++i) { @@ -372,10 +385,12 @@ TEST_P(SplinesLinearProblemSizesFixture, 2x2Blocks) TEST_P(SplinesLinearProblemSizesFixture, 3x3Blocks) { auto const [N, k] = GetParam(); + std::size_t const bottom_size = 3; + std::size_t const top_size = 2; std::unique_ptr> splines_linear_problem = ddc::detail::SplinesLinearProblemMaker::make_new_block_matrix_with_band_main_block< - Kokkos::DefaultExecutionSpace>(N, k, k, false, 3, 2); + Kokkos::DefaultExecutionSpace>(N, k, k, false, bottom_size, top_size); // Build a non-symmetric full-rank band matrix for (std::size_t i(0); i < N; ++i) { From 1c7389e34c5b21e2715a4e4d63ce80d75039628c Mon Sep 17 00:00:00 2001 From: Baptiste Legouix Date: Tue, 16 Jul 2024 13:40:25 +0200 Subject: [PATCH 4/6] Rename mesh -> discrete_dimension in splines (#539) --- benchmarks/splines.cpp | 5 +- examples/characteristics_advection.cpp | 2 +- .../kernels/splines/bsplines_non_uniform.hpp | 51 ++++++----- .../ddc/kernels/splines/bsplines_uniform.hpp | 48 +++++----- .../splines/greville_interpolation_points.hpp | 10 +-- .../splines/knots_as_interpolation_points.hpp | 6 +- .../ddc/kernels/splines/spline_builder.hpp | 65 +++++++------- .../ddc/kernels/splines/spline_builder_2d.hpp | 46 ++++++---- .../ddc/kernels/splines/spline_evaluator.hpp | 28 +++--- .../kernels/splines/spline_evaluator_2d.hpp | 87 ++++++++++--------- tests/splines/batched_2d_spline_builder.cpp | 2 +- tests/splines/batched_spline_builder.cpp | 2 +- tests/splines/extrapolation_rule.cpp | 4 +- tests/splines/non_periodic_spline_builder.cpp | 2 +- tests/splines/periodic_spline_builder.cpp | 2 +- ...periodic_spline_builder_ordered_points.cpp | 2 +- tests/splines/periodicity_spline_builder.cpp | 4 +- 17 files changed, 200 insertions(+), 166 deletions(-) diff --git a/benchmarks/splines.cpp b/benchmarks/splines.cpp index 63050cc9f..49ef44302 100644 --- a/benchmarks/splines.cpp +++ b/benchmarks/splines.cpp @@ -30,7 +30,7 @@ namespace DDC_HIP_5_7_ANONYMOUS_NAMESPACE_WORKAROUND(SPLINES_CPP) BSplinesX, ddc::BoundCond::PERIODIC, ddc::BoundCond::PERIODIC>; - struct DDimX : GrevillePoints::interpolation_mesh_type + struct DDimX : GrevillePoints::interpolation_discrete_dimension_type { }; @@ -148,7 +148,8 @@ static void characteristics_advection(benchmark::State& state) ddc::parallel_for_each( feet_coords.domain(), KOKKOS_LAMBDA(ddc::DiscreteElement const e) { - feet_coords(e) = ddc::coordinate(ddc::select(e)) - ddc::Coordinate(0.0176429863); + feet_coords(e) = ddc::coordinate(ddc::select(e)) + - ddc::Coordinate(0.0176429863); }); Kokkos::Profiling::popRegion(); Kokkos::Profiling::pushRegion("SplineBuilder"); diff --git a/examples/characteristics_advection.cpp b/examples/characteristics_advection.cpp index f5775d538..6f51ba9ec 100644 --- a/examples/characteristics_advection.cpp +++ b/examples/characteristics_advection.cpp @@ -47,7 +47,7 @@ struct BSplinesX : ddc::UniformBSplines }; using GrevillePoints = ddc:: GrevilleInterpolationPoints; -struct DDimX : GrevillePoints::interpolation_mesh_type +struct DDimX : GrevillePoints::interpolation_discrete_dimension_type { }; //! [X-discretization] diff --git a/include/ddc/kernels/splines/bsplines_non_uniform.hpp b/include/ddc/kernels/splines/bsplines_non_uniform.hpp index 962a77af0..cfa47b6e9 100644 --- a/include/ddc/kernels/splines/bsplines_non_uniform.hpp +++ b/include/ddc/kernels/splines/bsplines_non_uniform.hpp @@ -89,7 +89,7 @@ class NonUniformBSplines : detail::NonUniformBSplinesBase public: /// @brief The type of the knots defining the B-splines. - using knot_mesh_type = NonUniformBsplinesKnots; + using knot_discrete_dimension_type = NonUniformBsplinesKnots; /// @brief The type of the discrete dimension representing the B-splines. using discrete_dimension_type = NonUniformBSplines; @@ -104,8 +104,8 @@ class NonUniformBSplines : detail::NonUniformBSplinesBase using discrete_vector_type = DiscreteVector; private: - ddc::DiscreteDomain m_knot_domain; - ddc::DiscreteDomain m_break_point_domain; + ddc::DiscreteDomain m_knot_domain; + ddc::DiscreteDomain m_break_point_domain; public: Impl() = default; @@ -265,10 +265,11 @@ class NonUniformBSplines : detail::NonUniformBSplinesBase * @param[in] ix DiscreteElement identifying the B-spline. * @return DiscreteElement of the lower bound of the support of the B-spline. */ - KOKKOS_INLINE_FUNCTION ddc::DiscreteElement get_first_support_knot( - discrete_element_type const& ix) const + KOKKOS_INLINE_FUNCTION ddc::DiscreteElement + get_first_support_knot(discrete_element_type const& ix) const { - return ddc::DiscreteElement((ix - discrete_element_type(0)).value()); + return ddc::DiscreteElement( + (ix - discrete_element_type(0)).value()); } /** @brief Returns the coordinate of the last support knot associated to a DiscreteElement identifying a B-spline. @@ -280,10 +281,11 @@ class NonUniformBSplines : detail::NonUniformBSplinesBase * @param[in] ix DiscreteElement identifying the B-spline. * @return DiscreteElement of the upper bound of the support of the B-spline. */ - KOKKOS_INLINE_FUNCTION ddc::DiscreteElement get_last_support_knot( - discrete_element_type const& ix) const + KOKKOS_INLINE_FUNCTION ddc::DiscreteElement + get_last_support_knot(discrete_element_type const& ix) const { - return get_first_support_knot(ix) + ddc::DiscreteVector(degree() + 1); + return get_first_support_knot(ix) + + ddc::DiscreteVector(degree() + 1); } /** @brief Returns the coordinate of the first break point of the domain on which the B-splines are defined. @@ -339,7 +341,8 @@ class NonUniformBSplines : detail::NonUniformBSplinesBase * * @return The discrete domain describing the break points. */ - KOKKOS_INLINE_FUNCTION ddc::DiscreteDomain break_point_domain() const + KOKKOS_INLINE_FUNCTION ddc::DiscreteDomain + break_point_domain() const { return m_break_point_domain; } @@ -379,8 +382,8 @@ class NonUniformBSplines : detail::NonUniformBSplinesBase } private: - KOKKOS_INLINE_FUNCTION discrete_element_type - get_first_bspline_in_cell(ddc::DiscreteElement const& ic) const + KOKKOS_INLINE_FUNCTION discrete_element_type get_first_bspline_in_cell( + ddc::DiscreteElement const& ic) const { return discrete_element_type((ic - m_break_point_domain.front()).value()); } @@ -390,7 +393,7 @@ class NonUniformBSplines : detail::NonUniformBSplinesBase * @param x The point whose location must be determined. * @returns The DiscreteElement describing the knot at the lower bound of the cell of interest. */ - KOKKOS_INLINE_FUNCTION ddc::DiscreteElement find_cell_start( + KOKKOS_INLINE_FUNCTION ddc::DiscreteElement find_cell_start( ddc::Coordinate const& x) const; }; }; @@ -415,13 +418,13 @@ NonUniformBSplines::Impl::Impl( RandomIt const break_begin, RandomIt const break_end) : m_knot_domain( - ddc::DiscreteElement(0), - ddc::DiscreteVector( + ddc::DiscreteElement(0), + ddc::DiscreteVector( (break_end - break_begin) + 2 * degree())) // Create a mesh of knots including the eventual periodic point , m_break_point_domain( - ddc::DiscreteElement(degree()), - ddc::DiscreteVector( + ddc::DiscreteElement(degree()), + ddc::DiscreteVector( (break_end - break_begin))) // Create a mesh of break points { std::vector> knots((break_end - break_begin) + 2 * degree()); @@ -449,7 +452,7 @@ NonUniformBSplines::Impl::Impl( knots[degree() + npoints() - 1 + i] = rmax; } } - ddc::init_discrete_space(knots); + ddc::init_discrete_space(knots); } template @@ -467,7 +470,7 @@ KOKKOS_INLINE_FUNCTION ddc::DiscreteElement NonUniformBSplines:: assert(values.size() == degree() + 1); // 1. Compute cell index 'icell' - ddc::DiscreteElement const icell = find_cell_start(x); + ddc::DiscreteElement const icell = find_cell_start(x); assert(icell >= m_break_point_domain.front()); assert(icell <= m_break_point_domain.back()); @@ -505,7 +508,7 @@ KOKKOS_INLINE_FUNCTION ddc::DiscreteElement NonUniformBSplines:: assert(derivs.size() == degree() + 1); // 1. Compute cell index 'icell' - ddc::DiscreteElement const icell = find_cell_start(x); + ddc::DiscreteElement const icell = find_cell_start(x); assert(icell >= m_break_point_domain.front()); assert(icell <= m_break_point_domain.back()); @@ -581,7 +584,7 @@ KOKKOS_INLINE_FUNCTION ddc::DiscreteElement NonUniformBSplines:: assert(derivs.extent(1) == 1 + n); // 1. Compute cell index 'icell' and x_offset - ddc::DiscreteElement const icell = find_cell_start(x); + ddc::DiscreteElement const icell = find_cell_start(x); assert(icell >= m_break_point_domain.front()); assert(icell <= m_break_point_domain.back()); @@ -672,9 +675,9 @@ KOKKOS_INLINE_FUNCTION ddc::DiscreteElement> NonUn return m_break_point_domain.back() - 1; // Binary search - ddc::DiscreteElement low = m_break_point_domain.front(); - ddc::DiscreteElement high = m_break_point_domain.back(); - ddc::DiscreteElement icell = low + (high - low) / 2; + ddc::DiscreteElement low = m_break_point_domain.front(); + ddc::DiscreteElement high = m_break_point_domain.back(); + ddc::DiscreteElement icell = low + (high - low) / 2; while (x < ddc::coordinate(icell) || x >= ddc::coordinate(icell + 1)) { if (x < ddc::coordinate(icell)) { high = icell; diff --git a/include/ddc/kernels/splines/bsplines_uniform.hpp b/include/ddc/kernels/splines/bsplines_uniform.hpp index 5eef2a9e7..e895c8ba7 100644 --- a/include/ddc/kernels/splines/bsplines_uniform.hpp +++ b/include/ddc/kernels/splines/bsplines_uniform.hpp @@ -89,7 +89,7 @@ class UniformBSplines : detail::UniformBSplinesBase public: /// @brief The type of the knots defining the B-splines. - using knot_mesh_type = UniformBsplinesKnots; + using knot_discrete_dimension_type = UniformBsplinesKnots; /// @brief The type of the discrete dimension representing the B-splines. using discrete_dimension_type = UniformBSplines; @@ -105,8 +105,8 @@ class UniformBSplines : detail::UniformBSplinesBase private: // In the periodic case, they contain the periodic point twice!!! - ddc::DiscreteDomain m_knot_domain; - ddc::DiscreteDomain m_break_point_domain; + ddc::DiscreteDomain m_knot_domain; + ddc::DiscreteDomain m_break_point_domain; public: Impl() = default; @@ -120,16 +120,17 @@ class UniformBSplines : detail::UniformBSplinesBase explicit Impl(ddc::Coordinate rmin, ddc::Coordinate rmax, std::size_t ncells) { assert(ncells > 0); - ddc::DiscreteDomain pre_ghost; - ddc::DiscreteDomain post_ghost; + ddc::DiscreteDomain pre_ghost; + ddc::DiscreteDomain post_ghost; std::tie(m_break_point_domain, m_knot_domain, pre_ghost, post_ghost) - = ddc::init_discrete_space( - knot_mesh_type::template init_ghosted( + = ddc::init_discrete_space( + knot_discrete_dimension_type::template init_ghosted< + knot_discrete_dimension_type>( rmin, rmax, - ddc::DiscreteVector(ncells + 1), - ddc::DiscreteVector(degree()), - ddc::DiscreteVector(degree()))); + ddc::DiscreteVector(ncells + 1), + ddc::DiscreteVector(degree()), + ddc::DiscreteVector(degree()))); } /** @brief Copy-constructs from another Impl with a different Kokkos memory space. @@ -246,10 +247,11 @@ class UniformBSplines : detail::UniformBSplinesBase * @param[in] ix DiscreteElement identifying the B-spline. * @return DiscreteElement of the lower bound of the support of the B-spline. */ - KOKKOS_INLINE_FUNCTION ddc::DiscreteElement get_first_support_knot( - discrete_element_type const& ix) const + KOKKOS_INLINE_FUNCTION ddc::DiscreteElement + get_first_support_knot(discrete_element_type const& ix) const { - return ddc::DiscreteElement((ix - discrete_element_type(0)).value()); + return ddc::DiscreteElement( + (ix - discrete_element_type(0)).value()); } /** @brief Returns the coordinate of the last support knot associated to a DiscreteElement identifying a B-spline. @@ -261,10 +263,11 @@ class UniformBSplines : detail::UniformBSplinesBase * @param[in] ix DiscreteElement identifying the B-spline. * @return DiscreteElement of the upper bound of the support of the B-spline. */ - KOKKOS_INLINE_FUNCTION ddc::DiscreteElement get_last_support_knot( - discrete_element_type const& ix) const + KOKKOS_INLINE_FUNCTION ddc::DiscreteElement + get_last_support_knot(discrete_element_type const& ix) const { - return get_first_support_knot(ix) + ddc::DiscreteVector(degree() + 1); + return get_first_support_knot(ix) + + ddc::DiscreteVector(degree() + 1); } /** @brief Returns the coordinate of the lower bound of the domain on which the B-splines are defined. @@ -320,7 +323,8 @@ class UniformBSplines : detail::UniformBSplinesBase * * @return The discrete domain describing the break points. */ - KOKKOS_INLINE_FUNCTION ddc::DiscreteDomain break_point_domain() const + KOKKOS_INLINE_FUNCTION ddc::DiscreteDomain + break_point_domain() const { return m_break_point_domain; } @@ -351,7 +355,7 @@ class UniformBSplines : detail::UniformBSplinesBase private: KOKKOS_INLINE_FUNCTION double inv_step() const noexcept { - return 1.0 / ddc::step(); + return 1.0 / ddc::step(); } KOKKOS_INLINE_FUNCTION discrete_element_type @@ -427,7 +431,7 @@ KOKKOS_INLINE_FUNCTION ddc::DiscreteElement UniformBSplines:: // 3. Compute derivatives of aforementioned B-splines // Derivatives are normalized, hence they should be divided by dx double xx, temp, saved; - derivs(0) = 1.0 / ddc::step(); + derivs(0) = 1.0 / ddc::step(); for (std::size_t j = 1; j < degree(); ++j) { xx = -offset; saved = 0.0; @@ -596,7 +600,7 @@ UniformBSplines::Impl::integrals( discrete_domain_type const dom_bsplines( full_dom_splines.take_first(discrete_vector_type {nbasis()})); for (auto ix : dom_bsplines) { - int_vals(ix) = ddc::step(); + int_vals(ix) = ddc::step(); } if (int_vals.size() == size()) { discrete_domain_type const dom_bsplines_repeated( @@ -610,7 +614,7 @@ UniformBSplines::Impl::integrals( = full_dom_splines .remove(discrete_vector_type(degree()), discrete_vector_type(degree())); for (auto ix : dom_bspline_entirely_in_domain) { - int_vals(ix) = ddc::step(); + int_vals(ix) = ddc::step(); } std::array edge_vals_ptr; @@ -625,7 +629,7 @@ UniformBSplines::Impl::integrals( for (std::size_t i = 0; i < degree(); ++i) { double const c_eval = ddc::detail::sum(edge_vals, 0, degree() - i); - double const edge_value = ddc::step() * (d_eval - c_eval); + double const edge_value = ddc::step() * (d_eval - c_eval); int_vals(discrete_element_type(i)) = edge_value; int_vals(discrete_element_type(nbasis() - 1 - i)) = edge_value; diff --git a/include/ddc/kernels/splines/greville_interpolation_points.hpp b/include/ddc/kernels/splines/greville_interpolation_points.hpp index 19f2e9f82..fb8d3a9c5 100644 --- a/include/ddc/kernels/splines/greville_interpolation_points.hpp +++ b/include/ddc/kernels/splines/greville_interpolation_points.hpp @@ -117,7 +117,7 @@ class GrevilleInterpolationPoints static constexpr int N_BE_MIN = n_boundary_equations(BcLower, BSplines::degree()); static constexpr int N_BE_MAX = n_boundary_equations(BcUpper, BSplines::degree()); template - static constexpr bool is_uniform_mesh_v + static constexpr bool is_uniform_discrete_dimension_v = U::is_uniform() && ((N_BE_MIN != 0 && N_BE_MAX != 0) || U::is_periodic()); public: @@ -136,7 +136,7 @@ class GrevilleInterpolationPoints class Sampling, typename U = BSplines, std::enable_if_t< - is_uniform_mesh_v, + is_uniform_discrete_dimension_v, bool> = true> // U must be in condition for SFINAE static auto get_sampling() { @@ -154,7 +154,7 @@ class GrevilleInterpolationPoints class Sampling, typename U = BSplines, std::enable_if_t< - !is_uniform_mesh_v, + !is_uniform_discrete_dimension_v, bool> = true> // U must be in condition for SFINAE static auto get_sampling() { @@ -259,8 +259,8 @@ class GrevilleInterpolationPoints * * This is either NonUniformPointSampling or UniformPointSampling. */ - using interpolation_mesh_type = std::conditional_t< - is_uniform_mesh_v, + using interpolation_discrete_dimension_type = std::conditional_t< + is_uniform_discrete_dimension_v, ddc::UniformPointSampling, ddc::NonUniformPointSampling>; diff --git a/include/ddc/kernels/splines/knots_as_interpolation_points.hpp b/include/ddc/kernels/splines/knots_as_interpolation_points.hpp index 9bd2b2b81..eeee23ceb 100644 --- a/include/ddc/kernels/splines/knots_as_interpolation_points.hpp +++ b/include/ddc/kernels/splines/knots_as_interpolation_points.hpp @@ -55,11 +55,11 @@ class KnotsAsInterpolationPoints } else { using SamplingImpl = typename Sampling::template Impl; std::vector knots(ddc::discrete_space().npoints()); - ddc::DiscreteDomain break_point_domain( + ddc::DiscreteDomain break_point_domain( ddc::discrete_space().break_point_domain()); ddc::for_each( break_point_domain, - [&](ddc::DiscreteElement ik) { + [&](ddc::DiscreteElement ik) { knots[ik - break_point_domain.front()] = ddc::coordinate(ik); }); return SamplingImpl(knots); @@ -67,7 +67,7 @@ class KnotsAsInterpolationPoints } /// The DDC type of the sampling for the interpolation points. - using interpolation_mesh_type = std::conditional_t< + using interpolation_discrete_dimension_type = std::conditional_t< is_uniform_bsplines_v, ddc::UniformPointSampling, ddc::NonUniformPointSampling>; diff --git a/include/ddc/kernels/splines/spline_builder.hpp b/include/ddc/kernels/splines/spline_builder.hpp index 97bd41b01..bacfcd4b6 100644 --- a/include/ddc/kernels/splines/spline_builder.hpp +++ b/include/ddc/kernels/splines/spline_builder.hpp @@ -32,22 +32,22 @@ enum class SplineSolver { * A class which contains an operator () which can be used to build a spline approximation * of a function. A spline approximation is represented by coefficients stored in a Chunk * of B-splines. The spline is constructed such that it respects the boundary conditions - * BcLower and BcUpper, and it interpolates the function at the points on the interpolation_mesh - * associated with interpolation_mesh_type. + * BcLower and BcUpper, and it interpolates the function at the points on the interpolation_discrete_dimension + * associated with interpolation_discrete_dimension_type. * @tparam ExecSpace The Kokkos execution space on which the spline approximation is performed. * @tparam MemorySpace The Kokkos memory space on which the data (interpolation function and splines coefficients) is stored. * @tparam BSplines The discrete dimension representing the B-splines. - * @tparam InterpolationMesh The discrete dimension on which interpolation points are defined. + * @tparam InterpolationDDim The discrete dimension on which interpolation points are defined. * @tparam BcLower The lower boundary condition. * @tparam BcUpper The upper boundary condition. * @tparam Solver The SplineSolver giving the backend used to perform the spline approximation. - * @tparam IDimX A variadic template of all the discrete dimensions forming the full space (InterpolationMesh + batched dimensions). + * @tparam IDimX A variadic template of all the discrete dimensions forming the full space (InterpolationDDim + batched dimensions). */ template < class ExecSpace, class MemorySpace, class BSplines, - class InterpolationMesh, + class InterpolationDDim, ddc::BoundCond BcLower, ddc::BoundCond BcUpper, SplineSolver Solver, @@ -61,7 +61,7 @@ class SplineBuilder && (BcUpper != ddc::BoundCond::PERIODIC))); private: - using continuous_dimension_type = typename InterpolationMesh::continuous_dimension_type; + using continuous_dimension_type = typename InterpolationDDim::continuous_dimension_type; public: /// @brief The type of the Kokkos execution space used by this class. @@ -71,7 +71,7 @@ class SplineBuilder using memory_space = MemorySpace; /// @brief The type of the interpolation discrete dimension (discrete dimension of interest) used by this class. - using interpolation_mesh_type = InterpolationMesh; + using interpolation_discrete_dimension_type = InterpolationDDim; /// @brief The discrete dimension representing the B-splines. using bsplines_type = BSplines; @@ -80,7 +80,7 @@ class SplineBuilder using deriv_type = ddc::Deriv; /// @brief The type of the domain for the 1D interpolation mesh used by this class. - using interpolation_domain_type = ddc::DiscreteDomain; + using interpolation_domain_type = ddc::DiscreteDomain; /// @brief The type of the whole domain representing interpolation points. using batched_interpolation_domain_type = ddc::DiscreteDomain; @@ -95,7 +95,7 @@ class SplineBuilder using batch_domain_type = typename ddc::detail::convert_type_seq_to_discrete_domain, - ddc::detail::TypeSeq>>; + ddc::detail::TypeSeq>>; /** * @brief The type of the whole spline domain (cartesian product of 1D spline domain @@ -107,7 +107,7 @@ class SplineBuilder using batched_spline_domain_type = typename ddc::detail::convert_type_seq_to_discrete_domain, - ddc::detail::TypeSeq, + ddc::detail::TypeSeq, ddc::detail::TypeSeq>>; private: @@ -123,7 +123,7 @@ class SplineBuilder ddc::detail::TypeSeq, ddc::type_seq_remove_t< ddc::detail::TypeSeq, - ddc::detail::TypeSeq>>>; + ddc::detail::TypeSeq>>>; public: /** @@ -136,7 +136,7 @@ class SplineBuilder using batched_derivs_domain_type = typename ddc::detail::convert_type_seq_to_discrete_domain, - ddc::detail::TypeSeq, + ddc::detail::TypeSeq, ddc::detail::TypeSeq>>; /// @brief Indicates if the degree of the splines is odd or even. @@ -292,7 +292,7 @@ class SplineBuilder batched_spline_domain_type batched_spline_domain() const noexcept { return ddc::replace_dim_of< - interpolation_mesh_type, + interpolation_discrete_dimension_type, bsplines_type>(batched_interpolation_domain(), spline_domain()); } @@ -324,7 +324,7 @@ class SplineBuilder */ batched_derivs_domain_type batched_derivs_xmin_domain() const noexcept { - return ddc::replace_dim_of( + return ddc::replace_dim_of( batched_interpolation_domain(), ddc::DiscreteDomain( ddc::DiscreteElement(1), @@ -340,7 +340,7 @@ class SplineBuilder */ batched_derivs_domain_type batched_derivs_xmax_domain() const noexcept { - return ddc::replace_dim_of( + return ddc::replace_dim_of( batched_interpolation_domain(), ddc::DiscreteDomain( ddc::DiscreteElement(1), @@ -415,7 +415,7 @@ template < class ExecSpace, class MemorySpace, class BSplines, - class InterpolationMesh, + class InterpolationDDim, ddc::BoundCond BcLower, ddc::BoundCond BcUpper, SplineSolver Solver, @@ -424,7 +424,7 @@ int SplineBuilder< ExecSpace, MemorySpace, BSplines, - InterpolationMesh, + InterpolationDDim, BcLower, BcUpper, Solver, @@ -438,7 +438,8 @@ int SplineBuilder< double, std::experimental::extents> const values(values_ptr.data()); - ddc::DiscreteElement start(interpolation_domain.front()); + ddc::DiscreteElement start( + interpolation_domain.front()); auto jmin = ddc::discrete_space() .eval_basis(values, ddc::coordinate(start + BSplines::degree())); if constexpr (bsplines_type::degree() % 2 == 0) { @@ -458,7 +459,7 @@ template < class ExecSpace, class MemorySpace, class BSplines, - class InterpolationMesh, + class InterpolationDDim, ddc::BoundCond BcLower, ddc::BoundCond BcUpper, SplineSolver Solver, @@ -467,7 +468,7 @@ void SplineBuilder< ExecSpace, MemorySpace, BSplines, - InterpolationMesh, + InterpolationDDim, BcLower, BcUpper, Solver, @@ -505,7 +506,7 @@ template < class ExecSpace, class MemorySpace, class BSplines, - class InterpolationMesh, + class InterpolationDDim, ddc::BoundCond BcLower, ddc::BoundCond BcUpper, SplineSolver Solver, @@ -514,7 +515,7 @@ void SplineBuilder< ExecSpace, MemorySpace, BSplines, - InterpolationMesh, + InterpolationDDim, BcLower, BcUpper, Solver, @@ -553,7 +554,7 @@ template < class ExecSpace, class MemorySpace, class BSplines, - class InterpolationMesh, + class InterpolationDDim, ddc::BoundCond BcLower, ddc::BoundCond BcUpper, SplineSolver Solver, @@ -562,7 +563,7 @@ void SplineBuilder< ExecSpace, MemorySpace, BSplines, - InterpolationMesh, + InterpolationDDim, BcLower, BcUpper, Solver, @@ -621,7 +622,7 @@ template < class ExecSpace, class MemorySpace, class BSplines, - class InterpolationMesh, + class InterpolationDDim, ddc::BoundCond BcLower, ddc::BoundCond BcUpper, SplineSolver Solver, @@ -630,7 +631,7 @@ void SplineBuilder< ExecSpace, MemorySpace, BSplines, - InterpolationMesh, + InterpolationDDim, BcLower, BcUpper, Solver, @@ -676,7 +677,7 @@ void SplineBuilder< ddc::for_each(interpolation_domain(), [&](auto ix) { auto jmin = ddc::discrete_space().eval_basis( values, - ddc::coordinate(ddc::DiscreteElement(ix))); + ddc::coordinate(ddc::DiscreteElement(ix))); for (std::size_t s = 0; s < bsplines_type::degree() + 1; ++s) { int const j = ddc::detail:: modulo(int(jmin.uid() - m_offset + s), @@ -723,7 +724,7 @@ template < class ExecSpace, class MemorySpace, class BSplines, - class InterpolationMesh, + class InterpolationDDim, ddc::BoundCond BcLower, ddc::BoundCond BcUpper, SplineSolver Solver, @@ -733,7 +734,7 @@ void SplineBuilder< ExecSpace, MemorySpace, BSplines, - InterpolationMesh, + InterpolationDDim, BcLower, BcUpper, Solver, @@ -752,7 +753,7 @@ operator()( Layout, memory_space>> const derivs_xmax) const { - assert(vals.template extent() + assert(vals.template extent() == ddc::discrete_space().nbasis() - s_nbc_xmin - s_nbc_xmax); assert((BcLower == ddc::BoundCond::HERMITE) @@ -800,7 +801,9 @@ operator()( spline[ddc::DiscreteDomain( ddc::DiscreteElement(s_nbc_xmin + m_offset), ddc::DiscreteVector(static_cast( - vals.domain().template extent())))] + vals.domain() + .template extent< + interpolation_discrete_dimension_type>())))] .allocation_kokkos_view(), vals.allocation_kokkos_view()); diff --git a/include/ddc/kernels/splines/spline_builder_2d.hpp b/include/ddc/kernels/splines/spline_builder_2d.hpp index 6d2b64fdc..4408ece2a 100644 --- a/include/ddc/kernels/splines/spline_builder_2d.hpp +++ b/include/ddc/kernels/splines/spline_builder_2d.hpp @@ -99,20 +99,25 @@ class SplineBuilder2D using deriv_type2 = typename builder_type2::deriv_type; /// @brief The type of the interpolation mesh in the first dimension. - using interpolation_mesh_type1 = typename builder_type1::interpolation_mesh_type; + using interpolation_discrete_dimension_type1 = + typename builder_type1::interpolation_discrete_dimension_type; /// @brief The type of the interpolation mesh in the second dimension. - using interpolation_mesh_type2 = typename builder_type2::interpolation_mesh_type; + using interpolation_discrete_dimension_type2 = + typename builder_type2::interpolation_discrete_dimension_type; /// @brief The type of the domain for the interpolation mesh in the first dimension. - using interpolation_domain_type1 = typename builder_type1::interpolation_mesh_type; + using interpolation_domain_type1 = + typename builder_type1::interpolation_discrete_dimension_type; /// @brief The type of the domain for the interpolation mesh in the second dimension. - using interpolation_domain_type2 = typename builder_type2::interpolation_mesh_type; + using interpolation_domain_type2 = + typename builder_type2::interpolation_discrete_dimension_type; /// @brief The type of the domain for the interpolation mesh in the 2D dimension. - using interpolation_domain_type - = ddc::DiscreteDomain; + using interpolation_domain_type = ddc::DiscreteDomain< + interpolation_discrete_dimension_type1, + interpolation_discrete_dimension_type2>; /// @brief The type of the whole domain representing interpolation points. using batched_interpolation_domain_type = ddc::DiscreteDomain; @@ -127,7 +132,9 @@ class SplineBuilder2D using batch_domain_type = ddc::detail::convert_type_seq_to_discrete_domain, - ddc::detail::TypeSeq>>; + ddc::detail::TypeSeq< + interpolation_discrete_dimension_type1, + interpolation_discrete_dimension_type2>>>; /** * @brief The type of the whole spline domain (cartesian product of 2D spline domain @@ -139,7 +146,9 @@ class SplineBuilder2D using batched_spline_domain_type = ddc::detail::convert_type_seq_to_discrete_domain, - ddc::detail::TypeSeq, + ddc::detail::TypeSeq< + interpolation_discrete_dimension_type1, + interpolation_discrete_dimension_type2>, ddc::detail::TypeSeq>>; /** @@ -161,7 +170,7 @@ class SplineBuilder2D using batched_derivs_domain_type2 = ddc::detail::convert_type_seq_to_discrete_domain, - ddc::detail::TypeSeq, + ddc::detail::TypeSeq, ddc::detail::TypeSeq>>; /** @@ -174,7 +183,9 @@ class SplineBuilder2D using batched_derivs_domain_type = ddc::detail::convert_type_seq_to_discrete_domain, - ddc::detail::TypeSeq, + ddc::detail::TypeSeq< + interpolation_discrete_dimension_type1, + interpolation_discrete_dimension_type2>, ddc::detail::TypeSeq>>; private: @@ -207,11 +218,12 @@ class SplineBuilder2D batched_interpolation_domain, cols_per_chunk, preconditioner_max_block_size) - , m_spline_builder_deriv1(ddc::replace_dim_of( - m_spline_builder1.batched_interpolation_domain(), - ddc::DiscreteDomain( - ddc::DiscreteElement(1), - ddc::DiscreteVector(bsplines_type2::degree() / 2)))) + , m_spline_builder_deriv1( + ddc::replace_dim_of( + m_spline_builder1.batched_interpolation_domain(), + ddc::DiscreteDomain( + ddc::DiscreteElement(1), + ddc::DiscreteVector(bsplines_type2::degree() / 2)))) , m_spline_builder2( m_spline_builder1.batched_spline_domain(), cols_per_chunk, @@ -304,9 +316,9 @@ class SplineBuilder2D */ batched_spline_domain_type batched_spline_domain() const noexcept { - return ddc::replace_dim_of( + return ddc::replace_dim_of( ddc::replace_dim_of< - interpolation_mesh_type2, + interpolation_discrete_dimension_type2, bsplines_type2>(batched_interpolation_domain(), spline_domain()), spline_domain()); } diff --git a/include/ddc/kernels/splines/spline_evaluator.hpp b/include/ddc/kernels/splines/spline_evaluator.hpp index 82d051464..1716792b2 100644 --- a/include/ddc/kernels/splines/spline_evaluator.hpp +++ b/include/ddc/kernels/splines/spline_evaluator.hpp @@ -23,16 +23,16 @@ namespace ddc { * @tparam ExecSpace The Kokkos execution space on which the spline evaluation is performed. * @tparam MemorySpace The Kokkos memory space on which the data (spline coefficients and evaluation) is stored. * @tparam BSplines The discrete dimension representing the B-splines. - * @tparam EvaluationMesh The discrete dimension on which evaluation points are defined. + * @tparam EvaluationDDim The discrete dimension on which evaluation points are defined. * @tparam LowerExtrapolationRule The lower extrapolation rule type. * @tparam UpperExtrapolationRule The upper extrapolation rule type. - * @tparam IDimX A variadic template of all the discrete dimensions forming the full space (EvaluationMesh + batched dimensions). + * @tparam IDimX A variadic template of all the discrete dimensions forming the full space (EvaluationDDim + batched dimensions). */ template < class ExecSpace, class MemorySpace, class BSplines, - class EvaluationMesh, + class EvaluationDDim, class LowerExtrapolationRule, class UpperExtrapolationRule, class... IDimX> @@ -63,13 +63,13 @@ class SplineEvaluator using memory_space = MemorySpace; /// @brief The type of the evaluation discrete dimension (discrete dimension of interest) used by this class. - using evaluation_mesh_type = EvaluationMesh; + using evaluation_discrete_dimension_type = EvaluationDDim; /// @brief The discrete dimension representing the B-splines. using bsplines_type = BSplines; /// @brief The type of the domain for the 1D evaluation mesh used by this class. - using evaluation_domain_type = ddc::DiscreteDomain; + using evaluation_domain_type = ddc::DiscreteDomain; /// @brief The type of the whole domain representing evaluation points. using batched_evaluation_domain_type = ddc::DiscreteDomain; @@ -84,7 +84,7 @@ class SplineEvaluator using batch_domain_type = typename ddc::detail::convert_type_seq_to_discrete_domain, - ddc::detail::TypeSeq>>; + ddc::detail::TypeSeq>>; /** * @brief The type of the whole spline domain (cartesian product of 1D spline domain @@ -93,7 +93,7 @@ class SplineEvaluator using batched_spline_domain_type = typename ddc::detail::convert_type_seq_to_discrete_domain, - ddc::detail::TypeSeq, + ddc::detail::TypeSeq, ddc::detail::TypeSeq>>; /// @brief The type of the extrapolation rule at the lower boundary. @@ -405,9 +405,10 @@ class SplineEvaluator ddc::ChunkSpan const spline_coef) const { - ddc::Coordinate - coord_eval_interest - = ddc::select(coord_eval); + ddc::Coordinate + coord_eval_interest = ddc::select< + typename evaluation_discrete_dimension_type::continuous_dimension_type>( + coord_eval); if constexpr (bsplines_type::is_periodic()) { if (coord_eval_interest < ddc::discrete_space().rmin() || coord_eval_interest > ddc::discrete_space().rmax()) { @@ -442,9 +443,10 @@ class SplineEvaluator double, std::experimental::extents> const vals(vals_ptr.data()); - ddc::Coordinate - coord_eval_interest - = ddc::select(coord_eval); + ddc::Coordinate + coord_eval_interest = ddc::select< + typename evaluation_discrete_dimension_type::continuous_dimension_type>( + coord_eval); if constexpr (std::is_same_v) { jmin = ddc::discrete_space().eval_basis(vals, coord_eval_interest); } else if constexpr (std::is_same_v) { diff --git a/include/ddc/kernels/splines/spline_evaluator_2d.hpp b/include/ddc/kernels/splines/spline_evaluator_2d.hpp index fff76e40e..cb6c9198e 100644 --- a/include/ddc/kernels/splines/spline_evaluator_2d.hpp +++ b/include/ddc/kernels/splines/spline_evaluator_2d.hpp @@ -24,21 +24,21 @@ namespace ddc { * @tparam MemorySpace The Kokkos memory space on which the data (spline coefficients and evaluation) is stored. * @tparam BSplines1 The discrete dimension representing the B-splines along the first dimension of interest. * @tparam BSplines2 The discrete dimension representing the B-splines along the second dimension of interest. - * @tparam EvaluationMesh1 The first discrete dimension on which evaluation points are defined. - * @tparam EvaluationMesh2 The second discrete dimension on which evaluation points are defined. + * @tparam EvaluationDDim1 The first discrete dimension on which evaluation points are defined. + * @tparam EvaluationDDim2 The second discrete dimension on which evaluation points are defined. * @tparam LowerExtrapolationRule1 The lower extrapolation rule type along first dimension of interest. * @tparam UpperExtrapolationRule1 The upper extrapolation rule type along first dimension of interest. * @tparam LowerExtrapolationRule2 The lower extrapolation rule type along second dimension of interest. * @tparam UpperExtrapolationRule2 The upper extrapolation rule type along second dimension of interest. - * @tparam IDimX A variadic template of all the discrete dimensions forming the full space (EvaluationMesh1 + EvaluationMesh2 + batched dimensions). + * @tparam IDimX A variadic template of all the discrete dimensions forming the full space (EvaluationDDim1 + EvaluationDDim2 + batched dimensions). */ template < class ExecSpace, class MemorySpace, class BSplines1, class BSplines2, - class EvaluationMesh1, - class EvaluationMesh2, + class EvaluationDDim1, + class EvaluationDDim2, class LowerExtrapolationRule1, class UpperExtrapolationRule1, class LowerExtrapolationRule2, @@ -72,10 +72,10 @@ class SplineEvaluator2D using memory_space = MemorySpace; /// @brief The type of the first discrete dimension of interest used by this class. - using evaluation_mesh_type1 = EvaluationMesh1; + using evaluation_discrete_dimension_type1 = EvaluationDDim1; /// @brief The type of the second discrete dimension of interest used by this class. - using evaluation_mesh_type2 = EvaluationMesh2; + using evaluation_discrete_dimension_type2 = EvaluationDDim2; /// @brief The discrete dimension representing the B-splines along first dimension. using bsplines_type1 = BSplines1; @@ -84,14 +84,15 @@ class SplineEvaluator2D using bsplines_type2 = BSplines2; /// @brief The type of the domain for the 1D evaluation mesh along first dimension used by this class. - using evaluation_domain_type1 = ddc::DiscreteDomain; + using evaluation_domain_type1 = ddc::DiscreteDomain; /// @brief The type of the domain for the 1D evaluation mesh along second dimension used by this class. - using evaluation_domain_type2 = ddc::DiscreteDomain; + using evaluation_domain_type2 = ddc::DiscreteDomain; /// @brief The type of the domain for the 2D evaluation mesh used by this class. - using evaluation_domain_type - = ddc::DiscreteDomain; + using evaluation_domain_type = ddc::DiscreteDomain< + evaluation_discrete_dimension_type1, + evaluation_discrete_dimension_type2>; /// @brief The type of the whole domain representing evaluation points. using batched_evaluation_domain_type = ddc::DiscreteDomain; @@ -112,7 +113,9 @@ class SplineEvaluator2D using batch_domain_type = typename ddc::detail::convert_type_seq_to_discrete_domain, - ddc::detail::TypeSeq>>; + ddc::detail::TypeSeq< + evaluation_discrete_dimension_type1, + evaluation_discrete_dimension_type2>>>; /** * @brief The type of the whole spline domain (cartesian product of 2D spline domain @@ -121,7 +124,9 @@ class SplineEvaluator2D using batched_spline_domain_type = typename ddc::detail::convert_type_seq_to_discrete_domain, - ddc::detail::TypeSeq, + ddc::detail::TypeSeq< + evaluation_discrete_dimension_type1, + evaluation_discrete_dimension_type2>, ddc::detail::TypeSeq>>; /// @brief The type of the extrapolation rule at the lower boundary along the first dimension. @@ -477,15 +482,17 @@ class SplineEvaluator2D static_assert( std::is_same_v< InterestDim, - typename evaluation_mesh_type1:: - continuous_dimension_type> || std::is_same_v); + typename evaluation_discrete_dimension_type1:: + continuous_dimension_type> || std::is_same_v); if constexpr (std::is_same_v< InterestDim, - typename evaluation_mesh_type1::continuous_dimension_type>) { + typename evaluation_discrete_dimension_type1:: + continuous_dimension_type>) { return deriv_dim_1(coord_eval, spline_coef); } else if constexpr (std::is_same_v< InterestDim, - typename evaluation_mesh_type2::continuous_dimension_type>) { + typename evaluation_discrete_dimension_type2:: + continuous_dimension_type>) { return deriv_dim_2(coord_eval, spline_coef); } } @@ -515,12 +522,12 @@ class SplineEvaluator2D static_assert( (std::is_same_v< InterestDim1, - typename evaluation_mesh_type1:: - continuous_dimension_type> && std::is_same_v) + typename evaluation_discrete_dimension_type1:: + continuous_dimension_type> && std::is_same_v) || (std::is_same_v< InterestDim2, - typename evaluation_mesh_type1:: - continuous_dimension_type> && std::is_same_v)); + typename evaluation_discrete_dimension_type1:: + continuous_dimension_type> && std::is_same_v)); return deriv_1_and_2(coord_eval, spline_coef); } @@ -713,15 +720,17 @@ class SplineEvaluator2D static_assert( std::is_same_v< InterestDim, - typename evaluation_mesh_type1:: - continuous_dimension_type> || std::is_same_v); + typename evaluation_discrete_dimension_type1:: + continuous_dimension_type> || std::is_same_v); if constexpr (std::is_same_v< InterestDim, - typename evaluation_mesh_type1::continuous_dimension_type>) { + typename evaluation_discrete_dimension_type1:: + continuous_dimension_type>) { return deriv_dim_1(spline_eval, coords_eval, spline_coef); } else if constexpr (std::is_same_v< InterestDim, - typename evaluation_mesh_type2::continuous_dimension_type>) { + typename evaluation_discrete_dimension_type2:: + continuous_dimension_type>) { return deriv_dim_2(spline_eval, coords_eval, spline_coef); } } @@ -770,12 +779,12 @@ class SplineEvaluator2D static_assert( (std::is_same_v< InterestDim1, - typename evaluation_mesh_type1:: - continuous_dimension_type> && std::is_same_v) + typename evaluation_discrete_dimension_type1:: + continuous_dimension_type> && std::is_same_v) || (std::is_same_v< InterestDim2, - typename evaluation_mesh_type1:: - continuous_dimension_type> && std::is_same_v)); + typename evaluation_discrete_dimension_type1:: + continuous_dimension_type> && std::is_same_v)); return deriv_1_and_2(spline_eval, coords_eval, spline_coef); } @@ -857,8 +866,8 @@ class SplineEvaluator2D ddc::ChunkSpan const spline_coef) const { - using Dim1 = typename evaluation_mesh_type1::continuous_dimension_type; - using Dim2 = typename evaluation_mesh_type2::continuous_dimension_type; + using Dim1 = typename evaluation_discrete_dimension_type1::continuous_dimension_type; + using Dim2 = typename evaluation_discrete_dimension_type2::continuous_dimension_type; if constexpr (bsplines_type1::is_periodic()) { if (ddc::get(coord_eval) < ddc::discrete_space().rmin() || ddc::get(coord_eval) > ddc::discrete_space().rmax()) { @@ -899,8 +908,8 @@ class SplineEvaluator2D } return eval_no_bc( ddc::Coordinate< - typename evaluation_mesh_type1::continuous_dimension_type, - typename evaluation_mesh_type2::continuous_dimension_type>( + typename evaluation_discrete_dimension_type1::continuous_dimension_type, + typename evaluation_discrete_dimension_type2::continuous_dimension_type>( ddc::get(coord_eval), ddc::get(coord_eval)), spline_coef); @@ -943,13 +952,13 @@ class SplineEvaluator2D double, std::experimental::extents> const vals2(vals2_ptr.data()); - ddc::Coordinate - coord_eval_interest1 - = ddc::select( + ddc::Coordinate + coord_eval_interest1 = ddc::select< + typename evaluation_discrete_dimension_type1::continuous_dimension_type>( coord_eval); - ddc::Coordinate - coord_eval_interest2 - = ddc::select( + ddc::Coordinate + coord_eval_interest2 = ddc::select< + typename evaluation_discrete_dimension_type2::continuous_dimension_type>( coord_eval); if constexpr (std::is_same_v) { diff --git a/tests/splines/batched_2d_spline_builder.cpp b/tests/splines/batched_2d_spline_builder.cpp index 0e4c20c0e..dd4fa3817 100644 --- a/tests/splines/batched_2d_spline_builder.cpp +++ b/tests/splines/batched_2d_spline_builder.cpp @@ -96,7 +96,7 @@ template struct IDim_ : std::conditional_t< B, - typename GrevillePoints>::interpolation_mesh_type, + typename GrevillePoints>::interpolation_discrete_dimension_type, ddc::UniformPointSampling> { }; diff --git a/tests/splines/batched_spline_builder.cpp b/tests/splines/batched_spline_builder.cpp index 41d968a60..4404f4ddc 100644 --- a/tests/splines/batched_spline_builder.cpp +++ b/tests/splines/batched_spline_builder.cpp @@ -95,7 +95,7 @@ template struct IDim_ : std::conditional_t< B, - typename GrevillePoints>::interpolation_mesh_type, + typename GrevillePoints>::interpolation_discrete_dimension_type, ddc::UniformPointSampling> { }; diff --git a/tests/splines/extrapolation_rule.cpp b/tests/splines/extrapolation_rule.cpp index 57455305a..db5f117de 100644 --- a/tests/splines/extrapolation_rule.cpp +++ b/tests/splines/extrapolation_rule.cpp @@ -78,7 +78,7 @@ template struct IDim : std::conditional_t< std::is_same_v || std::is_same_v, - typename GrevillePoints>::interpolation_mesh_type, + typename GrevillePoints>::interpolation_discrete_dimension_type, ddc::UniformPointSampling> { }; @@ -93,7 +93,7 @@ template struct IDim : std::conditional_t< std::is_same_v || std::is_same_v, - typename GrevillePoints>::interpolation_mesh_type, + typename GrevillePoints>::interpolation_discrete_dimension_type, ddc::NonUniformPointSampling> { }; diff --git a/tests/splines/non_periodic_spline_builder.cpp b/tests/splines/non_periodic_spline_builder.cpp index 9f656d1e0..5af5c6d2b 100644 --- a/tests/splines/non_periodic_spline_builder.cpp +++ b/tests/splines/non_periodic_spline_builder.cpp @@ -50,7 +50,7 @@ struct BSplinesX : ddc::NonUniformBSplines using GrevillePoints = ddc::GrevilleInterpolationPoints; -struct IDimX : GrevillePoints::interpolation_mesh_type +struct IDimX : GrevillePoints::interpolation_discrete_dimension_type { }; diff --git a/tests/splines/periodic_spline_builder.cpp b/tests/splines/periodic_spline_builder.cpp index 794b8f562..78c4d5f70 100644 --- a/tests/splines/periodic_spline_builder.cpp +++ b/tests/splines/periodic_spline_builder.cpp @@ -39,7 +39,7 @@ struct BSplinesX : ddc::NonUniformBSplines using GrevillePoints = ddc:: GrevilleInterpolationPoints; -struct IDimX : GrevillePoints::interpolation_mesh_type +struct IDimX : GrevillePoints::interpolation_discrete_dimension_type { }; diff --git a/tests/splines/periodic_spline_builder_ordered_points.cpp b/tests/splines/periodic_spline_builder_ordered_points.cpp index b69311f58..05a9fec7b 100644 --- a/tests/splines/periodic_spline_builder_ordered_points.cpp +++ b/tests/splines/periodic_spline_builder_ordered_points.cpp @@ -29,7 +29,7 @@ struct BSplinesX : ddc::NonUniformBSplines using GrevillePoints = ddc:: GrevilleInterpolationPoints; -struct IDimX : GrevillePoints::interpolation_mesh_type +struct IDimX : GrevillePoints::interpolation_discrete_dimension_type { }; diff --git a/tests/splines/periodicity_spline_builder.cpp b/tests/splines/periodicity_spline_builder.cpp index b21093b0b..8898861ad 100644 --- a/tests/splines/periodicity_spline_builder.cpp +++ b/tests/splines/periodicity_spline_builder.cpp @@ -38,7 +38,7 @@ struct BSplines : ddc::UniformBSplines // Gives discrete dimension. In the dimension of interest, it is deduced from the BSplines type. In the other dimensions, it has to be newly defined. In practice both types coincide in the test, but it may not be the case. template -struct IDim : GrevillePoints>::interpolation_mesh_type +struct IDim : GrevillePoints>::interpolation_discrete_dimension_type { }; @@ -49,7 +49,7 @@ struct BSplines : ddc::NonUniformBSplines }; template -struct IDim : GrevillePoints>::interpolation_mesh_type +struct IDim : GrevillePoints>::interpolation_discrete_dimension_type { }; From 78c60831e279983aaa3b69316c7730b0366d983f Mon Sep 17 00:00:00 2001 From: Thomas Padioleau Date: Tue, 16 Jul 2024 13:47:59 +0200 Subject: [PATCH 5/6] Update formatting check (#546) --- .github/workflows/tests.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml index c15980bf3..bc9e0a0d7 100644 --- a/.github/workflows/tests.yml +++ b/.github/workflows/tests.yml @@ -16,7 +16,7 @@ jobs: - uses: actions/checkout@v4 - uses: DoozyX/clang-format-lint-action@v0.17 with: - source: 'include/ddc/ tests/ examples/' + source: 'benchmarks/ examples/ include/ddc/ tests/' exclude: '' extensions: 'hpp,cpp' clangFormatVersion: 12 From 2e705cb683681e5889dd40f68707c4c089f9e5fc Mon Sep 17 00:00:00 2001 From: Baptiste Legouix Date: Tue, 16 Jul 2024 14:58:40 +0200 Subject: [PATCH 6/6] Clean splines continuous dimension (#540) --- .../ddc/kernels/splines/spline_builder.hpp | 6 +-- .../ddc/kernels/splines/spline_builder_2d.hpp | 28 ++++++-------- .../ddc/kernels/splines/spline_evaluator.hpp | 17 ++++----- .../kernels/splines/spline_evaluator_2d.hpp | 38 +++++++++---------- 4 files changed, 39 insertions(+), 50 deletions(-) diff --git a/include/ddc/kernels/splines/spline_builder.hpp b/include/ddc/kernels/splines/spline_builder.hpp index bacfcd4b6..4e2e1fb8e 100644 --- a/include/ddc/kernels/splines/spline_builder.hpp +++ b/include/ddc/kernels/splines/spline_builder.hpp @@ -60,9 +60,6 @@ class SplineBuilder || (!BSplines::is_periodic() && (BcLower != ddc::BoundCond::PERIODIC) && (BcUpper != ddc::BoundCond::PERIODIC))); -private: - using continuous_dimension_type = typename InterpolationDDim::continuous_dimension_type; - public: /// @brief The type of the Kokkos execution space used by this class. using exec_space = ExecSpace; @@ -70,6 +67,9 @@ class SplineBuilder /// @brief The type of the Kokkos memory space used by this class. using memory_space = MemorySpace; + /// @brief The type of the interpolation continuous dimension (continuous dimension of interest) used by this class. + using continuous_dimension_type = typename InterpolationDDim::continuous_dimension_type; + /// @brief The type of the interpolation discrete dimension (discrete dimension of interest) used by this class. using interpolation_discrete_dimension_type = InterpolationDDim; diff --git a/include/ddc/kernels/splines/spline_builder_2d.hpp b/include/ddc/kernels/splines/spline_builder_2d.hpp index 4408ece2a..d9a9c0964 100644 --- a/include/ddc/kernels/splines/spline_builder_2d.hpp +++ b/include/ddc/kernels/splines/spline_builder_2d.hpp @@ -76,16 +76,20 @@ class SplineBuilder2D typename builder_type2::deriv_type, IDimX>...>; -private: - /// @brief The tag of the dimension of the first 1D SplineBuilder. - using continuous_dimension_type1 = - typename builder_type1::bsplines_type::continuous_dimension_type; + /// @brief The type of the first interpolation continuous dimension. + using continuous_dimension_type1 = typename builder_type1::continuous_dimension_type; - /// @brief The tag of the dimension of the second 1D SplineBuilder. - using continuous_dimension_type2 = - typename builder_type2::bsplines_type::continuous_dimension_type; + /// @brief The type of the second interpolation continuous dimension. + using continuous_dimension_type2 = typename builder_type2::continuous_dimension_type; + + /// @brief The type of the first interpolation discrete dimension. + using interpolation_discrete_dimension_type1 = + typename builder_type1::interpolation_discrete_dimension_type; + + /// @brief The type of the second interpolation discrete dimension. + using interpolation_discrete_dimension_type2 = + typename builder_type2::interpolation_discrete_dimension_type; -public: /// @brief The type of the B-splines in the first dimension. using bsplines_type1 = typename builder_type1::bsplines_type; @@ -98,14 +102,6 @@ class SplineBuilder2D /// @brief The type of the Deriv domain on boundaries in the second dimension. using deriv_type2 = typename builder_type2::deriv_type; - /// @brief The type of the interpolation mesh in the first dimension. - using interpolation_discrete_dimension_type1 = - typename builder_type1::interpolation_discrete_dimension_type; - - /// @brief The type of the interpolation mesh in the second dimension. - using interpolation_discrete_dimension_type2 = - typename builder_type2::interpolation_discrete_dimension_type; - /// @brief The type of the domain for the interpolation mesh in the first dimension. using interpolation_domain_type1 = typename builder_type1::interpolation_discrete_dimension_type; diff --git a/include/ddc/kernels/splines/spline_evaluator.hpp b/include/ddc/kernels/splines/spline_evaluator.hpp index 1716792b2..144e884f0 100644 --- a/include/ddc/kernels/splines/spline_evaluator.hpp +++ b/include/ddc/kernels/splines/spline_evaluator.hpp @@ -53,8 +53,6 @@ class SplineEvaluator { }; - using continuous_dimension_type = typename BSplines::continuous_dimension_type; - public: /// @brief The type of the Kokkos execution space used by this class. using exec_space = ExecSpace; @@ -62,6 +60,9 @@ class SplineEvaluator /// @brief The type of the Kokkos memory space used by this class. using memory_space = MemorySpace; + /// @brief The type of the evaluation continuous dimension (continuous dimension of interest) used by this class. + using continuous_dimension_type = typename BSplines::continuous_dimension_type; + /// @brief The type of the evaluation discrete dimension (discrete dimension of interest) used by this class. using evaluation_discrete_dimension_type = EvaluationDDim; @@ -405,10 +406,8 @@ class SplineEvaluator ddc::ChunkSpan const spline_coef) const { - ddc::Coordinate - coord_eval_interest = ddc::select< - typename evaluation_discrete_dimension_type::continuous_dimension_type>( - coord_eval); + ddc::Coordinate coord_eval_interest + = ddc::select(coord_eval); if constexpr (bsplines_type::is_periodic()) { if (coord_eval_interest < ddc::discrete_space().rmin() || coord_eval_interest > ddc::discrete_space().rmax()) { @@ -443,10 +442,8 @@ class SplineEvaluator double, std::experimental::extents> const vals(vals_ptr.data()); - ddc::Coordinate - coord_eval_interest = ddc::select< - typename evaluation_discrete_dimension_type::continuous_dimension_type>( - coord_eval); + ddc::Coordinate coord_eval_interest + = ddc::select(coord_eval); if constexpr (std::is_same_v) { jmin = ddc::discrete_space().eval_basis(vals, coord_eval_interest); } else if constexpr (std::is_same_v) { diff --git a/include/ddc/kernels/splines/spline_evaluator_2d.hpp b/include/ddc/kernels/splines/spline_evaluator_2d.hpp index cb6c9198e..818f6685f 100644 --- a/include/ddc/kernels/splines/spline_evaluator_2d.hpp +++ b/include/ddc/kernels/splines/spline_evaluator_2d.hpp @@ -61,10 +61,13 @@ class SplineEvaluator2D { }; +public: + /// @brief The type of the first evaluation continuous dimension used by this class. using continuous_dimension_type1 = typename BSplines1::continuous_dimension_type; + + /// @brief The type of the second evaluation continuous dimension used by this class. using continuous_dimension_type2 = typename BSplines2::continuous_dimension_type; -public: /// @brief The type of the Kokkos execution space used by this class. using exec_space = ExecSpace; @@ -482,8 +485,7 @@ class SplineEvaluator2D static_assert( std::is_same_v< InterestDim, - typename evaluation_discrete_dimension_type1:: - continuous_dimension_type> || std::is_same_v); + continuous_dimension_type1> || std::is_same_v); if constexpr (std::is_same_v< InterestDim, typename evaluation_discrete_dimension_type1:: @@ -523,11 +525,11 @@ class SplineEvaluator2D (std::is_same_v< InterestDim1, typename evaluation_discrete_dimension_type1:: - continuous_dimension_type> && std::is_same_v) + continuous_dimension_type> && std::is_same_v) || (std::is_same_v< InterestDim2, typename evaluation_discrete_dimension_type1:: - continuous_dimension_type> && std::is_same_v)); + continuous_dimension_type> && std::is_same_v)); return deriv_1_and_2(coord_eval, spline_coef); } @@ -721,7 +723,7 @@ class SplineEvaluator2D std::is_same_v< InterestDim, typename evaluation_discrete_dimension_type1:: - continuous_dimension_type> || std::is_same_v); + continuous_dimension_type> || std::is_same_v); if constexpr (std::is_same_v< InterestDim, typename evaluation_discrete_dimension_type1:: @@ -780,11 +782,11 @@ class SplineEvaluator2D (std::is_same_v< InterestDim1, typename evaluation_discrete_dimension_type1:: - continuous_dimension_type> && std::is_same_v) + continuous_dimension_type> && std::is_same_v) || (std::is_same_v< InterestDim2, typename evaluation_discrete_dimension_type1:: - continuous_dimension_type> && std::is_same_v)); + continuous_dimension_type> && std::is_same_v)); return deriv_1_and_2(spline_eval, coords_eval, spline_coef); } @@ -866,8 +868,8 @@ class SplineEvaluator2D ddc::ChunkSpan const spline_coef) const { - using Dim1 = typename evaluation_discrete_dimension_type1::continuous_dimension_type; - using Dim2 = typename evaluation_discrete_dimension_type2::continuous_dimension_type; + using Dim1 = continuous_dimension_type1; + using Dim2 = continuous_dimension_type2; if constexpr (bsplines_type1::is_periodic()) { if (ddc::get(coord_eval) < ddc::discrete_space().rmin() || ddc::get(coord_eval) > ddc::discrete_space().rmax()) { @@ -907,9 +909,7 @@ class SplineEvaluator2D } } return eval_no_bc( - ddc::Coordinate< - typename evaluation_discrete_dimension_type1::continuous_dimension_type, - typename evaluation_discrete_dimension_type2::continuous_dimension_type>( + ddc::Coordinate( ddc::get(coord_eval), ddc::get(coord_eval)), spline_coef); @@ -952,14 +952,10 @@ class SplineEvaluator2D double, std::experimental::extents> const vals2(vals2_ptr.data()); - ddc::Coordinate - coord_eval_interest1 = ddc::select< - typename evaluation_discrete_dimension_type1::continuous_dimension_type>( - coord_eval); - ddc::Coordinate - coord_eval_interest2 = ddc::select< - typename evaluation_discrete_dimension_type2::continuous_dimension_type>( - coord_eval); + ddc::Coordinate coord_eval_interest1 + = ddc::select(coord_eval); + ddc::Coordinate coord_eval_interest2 + = ddc::select(coord_eval); if constexpr (std::is_same_v) { jmin1 = ddc::discrete_space().eval_basis(vals1, coord_eval_interest1);