diff --git a/src/geometryRTheta/advection_field/advection_field_rp.hpp b/src/geometryRTheta/advection_field/advection_field_rp.hpp
index f9b152805..578dce584 100644
--- a/src/geometryRTheta/advection_field/advection_field_rp.hpp
+++ b/src/geometryRTheta/advection_field/advection_field_rp.hpp
@@ -98,7 +98,7 @@ class AdvectionFieldFinder
 private:
     Mapping const& m_mapping;
 
-    PolarSplineEvaluator<PolarBSplinesRTheta, ddc::NullExtrapolationRule> const
+    PolarSplineEvaluator<PolarBSplinesRTheta, ddc::NullExtrapolationRule, Kokkos::HostSpace> const
             m_polar_spline_evaluator;
 
     SplineRThetaEvaluatorNullBound const m_spline_evaluator;
@@ -188,7 +188,7 @@ class AdvectionFieldFinder
      *      The advection field on the physical axis. 
      */
     void operator()(
-            SplinePolar& electrostatic_potential_coef,
+            host_t<SplinePolar>& electrostatic_potential_coef,
             host_t<DVectorFieldRTheta<X, Y>> advection_field_xy) const
     {
         compute_advection_field_XY(
@@ -223,7 +223,8 @@ class AdvectionFieldFinder
                             Evaluator,
                             PolarSplineEvaluator<
                                     PolarBSplinesRTheta,
-                                    ddc::NullExtrapolationRule>> && std::is_same_v<SplineType, SplinePolar>));
+                                    ddc::NullExtrapolationRule,
+                                    Kokkos::HostSpace>> && std::is_same_v<SplineType, host_t<SplinePolar>>));
 
         IdxRangeRTheta const grid = get_idx_range(advection_field_xy);
         host_t<DVectorFieldMemRTheta<X, Y>> electric_field(grid);
@@ -412,7 +413,7 @@ class AdvectionFieldFinder
      *      The advection field on the physical axis at the O-point. 
      */
     void operator()(
-            SplinePolar& electrostatic_potential_coef,
+            host_t<SplinePolar>& electrostatic_potential_coef,
             host_t<DVectorFieldRTheta<R, Theta>> advection_field_rp,
             CoordXY& advection_field_xy_center) const
     {
@@ -453,7 +454,8 @@ class AdvectionFieldFinder
                             Evaluator,
                             PolarSplineEvaluator<
                                     PolarBSplinesRTheta,
-                                    ddc::NullExtrapolationRule>> && std::is_same_v<SplineType, SplinePolar>));
+                                    ddc::NullExtrapolationRule,
+                                    Kokkos::HostSpace>> && std::is_same_v<SplineType, host_t<SplinePolar>>));
 
         IdxRangeRTheta const grid_without_Opoint = get_idx_range(advection_field_rp);
 
diff --git a/src/geometryRTheta/poisson/polarpoissonlikesolver.hpp b/src/geometryRTheta/poisson/polarpoissonlikesolver.hpp
index 17396d557..c96acf6ab 100644
--- a/src/geometryRTheta/poisson/polarpoissonlikesolver.hpp
+++ b/src/geometryRTheta/poisson/polarpoissonlikesolver.hpp
@@ -186,7 +186,8 @@ class PolarSplineFEMPoissonLikeSolver
 
     host_t<FieldMem<double, IdxRangeQuadratureRTheta>> int_volume;
 
-    PolarSplineEvaluator<PolarBSplinesRTheta, ddc::NullExtrapolationRule> m_polar_spline_evaluator;
+    PolarSplineEvaluator<PolarBSplinesRTheta, ddc::NullExtrapolationRule, Kokkos::HostSpace>
+            m_polar_spline_evaluator;
     std::unique_ptr<MatrixBatchCsr<Kokkos::DefaultExecutionSpace, MatrixBatchCsrSolver::CG>>
             m_gko_matrix;
     const int batch_idx {0};
@@ -653,7 +654,7 @@ class PolarSplineFEMPoissonLikeSolver
      *      The spline representation of the solution @f$\phi@f$.
      */
     template <class RHSFunction>
-    void operator()(RHSFunction const& rhs, SplinePolar& spline) const
+    void operator()(RHSFunction const& rhs, host_t<SplinePolar>& spline) const
     {
         const int b_size = ddc::discrete_space<PolarBSplinesRTheta>().nbasis()
                            - ddc::discrete_space<BSplinesTheta_Polar>().nbasis();
@@ -806,7 +807,7 @@ class PolarSplineFEMPoissonLikeSolver
     {
         IdxRangeBSTheta_Polar polar_idx_range(
                 ddc::discrete_space<BSplinesTheta_Polar>().full_domain());
-        SplinePolar
+        host_t<SplinePolar>
                 spline(PolarBSplinesRTheta::singular_idx_range<PolarBSplinesRTheta>(),
                        IdxRangeBSRTheta(radial_bsplines, polar_idx_range));
 
diff --git a/src/geometryRTheta/time_solver/bsl_predcorr.hpp b/src/geometryRTheta/time_solver/bsl_predcorr.hpp
index e9bba40d4..8a52a9f04 100644
--- a/src/geometryRTheta/time_solver/bsl_predcorr.hpp
+++ b/src/geometryRTheta/time_solver/bsl_predcorr.hpp
@@ -123,11 +123,11 @@ class BslPredCorrRTheta : public ITimeSolverRTheta
                 IdxStep<BSplinesR> {PolarBSplinesRTheta::continuity + 1}));
         IdxRangeBSTheta polar_idx_range(ddc::discrete_space<BSplinesTheta>().full_domain());
 
