WIP: Piecewise cubic monotone interpolation

[BenWibking]

PR Summary

Adds monotone interpolation using cubic Hermite polynomials that can be used on device. However, the MonotoneInterpolator itself must be constructed on the host, with the contents of the views passed to it accessible from both host and device (if used).

This depends on Kokkos 4 support in order to use SharedHostPinnedSpace (kokkos/kokkos#5405).

Reference: "Monotone piecewise cubic interpolation," SIAM J. Numer. Anal., Vol. 17, No. 2, April 1980. https://doi.org/10.1137/0717021

PR Checklist

• Code passes cpplint
• New features are documented.
• Adds a test for any bugs fixed. Adds tests for new features.
• Code is formatted
• Changes are summarized in CHANGELOG.md
• CI has been triggered on Darwin for performance regression tests.
• Docs build
• (@lanl.gov employees) Update copyright on changed files

[brryan]

This is neat! Two quick thoughts now (mostly because I saw the include guard and got to thinking about that more generally)

[Yurlungur]

I think this is useful overall... I assuem a fleshed out PR would have tests. My one big concern is you show using it for both std::vector and a kokkos view... those types have different kinds of semantics: vectors are value semantics and views reference semantics. I worry that could get one in to trouble... but maybe it's fine.

[Yurlungur]

Are these assumed to have reference semantics?

[BenWibking]

I don't think so... Can you provide an example of how this could cause trouble?

[Yurlungur]

you're constructing x_vec_ by copy when you call the interpolator, so you're copying around the data unless VectorContainer has refernece semantics. That may be a feature or a bug. If VectoContainer is always a Kokkos::View, then we're fine though.

[Yurlungur]

Linear search seems like a bad idea IMO. Probably better to just branch.

[BenWibking]

I can replace this with a binary search.

[Yurlungur]

I wonder if this can be generalized using @lroberts36 concepts so all view types could be specialized for at once.

[BenWibking]

That would definitely be nice to have.

[BenWibking]

Currently PinnedArray1D is a typedef for SharedHostPinnedSpace in Kokkos 4, so this adds a requirement for Kokkos 4 if we have to explicitly specialize for common vector containers/views.

[Yurlungur]

I think this should just work, but Luke should comment

template<typename VectorContainer, class = ENABLEIF(implements<kokkos_view(VectorContainer)>::value)>
auto MonotoneInterpolator<VectorContainer>::ConstructVectorContainer(size_t size) {
  return VectorContainer("d", size);
}

you'll need to include utils/concepts_lite.hpp

[lroberts36]

I think that should work, but it is probably more consistent to use REQUIRES instead

template<typename VectorContainer, REQUIRES(implements<kokkos_view(VectorContainer)>::value)>
auto MonotoneInterpolator<VectorContainer>::ConstructVectorContainer(size_t size) {
  return VectorContainer("d", size);
}

[BenWibking]

I can't get this to compile on my machine. I don't fully understand concepts, so I'm not sure why it's complaining:

/Users/benwibking/parthenon/src/utils/interpolate.hpp:54:45: error: nested name specifier 'MonotoneInterpolator<VectorContainer>::' for declaration does not refer into a class, class template or class template partial specialization
auto MonotoneInterpolator<VectorContainer>::ConstructVectorContainer(size_t size) {
     ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

[BenWibking]

@Yurlungur Yes, planning to add tests.

[Yurlungur]

Tests look good, but it would be good to also test on device if you plan to use this tool on GPUs.

[Yurlungur]

you're constructing x_vec_ by copy when you call the interpolator, so you're copying around the data unless VectorContainer has refernece semantics. That may be a feature or a bug. If VectoContainer is always a Kokkos::View, then we're fine though.