Skip to content

Commit

Permalink
Alternative implementation to workaround llvm <17
Browse files Browse the repository at this point in the history
  • Loading branch information
tpadioleau committed Nov 2, 2024
1 parent 81cf2db commit bb724ee
Show file tree
Hide file tree
Showing 7 changed files with 64 additions and 75 deletions.
4 changes: 2 additions & 2 deletions include/ddc/chunk.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -193,7 +193,7 @@ class Chunk<ElementType, DiscreteDomain<DDims...>, Allocator>
static_assert((is_discrete_element_v<DElems> && ...), "Expected DiscreteElements");
assert(((select<DDims>(take<DDims>(delems...)) >= front<DDims>(this->m_domain)) && ...));
assert(((select<DDims>(take<DDims>(delems...)) <= back<DDims>(this->m_domain)) && ...));
return detail::mdspan_call(this->m_internal_mdspan, uid<DDims>(take<DDims>(delems...))...);
return DDC_MDSPAN_ACCESS_OP(this->m_internal_mdspan, uid<DDims>(take<DDims>(delems...))...);
}

/** Element access using a list of DiscreteElement
Expand All @@ -209,7 +209,7 @@ class Chunk<ElementType, DiscreteDomain<DDims...>, Allocator>
static_assert((is_discrete_element_v<DElems> && ...), "Expected DiscreteElements");
assert(((select<DDims>(take<DDims>(delems...)) >= front<DDims>(this->m_domain)) && ...));
assert(((select<DDims>(take<DDims>(delems...)) <= back<DDims>(this->m_domain)) && ...));
return detail::mdspan_call(this->m_internal_mdspan, uid<DDims>(take<DDims>(delems...))...);
return DDC_MDSPAN_ACCESS_OP(this->m_internal_mdspan, uid<DDims>(take<DDims>(delems...))...);
}

/** Returns the label of the Chunk
Expand Down
2 changes: 1 addition & 1 deletion include/ddc/chunk_span.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -330,7 +330,7 @@ class ChunkSpan<ElementType, DiscreteDomain<DDims...>, LayoutStridedPolicy, Memo
static_assert((is_discrete_element_v<DElems> && ...), "Expected DiscreteElements");
assert(((select<DDims>(take<DDims>(delems...)) >= front<DDims>(this->m_domain)) && ...));
assert(((select<DDims>(take<DDims>(delems...)) <= back<DDims>(this->m_domain)) && ...));
return detail::mdspan_call(this->m_internal_mdspan, uid<DDims>(take<DDims>(delems...))...);
return DDC_MDSPAN_ACCESS_OP(this->m_internal_mdspan, uid<DDims>(take<DDims>(delems...))...);
}

/** Access to the underlying allocation pointer
Expand Down
25 changes: 7 additions & 18 deletions include/ddc/detail/kokkos.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -13,28 +13,17 @@

#include "macros.hpp"

namespace ddc::detail {

template <
class ElementType,
class Extents,
class LayoutPolicy,
class AccessorPolicy,
class... SizeTypes>
KOKKOS_FORCEINLINE_FUNCTION constexpr decltype(auto) mdspan_call(
Kokkos::mdspan<ElementType, Extents, LayoutPolicy, AccessorPolicy> const& mdspan,
SizeTypes const&... indices)
{
#if defined(MDSPAN_USE_BRACKET_OPERATOR)
#if MDSPAN_USE_BRACKET_OPERATOR
return mdspan[indices...];
#define DDC_MDSPAN_ACCESS_OP(mds, ...) mds[__VA_ARGS__]
#else
return mdspan(indices...);
#define DDC_MDSPAN_ACCESS_OP(mds, ...) mds(__VA_ARGS__)
#endif
#else
return mdspan[indices...];
#define DDC_MDSPAN_ACCESS_OP(mds, ...) mds[__VA_ARGS__]
#endif
}

namespace ddc::detail {

template <class T>
struct type_holder
Expand Down Expand Up @@ -154,8 +143,8 @@ KOKKOS_FUNCTION mdspan_to_kokkos_layout_t<typename MP::layout_type> build_kokkos
std::size_t,
Kokkos::extents<std::size_t, sizeof...(Is), 2>,
Kokkos::layout_right> const interleaved_extents_strides(storage.data());
((mdspan_call(interleaved_extents_strides, Is, 0) = ep.extent(Is),
mdspan_call(interleaved_extents_strides, Is, 1) = mapping.stride(Is)),
((DDC_MDSPAN_ACCESS_OP(interleaved_extents_strides, Is, 0) = ep.extent(Is),
DDC_MDSPAN_ACCESS_OP(interleaved_extents_strides, Is, 1) = mapping.stride(Is)),
...);
return make_layout_stride(storage, std::make_index_sequence<sizeof...(Is) * 2> {});
} else {
Expand Down
38 changes: 19 additions & 19 deletions include/ddc/kernels/splines/bsplines_non_uniform.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -600,63 +600,63 @@ KOKKOS_INLINE_FUNCTION ddc::DiscreteElement<DDim> NonUniformBSplines<CDim, D>::

double saved;
double temp;
detail::mdspan_call(ndu, 0, 0) = 1.0;
DDC_MDSPAN_ACCESS_OP(ndu, 0, 0) = 1.0;
for (std::size_t j = 0; j < degree(); ++j) {
left[j] = x - ddc::coordinate(icell - j);
right[j] = ddc::coordinate(icell + j + 1) - x;
saved = 0.0;
for (std::size_t r = 0; r < j + 1; ++r) {
// compute inverse of knot differences and save them into lower
// triangular part of ndu
detail::mdspan_call(ndu, r, j + 1) = 1.0 / (right[r] + left[j - r]);
DDC_MDSPAN_ACCESS_OP(ndu, r, j + 1) = 1.0 / (right[r] + left[j - r]);
// compute basis functions and save them into upper triangular part
// of ndu
temp = detail::mdspan_call(ndu, j, r) * detail::mdspan_call(ndu, r, j + 1);
detail::mdspan_call(ndu, j + 1, r) = saved + right[r] * temp;
temp = DDC_MDSPAN_ACCESS_OP(ndu, j, r) * DDC_MDSPAN_ACCESS_OP(ndu, r, j + 1);
DDC_MDSPAN_ACCESS_OP(ndu, j + 1, r) = saved + right[r] * temp;
saved = left[j - r] * temp;
}
detail::mdspan_call(ndu, j + 1, j + 1) = saved;
DDC_MDSPAN_ACCESS_OP(ndu, j + 1, j + 1) = saved;
}
// Save 0-th derivative
for (std::size_t j = 0; j < degree() + 1; ++j) {
detail::mdspan_call(derivs, j, 0) = detail::mdspan_call(ndu, degree(), j);
DDC_MDSPAN_ACCESS_OP(derivs, j, 0) = DDC_MDSPAN_ACCESS_OP(ndu, degree(), j);
}

for (int r = 0; r < int(degree() + 1); ++r) {
int s1 = 0;
int s2 = 1;
detail::mdspan_call(a, 0, 0) = 1.0;
DDC_MDSPAN_ACCESS_OP(a, 0, 0) = 1.0;
for (int k = 1; k < int(n + 1); ++k) {
double d = 0.0;
int const rk = r - k;
int const pk = degree() - k;
if (r >= k) {
detail::mdspan_call(a, 0, s2)
= detail::mdspan_call(a, 0, s1) * detail::mdspan_call(ndu, rk, pk + 1);
d = detail::mdspan_call(a, 0, s2) * detail::mdspan_call(ndu, pk, rk);
DDC_MDSPAN_ACCESS_OP(a, 0, s2)
= DDC_MDSPAN_ACCESS_OP(a, 0, s1) * DDC_MDSPAN_ACCESS_OP(ndu, rk, pk + 1);
d = DDC_MDSPAN_ACCESS_OP(a, 0, s2) * DDC_MDSPAN_ACCESS_OP(ndu, pk, rk);
}
int const j1 = rk > -1 ? 1 : (-rk);
int const j2 = (r - 1) <= pk ? k : (degree() - r + 1);
for (int j = j1; j < j2; ++j) {
detail::mdspan_call(a, j, s2)
= (detail::mdspan_call(a, j, s1) - detail::mdspan_call(a, j - 1, s1))
* detail::mdspan_call(ndu, rk + j, pk + 1);
d += detail::mdspan_call(a, j, s2) * detail::mdspan_call(ndu, pk, rk + j);
DDC_MDSPAN_ACCESS_OP(a, j, s2)
= (DDC_MDSPAN_ACCESS_OP(a, j, s1) - DDC_MDSPAN_ACCESS_OP(a, j - 1, s1))
* DDC_MDSPAN_ACCESS_OP(ndu, rk + j, pk + 1);
d += DDC_MDSPAN_ACCESS_OP(a, j, s2) * DDC_MDSPAN_ACCESS_OP(ndu, pk, rk + j);
}
if (r <= pk) {
detail::mdspan_call(a, k, s2)
= -detail::mdspan_call(a, k - 1, s1) * detail::mdspan_call(ndu, r, pk + 1);
d += detail::mdspan_call(a, k, s2) * detail::mdspan_call(ndu, pk, r);
DDC_MDSPAN_ACCESS_OP(a, k, s2) = -DDC_MDSPAN_ACCESS_OP(a, k - 1, s1)
* DDC_MDSPAN_ACCESS_OP(ndu, r, pk + 1);
d += DDC_MDSPAN_ACCESS_OP(a, k, s2) * DDC_MDSPAN_ACCESS_OP(ndu, pk, r);
}
detail::mdspan_call(derivs, r, k) = d;
DDC_MDSPAN_ACCESS_OP(derivs, r, k) = d;
Kokkos::kokkos_swap(s1, s2);
}
}

int r = degree();
for (int k = 1; k < int(n + 1); ++k) {
for (std::size_t i = 0; i < derivs.extent(0); ++i) {
detail::mdspan_call(derivs, i, k) *= r;
DDC_MDSPAN_ACCESS_OP(derivs, i, k) *= r;
}
r *= degree() - k;
}
Expand Down
54 changes: 27 additions & 27 deletions include/ddc/kernels/splines/bsplines_uniform.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -413,17 +413,17 @@ KOKKOS_INLINE_FUNCTION ddc::DiscreteElement<DDim> UniformBSplines<CDim, D>::
double xx;
double temp;
double saved;
detail::mdspan_call(values, 0) = 1.0;
DDC_MDSPAN_ACCESS_OP(values, 0) = 1.0;
for (std::size_t j = 1; j < values.size(); ++j) {
xx = -offset;
saved = 0.0;
for (std::size_t r = 0; r < j; ++r) {
xx += 1;
temp = detail::mdspan_call(values, r) / j;
detail::mdspan_call(values, r) = saved + xx * temp;
temp = DDC_MDSPAN_ACCESS_OP(values, r) / j;
DDC_MDSPAN_ACCESS_OP(values, r) = saved + xx * temp;
saved = (j - xx) * temp;
}
detail::mdspan_call(values, j) = saved;
DDC_MDSPAN_ACCESS_OP(values, j) = saved;
}

return discrete_element_type(jmin);
Expand All @@ -447,29 +447,29 @@ KOKKOS_INLINE_FUNCTION ddc::DiscreteElement<DDim> UniformBSplines<CDim, D>::
double xx;
double temp;
double saved;
detail::mdspan_call(derivs, 0) = 1.0 / ddc::step<knot_discrete_dimension_type>();
DDC_MDSPAN_ACCESS_OP(derivs, 0) = 1.0 / ddc::step<knot_discrete_dimension_type>();
for (std::size_t j = 1; j < degree(); ++j) {
xx = -offset;
saved = 0.0;
for (std::size_t r = 0; r < j; ++r) {
xx += 1.0;
temp = detail::mdspan_call(derivs, r) / j;
detail::mdspan_call(derivs, r) = saved + xx * temp;
temp = DDC_MDSPAN_ACCESS_OP(derivs, r) / j;
DDC_MDSPAN_ACCESS_OP(derivs, r) = saved + xx * temp;
saved = (j - xx) * temp;
}
detail::mdspan_call(derivs, j) = saved;
DDC_MDSPAN_ACCESS_OP(derivs, j) = saved;
}

// Compute derivatives
double bjm1 = derivs[0];
double bj = bjm1;
detail::mdspan_call(derivs, 0) = -bjm1;
DDC_MDSPAN_ACCESS_OP(derivs, 0) = -bjm1;
for (std::size_t j = 1; j < degree(); ++j) {
bj = detail::mdspan_call(derivs, j);
detail::mdspan_call(derivs, j) = bjm1 - bj;
bj = DDC_MDSPAN_ACCESS_OP(derivs, j);
DDC_MDSPAN_ACCESS_OP(derivs, j) = bjm1 - bj;
bjm1 = bj;
}
detail::mdspan_call(derivs, degree()) = bj;
DDC_MDSPAN_ACCESS_OP(derivs, degree()) = bj;

return discrete_element_type(jmin);
}
Expand Down Expand Up @@ -507,47 +507,47 @@ KOKKOS_INLINE_FUNCTION ddc::DiscreteElement<DDim> UniformBSplines<CDim, D>::
double xx;
double temp;
double saved;
detail::mdspan_call(ndu, 0, 0) = 1.0;
DDC_MDSPAN_ACCESS_OP(ndu, 0, 0) = 1.0;
for (std::size_t j = 1; j < degree() + 1; ++j) {
xx = -offset;
saved = 0.0;
for (std::size_t r = 0; r < j; ++r) {
xx += 1.0;
temp = detail::mdspan_call(ndu, j - 1, r) / j;
detail::mdspan_call(ndu, j, r) = saved + xx * temp;
temp = DDC_MDSPAN_ACCESS_OP(ndu, j - 1, r) / j;
DDC_MDSPAN_ACCESS_OP(ndu, j, r) = saved + xx * temp;
saved = (j - xx) * temp;
}
detail::mdspan_call(ndu, j, j) = saved;
DDC_MDSPAN_ACCESS_OP(ndu, j, j) = saved;
}
for (std::size_t i = 0; i < ndu.extent(1); ++i) {
detail::mdspan_call(derivs, i, 0) = detail::mdspan_call(ndu, degree(), i);
DDC_MDSPAN_ACCESS_OP(derivs, i, 0) = DDC_MDSPAN_ACCESS_OP(ndu, degree(), i);
}

for (int r = 0; r < int(degree() + 1); ++r) {
int s1 = 0;
int s2 = 1;
detail::mdspan_call(a, 0, 0) = 1.0;
DDC_MDSPAN_ACCESS_OP(a, 0, 0) = 1.0;
for (int k = 1; k < int(n + 1); ++k) {
double d = 0.0;
int const rk = r - k;
int const pk = degree() - k;
if (r >= k) {
detail::mdspan_call(a, 0, s2) = detail::mdspan_call(a, 0, s1) / (pk + 1);
d = detail::mdspan_call(a, 0, s2) * detail::mdspan_call(ndu, pk, rk);
DDC_MDSPAN_ACCESS_OP(a, 0, s2) = DDC_MDSPAN_ACCESS_OP(a, 0, s1) / (pk + 1);
d = DDC_MDSPAN_ACCESS_OP(a, 0, s2) * DDC_MDSPAN_ACCESS_OP(ndu, pk, rk);
}
int const j1 = rk > -1 ? 1 : (-rk);
int const j2 = (r - 1) <= pk ? k : (degree() - r + 1);
for (int j = j1; j < j2; ++j) {
detail::mdspan_call(a, j, s2)
= (detail::mdspan_call(a, j, s1) - detail::mdspan_call(a, j - 1, s1))
DDC_MDSPAN_ACCESS_OP(a, j, s2)
= (DDC_MDSPAN_ACCESS_OP(a, j, s1) - DDC_MDSPAN_ACCESS_OP(a, j - 1, s1))
/ (pk + 1);
d += detail::mdspan_call(a, j, s2) * detail::mdspan_call(ndu, pk, rk + j);
d += DDC_MDSPAN_ACCESS_OP(a, j, s2) * DDC_MDSPAN_ACCESS_OP(ndu, pk, rk + j);
}
if (r <= pk) {
detail::mdspan_call(a, k, s2) = -detail::mdspan_call(a, k - 1, s1) / (pk + 1);
d += detail::mdspan_call(a, k, s2) * detail::mdspan_call(ndu, pk, r);
DDC_MDSPAN_ACCESS_OP(a, k, s2) = -DDC_MDSPAN_ACCESS_OP(a, k - 1, s1) / (pk + 1);
d += DDC_MDSPAN_ACCESS_OP(a, k, s2) * DDC_MDSPAN_ACCESS_OP(ndu, pk, r);
}
detail::mdspan_call(derivs, r, k) = d;
DDC_MDSPAN_ACCESS_OP(derivs, r, k) = d;
Kokkos::kokkos_swap(s1, s2);
}
}
Expand All @@ -559,7 +559,7 @@ KOKKOS_INLINE_FUNCTION ddc::DiscreteElement<DDim> UniformBSplines<CDim, D>::
double d = degree() * inv_dx;
for (int k = 1; k < int(n + 1); ++k) {
for (std::size_t i = 0; i < derivs.extent(0); ++i) {
detail::mdspan_call(derivs, i, k) *= d;
DDC_MDSPAN_ACCESS_OP(derivs, i, k) *= d;
}
d *= (degree() - k) * inv_dx;
}
Expand Down
12 changes: 6 additions & 6 deletions include/ddc/kernels/splines/spline_builder.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -501,7 +501,7 @@ int SplineBuilder<
} else {
int const mid = bsplines_type::degree() / 2;
offset = jmin.uid() - start.uid()
+ (detail::mdspan_call(values, mid) > detail::mdspan_call(values, mid + 1)
+ (DDC_MDSPAN_ACCESS_OP(values, mid) > DDC_MDSPAN_ACCESS_OP(values, mid + 1)
? mid
: mid + 1)
- BSplines::degree();
Expand Down Expand Up @@ -680,7 +680,7 @@ void SplineBuilder<
// all derivatives by multiplying the i-th derivative by dx^i
for (std::size_t i = 0; i < bsplines_type::degree() + 1; ++i) {
for (std::size_t j = 1; j < bsplines_type::degree() / 2 + 1; ++j) {
detail::mdspan_call(derivs, i, j) *= ddc::detail::ipow(m_dx, j);
DDC_MDSPAN_ACCESS_OP(derivs, i, j) *= ddc::detail::ipow(m_dx, j);
}
}

Expand All @@ -690,7 +690,7 @@ void SplineBuilder<
matrix->set_element(
i,
j,
detail::mdspan_call(derivs, j, s_nbc_xmin - i - 1 + s_odd));
DDC_MDSPAN_ACCESS_OP(derivs, j, s_nbc_xmin - i - 1 + s_odd));
}
}
}
Expand All @@ -709,7 +709,7 @@ void SplineBuilder<
int const j = ddc::detail::
modulo(int(jmin.uid() - m_offset + s),
static_cast<int>(ddc::discrete_space<BSplines>().nbasis()));
matrix->set_element(ix.uid() - start + s_nbc_xmin, j, detail::mdspan_call(values, s));
matrix->set_element(ix.uid() - start + s_nbc_xmin, j, DDC_MDSPAN_ACCESS_OP(values, s));
}
});

Expand All @@ -733,15 +733,15 @@ void SplineBuilder<
// all derivatives by multiplying the i-th derivative by dx^i
for (std::size_t i = 0; i < bsplines_type::degree() + 1; ++i) {
for (std::size_t j = 1; j < bsplines_type::degree() / 2 + 1; ++j) {
detail::mdspan_call(derivs, i, j) *= ddc::detail::ipow(m_dx, j);
DDC_MDSPAN_ACCESS_OP(derivs, i, j) *= ddc::detail::ipow(m_dx, j);
}
}

int const i0 = ddc::discrete_space<BSplines>().nbasis() - s_nbc_xmax;
int const j0 = ddc::discrete_space<BSplines>().nbasis() - bsplines_type::degree();
for (std::size_t j = 0; j < bsplines_type::degree(); ++j) {
for (std::size_t i = 0; i < s_nbc_xmax; ++i) {
matrix->set_element(i0 + i, j0 + j, detail::mdspan_call(derivs, j + 1, i + s_odd));
matrix->set_element(i0 + i, j0 + j, DDC_MDSPAN_ACCESS_OP(derivs, j + 1, i + s_odd));
}
}
}
Expand Down
4 changes: 2 additions & 2 deletions tests/splines/bsplines.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -75,7 +75,7 @@ TYPED_TEST(BSplinesFixture, PartitionOfUnityUniform)
ddc::discrete_space<BSplinesX>().eval_basis(values, test_point);
double sum = 0.0;
for (std::size_t j(0); j < degree + 1; ++j) {
sum += ddc::detail::mdspan_call(values, j);
sum += DDC_MDSPAN_ACCESS_OP(values, j);
}
EXPECT_LE(fabs(sum - 1.0), 1.0e-15);
}
Expand Down Expand Up @@ -110,7 +110,7 @@ TYPED_TEST(BSplinesFixture, PartitionOfUnityNonUniform)
ddc::discrete_space<BSplinesX>().eval_basis(values, test_point);
double sum = 0.0;
for (std::size_t j(0); j < degree + 1; ++j) {
sum += ddc::detail::mdspan_call(values, j);
sum += DDC_MDSPAN_ACCESS_OP(values, j);
}
EXPECT_LE(fabs(sum - 1.0), 1.0e-15);
}
Expand Down

0 comments on commit bb724ee

Please sign in to comment.