-        SplinePolar electrostatic_potential_coef(
+        host_t<SplinePolar> electrostatic_potential_coef(
                 PolarBSplinesRTheta::singular_idx_range<PolarBSplinesRTheta>(),
                 IdxRangeBSRTheta(radial_bsplines, polar_idx_range));
         ddc::NullExtrapolationRule extrapolation_rule;
-        PolarSplineEvaluator<PolarBSplinesRTheta, ddc::NullExtrapolationRule>
+        PolarSplineEvaluator<PolarBSplinesRTheta, ddc::NullExtrapolationRule, Kokkos::HostSpace>
                 polar_spline_evaluator(extrapolation_rule);
 
 
diff --git a/src/geometryRTheta/time_solver/bsl_predcorr_second_order_explicit.hpp b/src/geometryRTheta/time_solver/bsl_predcorr_second_order_explicit.hpp
index b715c8e15..b05802784 100644
--- a/src/geometryRTheta/time_solver/bsl_predcorr_second_order_explicit.hpp
+++ b/src/geometryRTheta/time_solver/bsl_predcorr_second_order_explicit.hpp
@@ -158,12 +158,12 @@ class BslExplicitPredCorrRTheta : public ITimeSolverRTheta
         // --- Electrostatic potential (phi). -------------------------------------------------------------
         host_t<DFieldMemRTheta> electrical_potential(grid);
 
-        SplinePolar electrostatic_potential_coef(
+        host_t<SplinePolar> electrostatic_potential_coef(
                 PolarBSplinesRTheta::singular_idx_range<PolarBSplinesRTheta>(),
                 IdxRangeBSRTheta(radial_bsplines, polar_idx_range));
 
         ddc::NullExtrapolationRule extrapolation_rule;
-        PolarSplineEvaluator<PolarBSplinesRTheta, ddc::NullExtrapolationRule>
+        PolarSplineEvaluator<PolarBSplinesRTheta, ddc::NullExtrapolationRule, Kokkos::HostSpace>
                 polar_spline_evaluator(extrapolation_rule);
 
         // --- For the computation of advection field from the electrostatic potential (phi): -------------
diff --git a/src/geometryRTheta/time_solver/bsl_predcorr_second_order_implicit.hpp b/src/geometryRTheta/time_solver/bsl_predcorr_second_order_implicit.hpp
index 117557349..73a45c0a7 100644
--- a/src/geometryRTheta/time_solver/bsl_predcorr_second_order_implicit.hpp
+++ b/src/geometryRTheta/time_solver/bsl_predcorr_second_order_implicit.hpp
@@ -158,12 +158,12 @@ class BslImplicitPredCorrRTheta : public ITimeSolverRTheta
         // --- Electrostatic potential (phi). -------------------------------------------------------------
         host_t<DFieldMemRTheta> electrical_potential(grid);
 
-        SplinePolar electrostatic_potential_coef(
+        host_t<SplinePolar> electrostatic_potential_coef(
                 PolarBSplinesRTheta::singular_idx_range<PolarBSplinesRTheta>(),
                 IdxRangeBSRTheta(radial_bsplines, polar_idx_range));
 
         ddc::NullExtrapolationRule extrapolation_rule;
-        PolarSplineEvaluator<PolarBSplinesRTheta, ddc::NullExtrapolationRule>
+        PolarSplineEvaluator<PolarBSplinesRTheta, ddc::NullExtrapolationRule, Kokkos::HostSpace>
                 polar_spline_evaluator(extrapolation_rule);
 
         // --- For the computation of advection field from the electrostatic potential (phi): -------------
diff --git a/src/utils/ddc_alias_inline_functions.hpp b/src/utils/ddc_alias_inline_functions.hpp
index 40b1d9838..1ecf36a53 100644
--- a/src/utils/ddc_alias_inline_functions.hpp
+++ b/src/utils/ddc_alias_inline_functions.hpp
@@ -82,17 +82,19 @@ inline constexpr bool
 template <class ElementType, class IdxRangeType, class MemoryType>
 inline constexpr bool enable_mem_type<FieldMem<ElementType, IdxRangeType, MemoryType>> = true;
 
-template <class PolarBSplines>
-inline constexpr bool enable_data_access_methods<PolarSpline<PolarBSplines>> = true;
+template <class PolarBSplines, class MemorySpace>
+inline constexpr bool enable_data_access_methods<PolarSpline<PolarBSplines, MemorySpace>> = true;
 
-template <class PolarBSplines>
-inline constexpr bool enable_data_access_methods<PolarSplineView<PolarBSplines>> = true;
+template <class PolarBSplines, class MemorySpace>
+inline constexpr bool
+        enable_data_access_methods<PolarSplineView<PolarBSplines, MemorySpace>> = true;
 
-template <class PolarBSplines>
-inline constexpr bool enable_data_access_methods<PolarSplineSpan<PolarBSplines>> = true;
+template <class PolarBSplines, class MemorySpace>
+inline constexpr bool
+        enable_data_access_methods<PolarSplineSpan<PolarBSplines, MemorySpace>> = true;
 
-template <class PolarBSplines>
-inline constexpr bool enable_mem_type<PolarSpline<PolarBSplines>> = true;
+template <class PolarBSplines, class MemorySpace>
+inline constexpr bool enable_mem_type<PolarSpline<PolarBSplines, MemorySpace>> = true;
 
 template <typename Type>
 static constexpr bool has_idx_range_v = detail::HasIdxRange<Type>::value;
diff --git a/src/utils/ddc_helper.hpp b/src/utils/ddc_helper.hpp
index dd06ba9b1..0be138765 100644
--- a/src/utils/ddc_helper.hpp
+++ b/src/utils/ddc_helper.hpp
@@ -7,6 +7,8 @@
 
 #include <ddc/ddc.hpp>
 
+#include <sll/polar_spline.hpp>
+
 #include "ddc_aliases.hpp"
 #include "directional_tag.hpp"
 #include "transpose.hpp"
@@ -313,6 +315,20 @@ struct OnMemorySpace<
     using type = VectorField<ElementType, SupportType, NDTag, NewMemorySpace, Layout>;
 };
 
+/**
+ * @brief Get a new `PolarSpline` type with the same parametrisation
+ * except in the memory space which is set to NewMemorySpace.
+ * @tparam NewMemorySpace The new memory space. 
+ * @tparam PolarSplineType Type of the B-splines.
+ * @tparam MemorySpace The original memory space of the chunk of the VectorFieldMem.
+ * @see VectorField
+ */
+template <class NewMemorySpace, class PolarSplineType, class MemorySpace>
+struct OnMemorySpace<NewMemorySpace, PolarSpline<PolarSplineType, MemorySpace>>
+{
+    using type = PolarSpline<PolarSplineType, NewMemorySpace>;
+};
+
 
 template <template <class Tag> class Templ, class TypeSeq>
 struct ApplyTemplateToTypeSeq;
diff --git a/tests/geometryRTheta/advection_field_rp/advection_field_gtest.cpp b/tests/geometryRTheta/advection_field_rp/advection_field_gtest.cpp
index d33ccc275..db3c8e9d5 100644
--- a/tests/geometryRTheta/advection_field_rp/advection_field_gtest.cpp
+++ b/tests/geometryRTheta/advection_field_rp/advection_field_gtest.cpp
@@ -112,8 +112,8 @@ TEST(AdvectionFieldRThetaComputation, TestAdvectionFieldFinder)
 
 
     ddc::NullExtrapolationRule r_extrapolation_rule;
-    PolarSplineEvaluator<PolarBSplinesRTheta, ddc::NullExtrapolationRule> polar_spline_evaluator(
-            r_extrapolation_rule);
+    PolarSplineEvaluator<PolarBSplinesRTheta, ddc::NullExtrapolationRule, Kokkos::HostSpace>
+            polar_spline_evaluator(r_extrapolation_rule);
 
     // --- Define the mapping. ------------------------------------------------------------------------
     const Mapping mapping;
diff --git a/vendor/sll/include/sll/polar_bsplines.hpp b/vendor/sll/include/sll/polar_bsplines.hpp
index 54bdbb20b..35a9f6e39 100644
--- a/vendor/sll/include/sll/polar_bsplines.hpp
+++ b/vendor/sll/include/sll/polar_bsplines.hpp
@@ -490,7 +490,7 @@ class PolarBSplines
          *
          * @param[out] int_vals The integrals of the basis splines.
          */
-        void integrals(PolarSplineSpan<DDim> int_vals) const;
+        void integrals(PolarSplineSpan<DDim, MemorySpace> int_vals) const;
 
         /**
          * Get the total number of basis functions.
@@ -662,7 +662,7 @@ ddc::DiscreteElement<BSplinesR, BSplinesTheta> PolarBSplines<BSplinesR, BSplines
 template <class BSplinesR, class BSplinesTheta, int C>
 template <class DDim, class MemorySpace>
 void PolarBSplines<BSplinesR, BSplinesTheta, C>::Impl<DDim, MemorySpace>::integrals(
-        PolarSplineSpan<DDim> int_vals) const
+        PolarSplineSpan<DDim, MemorySpace> int_vals) const
 {
     auto r_bspl_space = ddc::discrete_space<BSplinesR>();
     auto theta_bspl_space = ddc::discrete_space<BSplinesTheta>();
@@ -676,11 +676,14 @@ void PolarBSplines<BSplinesR, BSplinesTheta, C>::Impl<DDim, MemorySpace>::integr
            || int_vals.spline_coef.domain().template extent<BSplinesTheta>()
                       == theta_bspl_space.size());
 
-    ddc::Chunk<double, ddc::DiscreteDomain<BSplinesR>> r_integrals_alloc(
-            r_bspl_space.full_domain().take_first(
+    ddc::Chunk<double, ddc::DiscreteDomain<BSplinesR>, ddc::KokkosAllocator<double, MemorySpace>>
+            r_integrals_alloc(r_bspl_space.full_domain().take_first(
                     ddc::DiscreteVector<BSplinesR> {r_bspl_space.nbasis()}));
-    ddc::Chunk<double, ddc::DiscreteDomain<BSplinesTheta>> theta_integrals_alloc(
-            theta_bspl_space.full_domain().take_first(
+    ddc::Chunk<
+            double,
+            ddc::DiscreteDomain<BSplinesTheta>,
+            ddc::KokkosAllocator<double, MemorySpace>>
+            theta_integrals_alloc(theta_bspl_space.full_domain().take_first(
                     ddc::DiscreteVector<BSplinesTheta> {theta_bspl_space.size()}));
     ddc::ChunkSpan r_integrals = r_integrals_alloc.span_view();
     ddc::ChunkSpan theta_integrals = theta_integrals_alloc.span_view();
@@ -688,9 +691,13 @@ void PolarBSplines<BSplinesR, BSplinesTheta, C>::Impl<DDim, MemorySpace>::integr
     ddc::integrals(Kokkos::DefaultHostExecutionSpace(), r_integrals);
     ddc::integrals(Kokkos::DefaultHostExecutionSpace(), theta_integrals);
 
-    ddc::for_each(singular_idx_range<DDim>(), [&](auto k) {
-        int_vals.singular_spline_coef(k) = ddc::transform_reduce(
-                ddc::select<BSplinesR, BSplinesTheta>(m_singular_basis_elements.domain()),
+    auto singular_domain
+            = ddc::select<BSplinesR, BSplinesTheta>(m_singular_basis_elements.domain());
+    ddc::ChunkSpan singular_spline_integrals = int_vals.singular_spline_coef.span_view();
+
+    ddc::for_each(singular_idx_range<DDim>(), [&](ddc::DiscreteElement<DDim> k) {
+        singular_spline_integrals(k) = ddc::transform_reduce(
+                singular_domain,
                 0.0,
                 ddc::reducer::sum<double>(),
                 [&](tensor_product_index_type const idx) {
@@ -702,19 +709,20 @@ void PolarBSplines<BSplinesR, BSplinesTheta, C>::Impl<DDim, MemorySpace>::integr
 
     ddc::DiscreteDomain<BSplinesR> r_tensor_product_dom(
             ddc::select<BSplinesR>(int_vals.spline_coef.domain()));
-
     tensor_product_idx_range_type
             tensor_bspline_idx_range(r_tensor_product_dom, theta_integrals.domain());
+    ddc::ChunkSpan spline_integrals = int_vals.spline_coef.span_view();
 
-    ddc::for_each(tensor_bspline_idx_range, [&](auto idx) {
-        int_vals.spline_coef(idx) = r_integrals(ddc::select<BSplinesR>(idx))
-                                    * theta_integrals(ddc::select<BSplinesTheta>(idx));
+    ddc::for_each(tensor_bspline_idx_range, [&](tensor_product_index_type idx) {
+        spline_integrals(idx) = r_integrals(ddc::select<BSplinesR>(idx))
+                                * theta_integrals(ddc::select<BSplinesTheta>(idx));
     });
-
     if (int_vals.spline_coef.domain().template extent<BSplinesTheta>() == theta_bspl_space.size()) {
         ddc::DiscreteDomain<BSplinesTheta> periodic_points(theta_integrals.domain().take_last(
                 ddc::DiscreteVector<BSplinesTheta> {BSplinesTheta::degree()}));
         tensor_product_idx_range_type repeat_idx_range(r_tensor_product_dom, periodic_points);
-        ddc::for_each(repeat_idx_range, [&](auto idx) { int_vals.spline_coef(idx) = 0.0; });
+        ddc::for_each(repeat_idx_range, [&](tensor_product_index_type idx) {
+            spline_integrals(idx) = 0.0;
+        });
     }
 }
diff --git a/vendor/sll/include/sll/polar_spline.hpp b/vendor/sll/include/sll/polar_spline.hpp
index c7d079d47..bf79b4236 100644
--- a/vendor/sll/include/sll/polar_spline.hpp
+++ b/vendor/sll/include/sll/polar_spline.hpp
@@ -2,10 +2,10 @@
 
 #include <ddc/ddc.hpp>
 
-template <class PolarBSplinesType>
+template <class PolarBSplinesType, class MemSpace = Kokkos::DefaultExecutionSpace::memory_space>
 struct PolarSplineSpan;
 
-template <class PolarBSplinesType>
+template <class PolarBSplinesType, class MemSpace = Kokkos::DefaultExecutionSpace::memory_space>
 struct PolarSplineView;
 
 
@@ -15,7 +15,7 @@ struct PolarSplineView;
  *
  * @tparam PolarBSplinesType The type of the polar bsplines on which this spline is defined.
  */
-template <class PolarBSplinesType>
+template <class PolarBSplinesType, class MemSpace = Kokkos::DefaultExecutionSpace::memory_space>
 struct PolarSpline
 {
 public:
@@ -29,13 +29,21 @@ struct PolarSpline
      * A Chunk containing the coefficients in front of the b-spline elements which can be
      * expressed as a tensor product of 1D bsplines.
      */
-    ddc::Chunk<double, ddc::DiscreteDomain<BSplinesR, BSplinesTheta>> spline_coef;
+    ddc::Chunk<
+            double,
+            ddc::DiscreteDomain<BSplinesR, BSplinesTheta>,
+            ddc::KokkosAllocator<double, MemSpace>>
+            spline_coef;
 
     /**
      * A Chunk containing the coefficients in front of the b-spline elements near the
      * singular point which cannot be expressed as a tensor product of 1D bsplines.
      */
-    ddc::Chunk<double, ddc::DiscreteDomain<PolarBSplinesType>> singular_spline_coef;
+    ddc::Chunk<
+            double,
+            ddc::DiscreteDomain<PolarBSplinesType>,
+            ddc::KokkosAllocator<double, MemSpace>>
+            singular_spline_coef;
 
 public:
     /**
@@ -47,7 +55,7 @@ struct PolarSpline
      *
      * @param domain A 2D domain of bsplines.
      */
-    PolarSpline<PolarBSplinesType>(ddc::DiscreteDomain<BSplinesR, BSplinesTheta> domain)
+    PolarSpline<PolarBSplinesType, MemSpace>(ddc::DiscreteDomain<BSplinesR, BSplinesTheta> domain)
         : spline_coef(ddc::DiscreteDomain<BSplinesR, BSplinesTheta>(
                 ddc::select<BSplinesR>(domain).remove_first(
                         ddc::DiscreteVector<BSplinesR>(PolarBSplinesType::continuity + 1)),
@@ -67,7 +75,7 @@ struct PolarSpline
      * @param domain          The domain for the coefficients in front of the b-spline elements
      *                        which can be expressed as a tensor product of 1D bsplines.
      */
-    PolarSpline<PolarBSplinesType>(
+    PolarSpline<PolarBSplinesType, MemSpace>(
             ddc::DiscreteDomain<PolarBSplinesType> singular_domain,
             ddc::DiscreteDomain<BSplinesR, BSplinesTheta> domain)
         : spline_coef(domain)
@@ -81,9 +89,9 @@ struct PolarSpline
      *
      * @return A modifiable reference to this polar spline.
      */
-    PolarSplineSpan<PolarBSplinesType> span_view()
+    PolarSplineSpan<PolarBSplinesType, MemSpace> span_view()
     {
-        return PolarSplineSpan<PolarBSplinesType>(*this);
+        return PolarSplineSpan<PolarBSplinesType, MemSpace>(*this);
     }
 
     /**
@@ -91,9 +99,9 @@ struct PolarSpline
      *
      * @return A constant reference to this polar spline.
      */
-    PolarSplineView<PolarBSplinesType> span_cview() const
+    PolarSplineView<PolarBSplinesType, MemSpace> span_cview() const
     {
-        return PolarSplineView<PolarBSplinesType>(*this);
+        return PolarSplineView<PolarBSplinesType, MemSpace>(*this);
     }
 };
 
@@ -103,7 +111,7 @@ struct PolarSpline
  *
  * @tparam PolarBSplinesType The type of the polar bsplines on which this spline is defined.
  */
-template <class PolarBSplinesType>
+template <class PolarBSplinesType, class MemSpace>
 struct PolarSplineSpan
 {
 public:
@@ -117,13 +125,23 @@ struct PolarSplineSpan
      * A ChunkSpan containing the coefficients in front of the b-spline elements which can be
      * expressed as a tensor product of 1D bsplines.
      */
-    ddc::ChunkSpan<double, ddc::DiscreteDomain<BSplinesR, BSplinesTheta>> spline_coef;
+    ddc::ChunkSpan<
+            double,
+            ddc::DiscreteDomain<BSplinesR, BSplinesTheta>,
+            std::experimental::layout_right,
+            MemSpace>
+            spline_coef;
 
     /**
      * A ChunkSpan containing the coefficients in front of the b-spline elements near the
      * singular point which cannot be expressed as a tensor product of 1D bsplines.
      */
-    ddc::ChunkSpan<double, ddc::DiscreteDomain<PolarBSplinesType>> singular_spline_coef;
+    ddc::ChunkSpan<
+            double,
+            ddc::DiscreteDomain<PolarBSplinesType>,
+            std::experimental::layout_right,
+            MemSpace>
+            singular_spline_coef;
 
 public:
     /**
@@ -131,7 +149,7 @@ struct PolarSplineSpan
      *
      * @param spl The PolarSpline being referenced.
      */
-    PolarSplineSpan<PolarBSplinesType>(PolarSpline<PolarBSplinesType>& spl)
+    PolarSplineSpan<PolarBSplinesType, MemSpace>(PolarSpline<PolarBSplinesType, MemSpace>& spl)
         : spline_coef(spl.spline_coef.span_view())
         , singular_spline_coef(spl.singular_spline_coef.span_view())
     {
@@ -142,7 +160,7 @@ struct PolarSplineSpan
      *
      * @return A modifiable reference to a polar spline.
      */
-    PolarSplineSpan<PolarBSplinesType> span_view()
+    PolarSplineSpan<PolarBSplinesType, MemSpace> span_view()
     {
         return *this;
     }
@@ -152,9 +170,9 @@ struct PolarSplineSpan
      *
      * @return A constant reference to a polar spline.
      */
-    PolarSplineView<PolarBSplinesType> span_cview() const
+    PolarSplineView<PolarBSplinesType, MemSpace> span_cview() const
     {
-        return PolarSplineView<PolarBSplinesType>(*this);
+        return PolarSplineView<PolarBSplinesType, MemSpace>(*this);
     }
 };
 
@@ -164,7 +182,7 @@ struct PolarSplineSpan
  *
  * @tparam PolarBSplinesType The type of the polar bsplines on which this spline is defined.
  */
-template <class PolarBSplinesType>
+template <class PolarBSplinesType, class MemSpace>
 struct PolarSplineView
 {
 public:
@@ -178,13 +196,21 @@ struct PolarSplineView
      * A ChunkView containing the coefficients in front of the b-spline elements which can be
      * expressed as a tensor product of 1D bsplines.
      */
-    ddc::ChunkSpan<double const, ddc::DiscreteDomain<BSplinesR, BSplinesTheta>> const spline_coef;
+    ddc::ChunkSpan<
+            double const,
+            ddc::DiscreteDomain<BSplinesR, BSplinesTheta>,
+            std::experimental::layout_right,
+            MemSpace> const spline_coef;
 
     /**
      * A ChunkView containing the coefficients in front of the b-spline elements near the
      * singular point which cannot be expressed as a tensor product of 1D bsplines.
      */
-    ddc::ChunkSpan<double const, ddc::DiscreteDomain<PolarBSplinesType>> const singular_spline_coef;
+    ddc::ChunkSpan<
+            double const,
+            ddc::DiscreteDomain<PolarBSplinesType>,
+            std::experimental::layout_right,
+            MemSpace> const singular_spline_coef;
 
 public:
     /**
@@ -192,7 +218,8 @@ struct PolarSplineView
      *
      * @param spl The PolarSpline being referenced.
      */
-    PolarSplineView<PolarBSplinesType>(PolarSpline<PolarBSplinesType> const& spl)
+    PolarSplineView<PolarBSplinesType, MemSpace>(
+            PolarSpline<PolarBSplinesType, MemSpace> const& spl)
         : spline_coef(spl.spline_coef.span_cview())
         , singular_spline_coef(spl.singular_spline_coef.span_cview())
     {
@@ -203,7 +230,8 @@ struct PolarSplineView
      *
      * @param spl The PolarSpline being referenced.
      */
-    PolarSplineView<PolarBSplinesType>(PolarSplineSpan<PolarBSplinesType> const& spl)
+    PolarSplineView<PolarBSplinesType, MemSpace>(
+            PolarSplineSpan<PolarBSplinesType, MemSpace> const& spl)
         : spline_coef(spl.spline_coef.span_cview())
         , singular_spline_coef(spl.singular_spline_coef.span_cview())
     {
@@ -214,7 +242,7 @@ struct PolarSplineView
      *
      * @return A reference to a polar spline.
      */
-    PolarSplineSpan<PolarBSplinesType> span_view() const
+    PolarSplineSpan<PolarBSplinesType, MemSpace> span_view() const
     {
         return *this;
     }
@@ -224,7 +252,7 @@ struct PolarSplineView
      *
      * @return A constant reference to a polar spline.
      */
-    PolarSplineView<PolarBSplinesType> span_cview() const
+    PolarSplineView<PolarBSplinesType, MemSpace> span_cview() const
     {
         return *this;
     }
@@ -233,11 +261,11 @@ struct PolarSplineView
 template <class T>
 inline constexpr bool is_polar_spline_v = false;
 
-template <class PolarBSplinesType>
-inline constexpr bool is_polar_spline_v<PolarSpline<PolarBSplinesType>> = true;
+template <class PolarBSplinesType, class MemSpace>
+inline constexpr bool is_polar_spline_v<PolarSpline<PolarBSplinesType, MemSpace>> = true;
 
-template <class PolarBSplinesType>
-inline constexpr bool is_polar_spline_v<PolarSplineSpan<PolarBSplinesType>> = true;
+template <class PolarBSplinesType, class MemSpace>
+inline constexpr bool is_polar_spline_v<PolarSplineSpan<PolarBSplinesType, MemSpace>> = true;
 
-template <class PolarBSplinesType>
-inline constexpr bool is_polar_spline_v<PolarSplineView<PolarBSplinesType>> = true;
+template <class PolarBSplinesType, class MemSpace>
+inline constexpr bool is_polar_spline_v<PolarSplineView<PolarBSplinesType, MemSpace>> = true;
diff --git a/vendor/sll/include/sll/polar_spline_evaluator.hpp b/vendor/sll/include/sll/polar_spline_evaluator.hpp
index b37dadb70..fe7f3f66b 100644
--- a/vendor/sll/include/sll/polar_spline_evaluator.hpp
+++ b/vendor/sll/include/sll/polar_spline_evaluator.hpp
@@ -8,7 +8,7 @@
  *
  * @see PolarBSplines
  */
-template <class PolarBSplinesType, class OuterExtrapolationRule>
+template <class PolarBSplinesType, class OuterExtrapolationRule, class MemSpace>
 class PolarSplineEvaluator
 {
 private:
@@ -141,7 +141,7 @@ class PolarSplineEvaluator
      */
     double operator()(
             ddc::Coordinate<DimR, DimTheta> coord_eval,
-            PolarSplineView<PolarBSplinesType> const spline_coef) const
+            PolarSplineView<PolarBSplinesType, MemSpace> const spline_coef) const
     {
         return eval(coord_eval, spline_coef);
     }
@@ -158,9 +158,14 @@ class PolarSplineEvaluator
      */
     template <class Domain>
     void operator()(
-            ddc::ChunkSpan<double, Domain> const spline_eval,
-            ddc::ChunkSpan<ddc::Coordinate<DimR, DimTheta> const, Domain> const coords_eval,
-            PolarSplineView<PolarBSplinesType> const spline_coef) const
+            ddc::ChunkSpan<double, Domain, std::experimental::layout_right, MemSpace> const
+                    spline_eval,
+            ddc::ChunkSpan<
+                    ddc::Coordinate<DimR, DimTheta> const,
+                    Domain,
+                    std::experimental::layout_right,
+                    MemSpace> const coords_eval,
+            PolarSplineView<PolarBSplinesType, MemSpace> const spline_coef) const
     {
         ddc::for_each(coords_eval.domain(), [=](auto i) {
             spline_eval(i) = eval(coords_eval(i), spline_coef);
@@ -181,7 +186,7 @@ class PolarSplineEvaluator
      */
     double deriv_dim_1(
             ddc::Coordinate<DimR, DimTheta> coord_eval,
-            PolarSplineView<PolarBSplinesType> const spline_coef) const
+            PolarSplineView<PolarBSplinesType, MemSpace> const spline_coef) const
     {
         return eval_no_bc(coord_eval, spline_coef, eval_deriv_r_type());
     }
@@ -200,7 +205,7 @@ class PolarSplineEvaluator
      */
     double deriv_dim_2(
             ddc::Coordinate<DimR, DimTheta> coord_eval,
-            PolarSplineView<PolarBSplinesType> const spline_coef) const
+            PolarSplineView<PolarBSplinesType, MemSpace> const spline_coef) const
     {
         return eval_no_bc(coord_eval, spline_coef, eval_deriv_theta_type());
     }
@@ -219,7 +224,7 @@ class PolarSplineEvaluator
      */
     double deriv_1_and_2(
             ddc::Coordinate<DimR, DimTheta> coord_eval,
-            PolarSplineView<PolarBSplinesType> const spline_coef) const
+            PolarSplineView<PolarBSplinesType, MemSpace> const spline_coef) const
     {
         return eval_no_bc(coord_eval, spline_coef, eval_deriv_r_theta_type());
     }
@@ -239,7 +244,7 @@ class PolarSplineEvaluator
     void deriv_dim_1(
             ddc::ChunkSpan<double, Domain> const spline_eval,
             ddc::ChunkSpan<ddc::Coordinate<DimR, DimTheta> const, Domain> const coords_eval,
-            PolarSplineView<PolarBSplinesType> const spline_coef) const
+            PolarSplineView<PolarBSplinesType, MemSpace> const spline_coef) const
     {
         ddc::for_each(coords_eval.domain(), [=](auto i) {
             spline_eval(i) = eval_no_bc(coords_eval(i), spline_coef, eval_deriv_r_type());
@@ -261,7 +266,7 @@ class PolarSplineEvaluator
     void deriv_dim_2(
             ddc::ChunkSpan<double, Domain> const spline_eval,
             ddc::ChunkSpan<ddc::Coordinate<DimR, DimTheta> const, Domain> const coords_eval,
-            PolarSplineView<PolarBSplinesType> const spline_coef) const
+            PolarSplineView<PolarBSplinesType, MemSpace> const spline_coef) const
     {
         ddc::for_each(coords_eval.domain(), [=](auto i) {
             spline_eval(i) = eval_no_bc(coords_eval(i), spline_coef, eval_deriv_theta_type());
@@ -282,7 +287,7 @@ class PolarSplineEvaluator
     void deriv_dim_1_and_2(
             ddc::ChunkSpan<double, Domain> const spline_eval,
             ddc::ChunkSpan<ddc::Coordinate<DimR, DimTheta> const, Domain> const coords_eval,
-            PolarSplineView<PolarBSplinesType> const spline_coef) const
+            PolarSplineView<PolarBSplinesType, MemSpace> const spline_coef) const
     {
         ddc::for_each(coords_eval.domain(), [=](auto i) {
             spline_eval(i) = eval_no_bc(coords_eval(i), spline_coef, eval_deriv_r_theta_type());
@@ -300,8 +305,9 @@ class PolarSplineEvaluator
      * @return The integral of the spline function over the domain.
      */
     template <class Mapping>
-    double integrate(PolarSplineView<PolarBSplinesType> const spline_coef, Mapping const mapping)
-            const
+    double integrate(
+            PolarSplineView<PolarBSplinesType, MemSpace> const spline_coef,
+            Mapping const mapping) const
     {
         int constexpr nr = ddc::discrete_space<BSplinesR>().ncells() + BSplinesR::degree() - 2;
         int constexpr ntheta
@@ -333,7 +339,7 @@ class PolarSplineEvaluator
 private:
     double eval(
             ddc::Coordinate<DimR, DimTheta> coord_eval,
-            PolarSplineView<PolarBSplinesType> const spline_coef) const
+            PolarSplineView<PolarBSplinesType, MemSpace> const spline_coef) const
     {
         const double coord_eval1 = ddc::get<DimR>(coord_eval);
         double coord_eval2 = ddc::get<DimTheta>(coord_eval);
@@ -353,7 +359,7 @@ class PolarSplineEvaluator
     template <class EvalType>
     double eval_no_bc(
             ddc::Coordinate<DimR, DimTheta> coord_eval,
-            PolarSplineView<PolarBSplinesType> const spline_coef,
+            PolarSplineView<PolarBSplinesType, MemSpace> const spline_coef,
             EvalType const) const
     {
         static_assert(
diff --git a/vendor/sll/tests/polar_splines.cpp b/vendor/sll/tests/polar_splines.cpp
index cc79d86cc..420d8ebab 100644
--- a/vendor/sll/tests/polar_splines.cpp
+++ b/vendor/sll/tests/polar_splines.cpp
@@ -77,11 +77,11 @@ TEST(PolarSplineTest, ConstantEval)
     using PolarCoord = ddc::Coordinate<R, Theta>;
     using CoordR = ddc::Coordinate<R>;
     using CoordTheta = ddc::Coordinate<Theta>;
-    using Spline = PolarSpline<BSplines>;
-    using Evaluator = PolarSplineEvaluator<BSplines, ddc::NullExtrapolationRule>;
+    using Spline = PolarSpline<BSplines, Kokkos::HostSpace>;
+    using Evaluator = PolarSplineEvaluator<BSplines, ddc::NullExtrapolationRule, Kokkos::HostSpace>;
     using BuilderRTheta = ddc::SplineBuilder2D<
             Kokkos::DefaultHostExecutionSpace,
-            Kokkos::DefaultHostExecutionSpace::memory_space,
+            Kokkos::HostSpace,
             BSplinesR,
             BSplinesTheta,
             GridR,
@@ -96,7 +96,7 @@ TEST(PolarSplineTest, ConstantEval)
 
     using EvaluatorRTheta = ddc::SplineEvaluator2D<
             Kokkos::DefaultHostExecutionSpace,
-            Kokkos::DefaultHostExecutionSpace::memory_space,
+            Kokkos::HostSpace,
             BSplinesR,
             BSplinesTheta,
             GridR,
@@ -199,7 +199,7 @@ TEST(PolarSplineTest, ConstantEval)
     }
 
     Spline integrals(builder_rtheta.spline_domain());
-    ddc::discrete_space<BSplines>().integrals(integrals);
+    ddc::host_discrete_space<BSplines>().integrals(integrals.span_view());
     double area = ddc::transform_reduce(
                           integrals.singular_spline_coef.domain(),
                           0.0,