diff --git a/include/ddc/kernels/splines/spline_builder.hpp b/include/ddc/kernels/splines/spline_builder.hpp index 8188892cb..b5fc6ebb1 100644 --- a/include/ddc/kernels/splines/spline_builder.hpp +++ b/include/ddc/kernels/splines/spline_builder.hpp @@ -10,6 +10,7 @@ #include #include #include +#include #include #include @@ -88,49 +89,55 @@ class SplineBuilder /// @brief The type of the domain for the 1D interpolation mesh used by this class. using interpolation_domain_type = ddc::DiscreteDomain; - /// @brief The type of the whole domain representing interpolation points. - template - using batched_interpolation_domain_type = ddc::DiscreteDomain; + /** + * @brief The type of the whole domain representing interpolation points. + * + * @tparam The batched discrete domain on which the interpolation points are defined. + */ + template >> + using batched_interpolation_domain_type = DDom; /** * @brief The type of the batch domain (obtained by removing the dimension of interest * from the whole domain). * + * @tparam The batched discrete domain on which the interpolation points are defined. + * * Example: For batched_interpolation_domain_type = DiscreteDomain and a dimension of interest Y, * this is DiscreteDomain */ - template - using batch_domain_type = ddc::remove_dims_of_t< - batched_interpolation_domain_type, - interpolation_discrete_dimension_type>; + template >> + using batch_domain_type = ddc::remove_dims_of_t; /** * @brief The type of the whole spline domain (cartesian product of 1D spline domain * and batch domain) preserving the underlying memory layout (order of dimensions). * + * @tparam The batched discrete domain on which the interpolation points are defined. + * * Example: For batched_interpolation_domain_type = DiscreteDomain and a dimension of interest Y * (associated to a B-splines tag BSplinesY), this is DiscreteDomain. */ - template - using batched_spline_domain_type = ddc::replace_dim_of_t< - batched_interpolation_domain_type, - interpolation_discrete_dimension_type, - bsplines_type>; + template >> + using batched_spline_domain_type + = ddc::replace_dim_of_t; private: /** * @brief The type of the whole spline domain (cartesian product of the 1D spline domain * and the batch domain) with 1D spline dimension being the leading dimension. * + * @tparam The batched discrete domain on which the interpolation points are defined. + * * Example: For batched_interpolation_domain_type = DiscreteDomain and a dimension of interest Y * (associated to a B-splines tag BSplinesY), this is DiscreteDomain. */ - template + template >> using batched_spline_tr_domain_type = typename ddc::detail::convert_type_seq_to_discrete_domain_t, ddc::type_seq_remove_t< - ddc::detail::TypeSeq, + ddc::to_type_seq_t, ddc::detail::TypeSeq>>>; public: @@ -138,14 +145,14 @@ class SplineBuilder * @brief The type of the whole Deriv domain (cartesian product of 1D Deriv domain * and batch domain) preserving the underlying memory layout (order of dimensions). * + * @tparam The batched discrete domain on which the interpolation points are defined. + * * Example: For batched_interpolation_domain_type = DiscreteDomain and a dimension of interest Y, * this is DiscreteDomain,Z> */ - template - using batched_derivs_domain_type = ddc::replace_dim_of_t< - batched_interpolation_domain_type, - interpolation_discrete_dimension_type, - deriv_type>; + template >> + using batched_derivs_domain_type + = ddc::replace_dim_of_t; /// @brief Indicates if the degree of the splines is odd or even. static constexpr bool s_odd = BSplines::degree() % 2; @@ -243,9 +250,9 @@ class SplineBuilder * * @see MatrixSparse */ - template + template >> explicit SplineBuilder( - batched_interpolation_domain_type const& batched_interpolation_domain, + DDom const& batched_interpolation_domain, std::optional cols_per_chunk = std::nullopt, std::optional preconditioner_max_block_size = std::nullopt) : SplineBuilder( @@ -296,12 +303,12 @@ class SplineBuilder * Values of the function must be provided on this domain in order * to build a spline representation of the function (cartesian product of 1D interpolation_domain and batch_domain). * + * @param batched_interpolation_domain The whole domain on which the interpolation points are defined. + * * @return The domain for the interpolation mesh. */ - template - batched_interpolation_domain_type batched_interpolation_domain( - batched_interpolation_domain_type const& batched_interpolation_domain) - const noexcept + template >> + DDom batched_interpolation_domain(DDom const& batched_interpolation_domain) const noexcept { return batched_interpolation_domain; } @@ -311,11 +318,12 @@ class SplineBuilder * * Obtained by removing the dimension of interest from the whole interpolation domain. * + * @param batched_interpolation_domain The whole domain on which the interpolation points are defined. + * * @return The batch domain. */ - template - batch_domain_type batch_domain(batched_interpolation_domain_type const& - batched_interpolation_domain) const noexcept + template + batch_domain_type batch_domain(DDom const& batched_interpolation_domain) const noexcept { return ddc::remove_dims_of(batched_interpolation_domain, interpolation_domain()); } @@ -337,12 +345,13 @@ class SplineBuilder * * Spline approximations (spline-transformed functions) are computed on this domain. * + * @param batched_interpolation_domain The whole domain on which the interpolation points are defined. + * * @return The domain for the spline coefficients. */ - template - batched_spline_domain_type batched_spline_domain( - batched_interpolation_domain_type const& batched_interpolation_domain) - const noexcept + template + batched_spline_domain_type batched_spline_domain( + DDom const& batched_interpolation_domain) const noexcept { return ddc::replace_dim_of< interpolation_discrete_dimension_type, @@ -355,14 +364,15 @@ class SplineBuilder * * 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. * + * @param batched_interpolation_domain The whole domain on which the interpolation points are defined. + * * @return The (transposed) domain for the spline coefficients. */ - template - batched_spline_tr_domain_type batched_spline_tr_domain( - batched_interpolation_domain_type const& batched_interpolation_domain) - const noexcept + template + batched_spline_tr_domain_type batched_spline_tr_domain( + DDom const& batched_interpolation_domain) const noexcept { - return batched_spline_tr_domain_type( + return batched_spline_tr_domain_type( ddc::replace_dim_of( batched_spline_domain(batched_interpolation_domain), ddc::DiscreteDomain( @@ -377,12 +387,13 @@ class SplineBuilder * * This is only used with BoundCond::HERMITE boundary conditions. * + * @param batched_interpolation_domain The whole domain on which the interpolation points are defined. + * * @return The domain for the Derivs values. */ - template - batched_derivs_domain_type batched_derivs_xmin_domain( - batched_interpolation_domain_type const& batched_interpolation_domain) - const noexcept + template + batched_derivs_domain_type batched_derivs_xmin_domain( + DDom const& batched_interpolation_domain) const noexcept { return ddc::replace_dim_of( batched_interpolation_domain, @@ -396,12 +407,13 @@ class SplineBuilder * * This is only used with BoundCond::HERMITE boundary conditions. * + * @param batched_interpolation_domain The whole domain on which the interpolation points are defined. + * * @return The domain for the Derivs values. */ - template - batched_derivs_domain_type batched_derivs_xmax_domain( - batched_interpolation_domain_type const& batched_interpolation_domain) - const noexcept + template + batched_derivs_domain_type batched_derivs_xmax_domain( + DDom const& batched_interpolation_domain) const noexcept { return ddc::replace_dim_of( batched_interpolation_domain, @@ -429,24 +441,19 @@ class SplineBuilder * @param[in] derivs_xmax The values of the derivatives at the upper boundary * (used only with BoundCond::HERMITE upper boundary condition). */ - template + template void operator()( - ddc::ChunkSpan, Layout, memory_space> - spline, - ddc::ChunkSpan< - double const, - batched_interpolation_domain_type, - Layout, - memory_space> vals, + ddc::ChunkSpan, Layout, memory_space> spline, + ddc::ChunkSpan vals, std::optional, + batched_derivs_domain_type, Layout, memory_space>> derivs_xmin = std::nullopt, std::optional, + batched_derivs_domain_type, Layout, memory_space>> derivs_xmax = std::nullopt) const; @@ -755,23 +762,19 @@ template < ddc::BoundCond BcLower, ddc::BoundCond BcUpper, SplineSolver Solver> -template +template void SplineBuilder:: operator()( - ddc::ChunkSpan, Layout, memory_space> spline, - ddc::ChunkSpan< - double const, - batched_interpolation_domain_type, - Layout, - memory_space> vals, + ddc::ChunkSpan, Layout, memory_space> spline, + ddc::ChunkSpan vals, std::optional, + batched_derivs_domain_type, Layout, memory_space>> const derivs_xmin, std::optional, + batched_derivs_domain_type, Layout, memory_space>> const derivs_xmax) const { @@ -793,8 +796,6 @@ operator()( assert(ddc::DiscreteElement(derivs_xmax->domain().front()).uid() == 1); } - using batch_domain_type = batch_domain_type; - // 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 @@ -806,7 +807,7 @@ operator()( "ddc_splines_hermite_compute_lower_coefficients", exec_space(), batch_domain(batched_interpolation_domain), - KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type j) { + KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type j) { for (int i = s_nbc_xmin; i > 0; --i) { spline(ddc::DiscreteElement(s_nbc_xmin - i), j) = derivs_xmin_values(ddc::DiscreteElement(i), j) @@ -849,7 +850,7 @@ operator()( "ddc_splines_hermite_compute_upper_coefficients", exec_space(), batch_domain(batched_interpolation_domain), - KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type j) { + KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type j) { for (int i = 0; i < s_nbc_xmax; ++i) { spline(ddc::DiscreteElement(nbasis_proxy - s_nbc_xmax - i), j) @@ -869,7 +870,7 @@ operator()( "ddc_splines_transpose_rhs", exec_space(), batch_domain(batched_interpolation_domain), - KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type const j) { + KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type const j) { for (std::size_t i = 0; i < nbasis_proxy; ++i) { spline_tr(ddc::DiscreteElement(i), j) = spline(ddc::DiscreteElement(i + offset_proxy), j); @@ -887,7 +888,7 @@ operator()( "ddc_splines_transpose_back_rhs", exec_space(), batch_domain(batched_interpolation_domain), - KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type const j) { + KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type const j) { for (std::size_t i = 0; i < nbasis_proxy; ++i) { spline(ddc::DiscreteElement(i + offset_proxy), j) = spline_tr(ddc::DiscreteElement(i), j); @@ -900,7 +901,7 @@ operator()( "ddc_splines_periodic_rows_duplicate_rhs", exec_space(), batch_domain(batched_interpolation_domain), - KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type const j) { + KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type const j) { if (offset_proxy != 0) { for (int i = 0; i < offset_proxy; ++i) { spline(ddc::DiscreteElement(i), j) = spline( diff --git a/include/ddc/kernels/splines/spline_builder_2d.hpp b/include/ddc/kernels/splines/spline_builder_2d.hpp index 8f7721b06..6d6ca6332 100644 --- a/include/ddc/kernels/splines/spline_builder_2d.hpp +++ b/include/ddc/kernels/splines/spline_builder_2d.hpp @@ -94,20 +94,26 @@ class SplineBuilder2D interpolation_discrete_dimension_type1, interpolation_discrete_dimension_type2>; - /// @brief The type of the whole domain representing interpolation points. - template - using batched_interpolation_domain_type = ddc::DiscreteDomain; + /** + * @brief The type of the whole domain representing interpolation points. + * + * @tparam The batched discrete domain on which the interpolation points are defined. + */ + template >> + using batched_interpolation_domain_type = DDom; /** * @brief The type of the batch domain (obtained by removing the dimensions of interest * from the whole domain). * + * @tparam The batched discrete domain on which the interpolation points are defined. + * * Example: For batched_interpolation_domain_type = DiscreteDomain and dimensions of interest X and Y, * this is DiscreteDomain. */ - template + template >> using batch_domain_type = ddc::remove_dims_of_t< - batched_interpolation_domain_type, + DDom, interpolation_discrete_dimension_type1, interpolation_discrete_dimension_type2>; @@ -115,13 +121,15 @@ class SplineBuilder2D * @brief The type of the whole spline domain (cartesian product of 2D spline domain * and batch domain) preserving the order of dimensions. * + * @tparam The batched discrete domain on which the interpolation points are defined. + * * Example: For batched_interpolation_domain_type = DiscreteDomain and dimensions of interest X and Y * (associated to B-splines tags BSplinesX and BSplinesY), this is DiscreteDomain */ - template + template >> using batched_spline_domain_type = ddc::detail::convert_type_seq_to_discrete_domain_t, + ddc::to_type_seq_t, ddc::detail::TypeSeq< interpolation_discrete_dimension_type1, interpolation_discrete_dimension_type2>, @@ -131,37 +139,41 @@ class SplineBuilder2D * @brief The type of the whole Derivs domain (cartesian product of the 1D Deriv domain * and the associated batch domain) in the first dimension, preserving the order of dimensions. * + * @tparam The batched discrete domain on which the interpolation points are defined. + * * Example: For batched_interpolation_domain_type = DiscreteDomain and dimensions of interest X and Y, * this is DiscreteDomain, Y, Z>. */ - template + template >> using batched_derivs_domain_type1 = - typename builder_type1::template batched_derivs_domain_type; + typename builder_type1::template batched_derivs_domain_type; /** * @brief The type of the whole Derivs domain (cartesian product of the 1D Deriv domain * and the associated batch domain) in the second dimension, preserving the order of dimensions. * + * @tparam The batched discrete domain on which the interpolation points are defined. + * * Example: For batched_interpolation_domain_type = DiscreteDomain and dimensions of interest X and Y, * this is DiscreteDomain, Z>. */ - template - using batched_derivs_domain_type2 = ddc::replace_dim_of_t< - batched_interpolation_domain_type, - interpolation_discrete_dimension_type2, - deriv_type2>; + template >> + using batched_derivs_domain_type2 + = ddc::replace_dim_of_t; /** * @brief The type of the whole Derivs domain (cartesian product of the 2D Deriv domain * and the batch domain) in the second dimension, preserving the order of dimensions. * + * @tparam The batched discrete domain on which the interpolation points are defined. + * * Example: For batched_interpolation_domain_type = DiscreteDomain and dimensions of interest X and Y, * this is DiscreteDomain, Deriv, Z>. */ - template + template >> using batched_derivs_domain_type = ddc::detail::convert_type_seq_to_discrete_domain_t, + ddc::to_type_seq_t, ddc::detail::TypeSeq< interpolation_discrete_dimension_type1, interpolation_discrete_dimension_type2>, @@ -215,9 +227,9 @@ class SplineBuilder2D * * @see SplinesLinearProblemSparse */ - template + template >> explicit SplineBuilder2D( - batched_interpolation_domain_type const& batched_interpolation_domain, + DDom const& batched_interpolation_domain, std::optional cols_per_chunk = std::nullopt, std::optional preconditioner_max_block_size = std::nullopt) : SplineBuilder2D( @@ -270,12 +282,12 @@ class SplineBuilder2D * Values of the function must be provided on this domain in order * to build a spline representation of the function (cartesian product of 2D interpolation_domain and batch_domain). * + * @param batched_interpolation_domain The whole domain on which the interpolation points are defined. + * * @return The domain for the interpolation mesh. */ - template - batched_interpolation_domain_type batched_interpolation_domain( - batched_interpolation_domain_type const& batched_interpolation_domain) - const noexcept + template >> + DDom batched_interpolation_domain(DDom const& batched_interpolation_domain) const noexcept { return batched_interpolation_domain; } @@ -285,11 +297,12 @@ class SplineBuilder2D * * Obtained by removing the dimensions of interest from the whole interpolation domain. * + * @param batched_interpolation_domain The whole domain on which the interpolation points are defined. + * * @return The batch domain. */ - template - batch_domain_type batch_domain(batched_interpolation_domain_type const& - batched_interpolation_domain) const noexcept + template >> + batch_domain_type batch_domain(DDom const& batched_interpolation_domain) const noexcept { return ddc::remove_dims_of(batched_interpolation_domain, interpolation_domain()); } @@ -313,12 +326,13 @@ class SplineBuilder2D * * Spline approximations (spline-transformed functions) are computed on this domain. * + * @param batched_interpolation_domain The whole domain on which the interpolation points are defined. + * * @return The domain for the spline coefficients. */ - template - batched_spline_domain_type batched_spline_domain( - batched_interpolation_domain_type const& batched_interpolation_domain) - const noexcept + template >> + batched_spline_domain_type batched_spline_domain( + DDom const& batched_interpolation_domain) const noexcept { return ddc::replace_dim_of( ddc::replace_dim_of< @@ -363,60 +377,55 @@ class SplineBuilder2D * The values of the the cross-derivatives at the upper boundary in the first dimension * and the upper boundary in the second dimension. */ - template + template void operator()( - ddc::ChunkSpan, Layout, memory_space> - spline, - ddc::ChunkSpan< - double const, - batched_interpolation_domain_type, - Layout, - memory_space> vals, + ddc::ChunkSpan, Layout, memory_space> spline, + ddc::ChunkSpan vals, std::optional, + batched_derivs_domain_type1, Layout, memory_space>> derivs_min1 = std::nullopt, std::optional, + batched_derivs_domain_type1, Layout, memory_space>> derivs_max1 = std::nullopt, std::optional, + batched_derivs_domain_type2, Layout, memory_space>> derivs_min2 = std::nullopt, std::optional, + batched_derivs_domain_type2, Layout, memory_space>> derivs_max2 = std::nullopt, std::optional, + batched_derivs_domain_type, Layout, memory_space>> mixed_derivs_min1_min2 = std::nullopt, std::optional, + batched_derivs_domain_type, Layout, memory_space>> mixed_derivs_max1_min2 = std::nullopt, std::optional, + batched_derivs_domain_type, Layout, memory_space>> mixed_derivs_min1_max2 = std::nullopt, std::optional, + batched_derivs_domain_type, Layout, memory_space>> mixed_derivs_max1_max2 = std::nullopt) const; @@ -435,7 +444,7 @@ template < ddc::BoundCond BcLower2, ddc::BoundCond BcUpper2, ddc::SplineSolver Solver> -template +template void SplineBuilder2D< ExecSpace, MemorySpace, @@ -449,50 +458,46 @@ void SplineBuilder2D< BcUpper2, Solver>:: operator()( - ddc::ChunkSpan, Layout, memory_space> spline, - ddc::ChunkSpan< - double const, - batched_interpolation_domain_type, - Layout, - memory_space> vals, + ddc::ChunkSpan, Layout, memory_space> spline, + ddc::ChunkSpan vals, std::optional, + batched_derivs_domain_type1, Layout, memory_space>> const derivs_min1, std::optional, + batched_derivs_domain_type1, Layout, memory_space>> const derivs_max1, std::optional, + batched_derivs_domain_type2, Layout, memory_space>> const derivs_min2, std::optional, + batched_derivs_domain_type2, Layout, memory_space>> const derivs_max2, std::optional, + batched_derivs_domain_type, Layout, memory_space>> const mixed_derivs_min1_min2, std::optional, + batched_derivs_domain_type, Layout, memory_space>> const mixed_derivs_max1_min2, std::optional, + batched_derivs_domain_type, Layout, memory_space>> const mixed_derivs_min1_max2, std::optional, + batched_derivs_domain_type, Layout, memory_space>> const mixed_derivs_max1_max2) const { diff --git a/include/ddc/kernels/splines/spline_evaluator.hpp b/include/ddc/kernels/splines/spline_evaluator.hpp index 105145090..7a924d7b5 100644 --- a/include/ddc/kernels/splines/spline_evaluator.hpp +++ b/include/ddc/kernels/splines/spline_evaluator.hpp @@ -73,9 +73,13 @@ class SplineEvaluator /// @brief The type of the domain for the 1D evaluation mesh used by this class. using evaluation_domain_type = ddc::DiscreteDomain; - /// @brief The type of the whole domain representing evaluation points. - template - using batched_evaluation_domain_type = ddc::DiscreteDomain; + /** + * @brief The type of the whole domain representing evaluation points. + * + * @tparam The batched discrete domain on which the interpolation points are defined. + */ + template >> + using batched_evaluation_domain_type = DDom; /// @brief The type of the 1D spline domain corresponding to the dimension of interest. using spline_domain_type = ddc::DiscreteDomain; @@ -83,21 +87,21 @@ class SplineEvaluator /** * @brief The type of the batch domain (obtained by removing the dimension of interest * from the whole domain). + * + * @tparam The batched discrete domain on which the interpolation points are defined. */ - template - using batch_domain_type = ddc::remove_dims_of_t< - batched_evaluation_domain_type, - evaluation_discrete_dimension_type>; + template >> + using batch_domain_type = ddc::remove_dims_of_t; /** * @brief The type of the whole spline domain (cartesian product of 1D spline domain * and batch domain) preserving the order of dimensions. + * + * @tparam The batched discrete domain on which the interpolation points are defined. */ - template - using batched_spline_domain_type = ddc::replace_dim_of_t< - batched_evaluation_domain_type, - evaluation_discrete_dimension_type, - bsplines_type>; + template >> + using batched_spline_domain_type + = ddc::replace_dim_of_t; /// @brief The type of the extrapolation rule at the lower boundary. using lower_extrapolation_rule_type = LowerExtrapolationRule; @@ -262,33 +266,26 @@ class SplineEvaluator * the set of 1D spline coefficients retained to perform the evaluation). * @param[in] spline_coef A ChunkSpan storing the spline coefficients. */ - template + template void operator()( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, - ddc::ChunkSpan< - ddc::Coordinate const, - batched_evaluation_domain_type, - Layout2, - memory_space> const coords_eval, + ddc::ChunkSpan const spline_eval, + ddc::ChunkSpan const, DDom, Layout2, memory_space> const + coords_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout3, memory_space> const spline_coef) const { evaluation_domain_type const evaluation_domain(spline_eval.domain()); - batch_domain_type const batch_domain(spline_eval.domain()); + batch_domain_type const batch_domain(spline_eval.domain()); ddc::parallel_for_each( "ddc_splines_evaluate", exec_space(), batch_domain, KOKKOS_CLASS_LAMBDA( - typename batch_domain_type::discrete_element_type const j) { + typename batch_domain_type::discrete_element_type const j) { const auto spline_eval_1D = spline_eval[j]; const auto coords_eval_1D = coords_eval[j]; const auto spline_coef_1D = spline_coef[j]; @@ -316,28 +313,24 @@ class SplineEvaluator * of the mesh. * @param[in] spline_coef A ChunkSpan storing the spline coefficients. */ - template + template void operator()( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, + ddc::ChunkSpan const spline_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout2, memory_space> const spline_coef) const { evaluation_domain_type const evaluation_domain(spline_eval.domain()); - batch_domain_type const batch_domain(spline_eval.domain()); + batch_domain_type const batch_domain(spline_eval.domain()); ddc::parallel_for_each( "ddc_splines_evaluate", exec_space(), batch_domain, KOKKOS_CLASS_LAMBDA( - typename batch_domain_type::discrete_element_type const j) { + typename batch_domain_type::discrete_element_type const j) { const auto spline_eval_1D = spline_eval[j]; const auto spline_coef_1D = spline_coef[j]; for (auto const i : evaluation_domain) { @@ -386,33 +379,26 @@ class SplineEvaluator * the set of 1D spline coefficients retained to perform the evaluation). * @param[in] spline_coef A ChunkSpan storing the spline coefficients. */ - template + template void deriv( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, - ddc::ChunkSpan< - ddc::Coordinate const, - batched_evaluation_domain_type, - Layout2, - memory_space> const coords_eval, + ddc::ChunkSpan const spline_eval, + ddc::ChunkSpan const, DDom, Layout2, memory_space> const + coords_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout3, memory_space> const spline_coef) const { evaluation_domain_type const evaluation_domain(spline_eval.domain()); - batch_domain_type const batch_domain(spline_eval.domain()); + batch_domain_type const batch_domain(spline_eval.domain()); ddc::parallel_for_each( "ddc_splines_differentiate", exec_space(), batch_domain, KOKKOS_CLASS_LAMBDA( - typename batch_domain_type::discrete_element_type const j) { + typename batch_domain_type::discrete_element_type const j) { const auto spline_eval_1D = spline_eval[j]; const auto coords_eval_1D = coords_eval[j]; const auto spline_coef_1D = spline_coef[j]; @@ -437,28 +423,24 @@ class SplineEvaluator * @param[out] spline_eval The derivatives of the spline function at the coordinates. * @param[in] spline_coef A ChunkSpan storing the spline coefficients. */ - template + template void deriv( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, + ddc::ChunkSpan const spline_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout2, memory_space> const spline_coef) const { evaluation_domain_type const evaluation_domain(spline_eval.domain()); - batch_domain_type const batch_domain(spline_eval.domain()); + batch_domain_type const batch_domain(spline_eval.domain()); ddc::parallel_for_each( "ddc_splines_differentiate", exec_space(), batch_domain, KOKKOS_CLASS_LAMBDA( - typename batch_domain_type::discrete_element_type const j) { + typename batch_domain_type::discrete_element_type const j) { const auto spline_eval_1D = spline_eval[j]; const auto spline_coef_1D = spline_coef[j]; for (auto const i : evaluation_domain) { @@ -483,17 +465,21 @@ class SplineEvaluator * points represented by this domain are unused and irrelevant. * @param[in] spline_coef A ChunkSpan storing the spline coefficients. */ - template + template void integrate( - ddc::ChunkSpan, Layout1, memory_space> const - integrals, - ddc::ChunkSpan< - double const, - batched_spline_domain_type, - Layout2, - memory_space> const spline_coef) const + ddc::ChunkSpan const integrals, + ddc::ChunkSpan const + spline_coef) const { - batch_domain_type const batch_domain(integrals.domain()); + static_assert( + ddc::in_tags_v>, + "The spline coefficients domain must contain the bsplines dimension"); + using batch_domain_type = ddc::remove_dims_of_t; + static_assert( + std::is_same_v, + "The integrals domain must only contain the batch dimensions"); + + BatchDDom const batch_domain(integrals.domain()); ddc::Chunk values_alloc( ddc::DiscreteDomain(spline_coef.domain()), ddc::KokkosAllocator()); @@ -504,7 +490,7 @@ class SplineEvaluator "ddc_splines_integrate", exec_space(), batch_domain, - KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type const j) { + KOKKOS_LAMBDA(typename BatchDDom::discrete_element_type const j) { integrals(j) = 0; for (typename spline_domain_type::discrete_element_type const i : values.domain()) { diff --git a/include/ddc/kernels/splines/spline_evaluator_2d.hpp b/include/ddc/kernels/splines/spline_evaluator_2d.hpp index deeef4f6e..494d228e1 100644 --- a/include/ddc/kernels/splines/spline_evaluator_2d.hpp +++ b/include/ddc/kernels/splines/spline_evaluator_2d.hpp @@ -97,9 +97,13 @@ class SplineEvaluator2D evaluation_discrete_dimension_type1, evaluation_discrete_dimension_type2>; - /// @brief The type of the whole domain representing evaluation points. - template - using batched_evaluation_domain_type = ddc::DiscreteDomain; + /** + * @brief The type of the whole domain representing evaluation points. + * + * @tparam The batched discrete domain on which the interpolation points are defined. + */ + template >> + using batched_evaluation_domain_type = DDom; /// @brief The type of the 1D spline domain corresponding to the first dimension of interest. using spline_domain_type1 = ddc::DiscreteDomain; @@ -113,21 +117,25 @@ class SplineEvaluator2D /** * @brief The type of the batch domain (obtained by removing the dimensions of interest * from the whole domain). + * + * @tparam The batched discrete domain on which the interpolation points are defined. */ - template + template >> using batch_domain_type = typename ddc::remove_dims_of_t< - batched_evaluation_domain_type, + DDom, evaluation_discrete_dimension_type1, evaluation_discrete_dimension_type2>; /** * @brief The type of the whole spline domain (cartesian product of 2D spline domain * and batch domain) preserving the underlying memory layout (order of dimensions). + * + * @tparam The batched discrete domain on which the interpolation points are defined. */ - template + template >> using batched_spline_domain_type = typename ddc::detail::convert_type_seq_to_discrete_domain_t, + ddc::to_type_seq_t, ddc::detail::TypeSeq< evaluation_discrete_dimension_type1, evaluation_discrete_dimension_type2>, @@ -373,25 +381,18 @@ class SplineEvaluator2D * the set of 2D spline coefficients retained to perform the evaluation). * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients. */ - template + template void operator()( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, - ddc::ChunkSpan< - ddc::Coordinate const, - batched_evaluation_domain_type, - Layout2, - memory_space> const coords_eval, + ddc::ChunkSpan const spline_eval, + ddc::ChunkSpan const, DDom, Layout2, memory_space> const + coords_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout3, memory_space> const spline_coef) const { - batch_domain_type const batch_domain(coords_eval.domain()); + batch_domain_type const batch_domain(coords_eval.domain()); evaluation_domain_type1 const evaluation_domain1(spline_eval.domain()); evaluation_domain_type2 const evaluation_domain2(spline_eval.domain()); ddc::parallel_for_each( @@ -399,7 +400,7 @@ class SplineEvaluator2D exec_space(), batch_domain, KOKKOS_CLASS_LAMBDA( - typename batch_domain_type::discrete_element_type const j) { + typename batch_domain_type::discrete_element_type const j) { const auto spline_eval_2D = spline_eval[j]; const auto coords_eval_2D = coords_eval[j]; const auto spline_coef_2D = spline_coef[j]; @@ -426,20 +427,16 @@ class SplineEvaluator2D * @param[out] spline_eval The values of the 2D spline function at their coordinates. * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients. */ - template + template void operator()( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, + ddc::ChunkSpan const spline_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout2, memory_space> const spline_coef) const { - batch_domain_type const batch_domain(spline_eval.domain()); + batch_domain_type const batch_domain(spline_eval.domain()); evaluation_domain_type1 const evaluation_domain1(spline_eval.domain()); evaluation_domain_type2 const evaluation_domain2(spline_eval.domain()); ddc::parallel_for_each( @@ -447,7 +444,7 @@ class SplineEvaluator2D exec_space(), batch_domain, KOKKOS_CLASS_LAMBDA( - typename batch_domain_type::discrete_element_type const j) { + typename batch_domain_type::discrete_element_type const j) { const auto spline_eval_2D = spline_eval[j]; const auto spline_coef_2D = spline_coef[j]; for (auto const i1 : evaluation_domain1) { @@ -609,25 +606,18 @@ class SplineEvaluator2D * the set of 2D spline coefficients retained to perform the evaluation). * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients. */ - template + template void deriv_dim_1( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, - ddc::ChunkSpan< - ddc::Coordinate const, - batched_evaluation_domain_type, - Layout2, - memory_space> const coords_eval, + ddc::ChunkSpan const spline_eval, + ddc::ChunkSpan const, DDom, Layout2, memory_space> const + coords_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout3, memory_space> const spline_coef) const { - batch_domain_type const batch_domain(coords_eval.domain()); + batch_domain_type const batch_domain(coords_eval.domain()); evaluation_domain_type1 const evaluation_domain1(spline_eval.domain()); evaluation_domain_type2 const evaluation_domain2(spline_eval.domain()); ddc::parallel_for_each( @@ -635,7 +625,7 @@ class SplineEvaluator2D exec_space(), batch_domain, KOKKOS_CLASS_LAMBDA( - typename batch_domain_type::discrete_element_type const j) { + typename batch_domain_type::discrete_element_type const j) { const auto spline_eval_2D = spline_eval[j]; const auto coords_eval_2D = coords_eval[j]; const auto spline_coef_2D = spline_coef[j]; @@ -662,20 +652,16 @@ class SplineEvaluator2D * @param[out] spline_eval The derivatives of the 2D spline function at the desired coordinates. * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients. */ - template + template void deriv_dim_1( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, + ddc::ChunkSpan const spline_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout2, memory_space> const spline_coef) const { - batch_domain_type const batch_domain(spline_eval.domain()); + batch_domain_type const batch_domain(spline_eval.domain()); evaluation_domain_type1 const evaluation_domain1(spline_eval.domain()); evaluation_domain_type2 const evaluation_domain2(spline_eval.domain()); ddc::parallel_for_each( @@ -683,7 +669,7 @@ class SplineEvaluator2D exec_space(), batch_domain, KOKKOS_CLASS_LAMBDA( - typename batch_domain_type::discrete_element_type const j) { + typename batch_domain_type::discrete_element_type const j) { const auto spline_eval_2D = spline_eval[j]; const auto spline_coef_2D = spline_coef[j]; for (auto const i1 : evaluation_domain1) { @@ -716,25 +702,18 @@ class SplineEvaluator2D * the set of 2D spline coefficients retained to perform the evaluation). * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients. */ - template + template void deriv_dim_2( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, - ddc::ChunkSpan< - ddc::Coordinate const, - batched_evaluation_domain_type, - Layout2, - memory_space> const coords_eval, + ddc::ChunkSpan const spline_eval, + ddc::ChunkSpan const, DDom, Layout2, memory_space> const + coords_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout3, memory_space> const spline_coef) const { - batch_domain_type const batch_domain(coords_eval.domain()); + batch_domain_type const batch_domain(coords_eval.domain()); evaluation_domain_type1 const evaluation_domain1(spline_eval.domain()); evaluation_domain_type2 const evaluation_domain2(spline_eval.domain()); ddc::parallel_for_each( @@ -742,7 +721,7 @@ class SplineEvaluator2D exec_space(), batch_domain, KOKKOS_CLASS_LAMBDA( - typename batch_domain_type::discrete_element_type const j) { + typename batch_domain_type::discrete_element_type const j) { const auto spline_eval_2D = spline_eval[j]; const auto coords_eval_2D = coords_eval[j]; const auto spline_coef_2D = spline_coef[j]; @@ -769,20 +748,16 @@ class SplineEvaluator2D * @param[out] spline_eval The derivatives of the 2D spline function at the desired coordinates. * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients. */ - template + template void deriv_dim_2( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, + ddc::ChunkSpan const spline_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout2, memory_space> const spline_coef) const { - batch_domain_type const batch_domain(spline_eval.domain()); + batch_domain_type const batch_domain(spline_eval.domain()); evaluation_domain_type1 const evaluation_domain1(spline_eval.domain()); evaluation_domain_type2 const evaluation_domain2(spline_eval.domain()); ddc::parallel_for_each( @@ -790,7 +765,7 @@ class SplineEvaluator2D exec_space(), batch_domain, KOKKOS_CLASS_LAMBDA( - typename batch_domain_type::discrete_element_type const j) { + typename batch_domain_type::discrete_element_type const j) { const auto spline_eval_2D = spline_eval[j]; const auto spline_coef_2D = spline_coef[j]; for (auto const i1 : evaluation_domain1) { @@ -823,25 +798,18 @@ class SplineEvaluator2D * the set of 2D spline coefficients retained to perform the evaluation). * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients. */ - template + template void deriv_1_and_2( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, - ddc::ChunkSpan< - ddc::Coordinate const, - batched_evaluation_domain_type, - Layout2, - memory_space> const coords_eval, + ddc::ChunkSpan const spline_eval, + ddc::ChunkSpan const, DDom, Layout2, memory_space> const + coords_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout3, memory_space> const spline_coef) const { - batch_domain_type const batch_domain(coords_eval.domain()); + batch_domain_type const batch_domain(coords_eval.domain()); evaluation_domain_type1 const evaluation_domain1(spline_eval.domain()); evaluation_domain_type2 const evaluation_domain2(spline_eval.domain()); ddc::parallel_for_each( @@ -849,7 +817,7 @@ class SplineEvaluator2D exec_space(), batch_domain, KOKKOS_CLASS_LAMBDA( - typename batch_domain_type::discrete_element_type const j) { + typename batch_domain_type::discrete_element_type const j) { const auto spline_eval_2D = spline_eval[j]; const auto coords_eval_2D = coords_eval[j]; const auto spline_coef_2D = spline_coef[j]; @@ -876,20 +844,16 @@ class SplineEvaluator2D * @param[out] spline_eval The cross-derivatives of the 2D spline function at the desired coordinates. * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients. */ - template + template void deriv_1_and_2( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, + ddc::ChunkSpan const spline_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout2, memory_space> const spline_coef) const { - batch_domain_type const batch_domain(spline_eval.domain()); + batch_domain_type const batch_domain(spline_eval.domain()); evaluation_domain_type1 const evaluation_domain1(spline_eval.domain()); evaluation_domain_type2 const evaluation_domain2(spline_eval.domain()); ddc::parallel_for_each( @@ -897,7 +861,7 @@ class SplineEvaluator2D exec_space(), batch_domain, KOKKOS_CLASS_LAMBDA( - typename batch_domain_type::discrete_element_type const j) { + typename batch_domain_type::discrete_element_type const j) { const auto spline_eval_2D = spline_eval[j]; const auto spline_coef_2D = spline_coef[j]; for (auto const i1 : evaluation_domain1) { @@ -936,22 +900,15 @@ class SplineEvaluator2D class Layout1, class Layout2, class Layout3, - class... CoordsDims, - class... DDimX> + class DDom, + class... CoordsDims> void deriv( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, - ddc::ChunkSpan< - ddc::Coordinate const, - batched_evaluation_domain_type, - Layout2, - memory_space> const coords_eval, + ddc::ChunkSpan const spline_eval, + ddc::ChunkSpan const, DDom, Layout2, memory_space> const + coords_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout3, memory_space> const spline_coef) const { @@ -987,16 +944,12 @@ class SplineEvaluator2D * @param[out] spline_eval The derivatives of the 2D spline function at the desired coordinates. * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients. */ - template + template void deriv( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, + ddc::ChunkSpan const spline_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout2, memory_space> const spline_coef) const { @@ -1047,22 +1000,15 @@ class SplineEvaluator2D class Layout1, class Layout2, class Layout3, - class... CoordsDims, - class... DDimX> + class DDom, + class... CoordsDims> void deriv2( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, - ddc::ChunkSpan< - ddc::Coordinate const, - batched_evaluation_domain_type, - Layout2, - memory_space> const coords_eval, + ddc::ChunkSpan const spline_eval, + ddc::ChunkSpan const, DDom, Layout2, memory_space> const + coords_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout3, memory_space> const spline_coef) const { @@ -1096,16 +1042,12 @@ class SplineEvaluator2D * @param[out] spline_eval The derivatives of the 2D spline function at the desired coordinates. * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients. */ - template + template void deriv2( - ddc::ChunkSpan< - double, - batched_evaluation_domain_type, - Layout1, - memory_space> const spline_eval, + ddc::ChunkSpan const spline_eval, ddc::ChunkSpan< double const, - batched_spline_domain_type, + batched_spline_domain_type, Layout2, memory_space> const spline_coef) const { @@ -1134,17 +1076,24 @@ class SplineEvaluator2D * points represented by this domain are unused and irrelevant. * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients. */ - template + template void integrate( - ddc::ChunkSpan, Layout1, memory_space> const - integrals, - ddc::ChunkSpan< - double const, - batched_spline_domain_type, - Layout2, - memory_space> const spline_coef) const + ddc::ChunkSpan const integrals, + ddc::ChunkSpan const + spline_coef) const { - batch_domain_type batch_domain(integrals.domain()); + static_assert( + ddc::type_seq_contains_v< + ddc::detail::TypeSeq, + to_type_seq_t>, + "The spline coefficients domain must contain the bsplines dimensions"); + using batch_domain_type + = ddc::remove_dims_of_t; + static_assert( + std::is_same_v, + "The integrals domain must only contain the batch dimensions"); + + batch_domain_type batch_domain(integrals.domain()); ddc::Chunk values1_alloc( ddc::DiscreteDomain(spline_coef.domain()), ddc::KokkosAllocator()); @@ -1160,7 +1109,7 @@ class SplineEvaluator2D "ddc_splines_integrate_bsplines", exec_space(), batch_domain, - KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type const j) { + KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type const j) { integrals(j) = 0; for (typename spline_domain_type1::discrete_element_type const i1 : values1.domain()) { diff --git a/tests/splines/batched_spline_builder.cpp b/tests/splines/batched_spline_builder.cpp index b0fa0b520..5d7147eaa 100644 --- a/tests/splines/batched_spline_builder.cpp +++ b/tests/splines/batched_spline_builder.cpp @@ -308,7 +308,8 @@ void BatchedSplineTest() BSplines, DDimI, extrapolation_rule_type, - extrapolation_rule_type> const spline_evaluator_batched(extrapolation_rule, extrapolation_rule); + extrapolation_rule_type> const + spline_evaluator_batched(extrapolation_rule, extrapolation_rule); // Instantiate chunk of coordinates of dom_interpolation ddc::Chunk coords_eval_alloc(dom_vals, ddc::KokkosAllocator, MemorySpace>()); @@ -332,7 +333,7 @@ void BatchedSplineTest() // Call spline_evaluator on the same mesh we started with spline_evaluator_batched(spline_eval, coords_eval.span_cview(), coef.span_cview()); spline_evaluator_batched.deriv(spline_eval_deriv, coords_eval.span_cview(), coef.span_cview()); - spline_evaluator_batched.template integrate(spline_eval_integrals, coef.span_cview()); + spline_evaluator_batched.integrate(spline_eval_integrals, coef.span_cview()); // Checking errors (we recover the initial values) double const max_norm_error = ddc::parallel_transform_reduce( @@ -359,7 +360,7 @@ void BatchedSplineTest() 0., ddc::reducer::max(), KOKKOS_LAMBDA(typename decltype(spline_builder)::template batch_domain_type< - DDims...>::discrete_element_type const e) { + ddc::DiscreteDomain>::discrete_element_type const e) { return Kokkos::abs( spline_eval_integrals(e) - evaluator.deriv(xN(), -1) + evaluator.deriv(x0(), -1)); diff --git a/tests/splines/non_periodic_spline_builder.cpp b/tests/splines/non_periodic_spline_builder.cpp index 6c3972824..cead12210 100644 --- a/tests/splines/non_periodic_spline_builder.cpp +++ b/tests/splines/non_periodic_spline_builder.cpp @@ -201,7 +201,7 @@ void TestNonPeriodicSplineBuilderTestIdentity() ddc::Chunk integral( spline_builder.batch_domain(interpolation_domain), ddc::KokkosAllocator()); - spline_evaluator.integrate(integral.span_view(), coef.span_cview()); + spline_evaluator.integrate(integral.span_view(), coef.span_cview()); ddc::Chunk< double, diff --git a/tests/splines/periodic_spline_builder.cpp b/tests/splines/periodic_spline_builder.cpp index 38108472a..683f1587c 100644 --- a/tests/splines/periodic_spline_builder.cpp +++ b/tests/splines/periodic_spline_builder.cpp @@ -141,7 +141,7 @@ void TestPeriodicSplineBuilderTestIdentity() ddc::Chunk integral( spline_builder.batch_domain(interpolation_domain), ddc::KokkosAllocator()); - spline_evaluator.integrate(integral.span_view(), coef.span_cview()); + spline_evaluator.integrate(integral.span_view(), coef.span_cview()); ddc::Chunk, ddc::KokkosAllocator> quadrature_coefficients_alloc;