Add bi- and trilinear interpolation #1020
Conversation
|
@Yurlungur do you by chance already have a unit test stashed? I didn't find one in Phoebus and so I just want to double check before reinventing the wheel. |
Unfortunately no... I don't think I was careful about unit testing this piece of code. |
src/utils/interpolation.hpp
Outdated
| } // namespace Cent | ||
| } // namespace interpolation | ||
| } // namespace parthenon | ||
| #endif // UTILS_INTERPOLATION_HPP_ |
src/utils/robust.hpp
Outdated
| } | ||
| } // namespace robust | ||
| } // namespace parthenon | ||
| #endif // UTILS_ROBUST_HPP_ |
|
@pgrete tests in spiner: https://github.com/lanl/spiner/blob/main/test/test.cpp |
src/utils/interpolation.hpp
Outdated
| class Linear : public Interpolation { | ||
| public: | ||
| KOKKOS_FUNCTION | ||
| Linear(const int n_support, const Real dx, const Real shift) | ||
| : Interpolation(n_support, dx, shift) { | ||
| PARTHENON_FAIL("Not written yet!"); | ||
| } | ||
|
|
||
| KOKKOS_INLINE_FUNCTION | ||
| void GetIndicesAndWeights(const int i, int *idx, Real *wgt) const override { | ||
| idx[0] = std::floor(i + shift_); | ||
| idx[1] = idx[0] + 1; | ||
|
|
||
| wgt[0] = wgt[1] = 1. - wgt[0]; | ||
|
|
||
| for (int nsup = 0; nsup < 2; nsup++) { | ||
| if (idx[nsup] < 0 || idx[nsup] >= n_support_) { | ||
| idx[nsup] = 0; | ||
| wgt[nsup] = 0.; | ||
| } | ||
| } | ||
| } | ||
|
|
||
| KOKKOS_INLINE_FUNCTION | ||
| int StencilSize() const override { return 2; } | ||
| }; |
There was a problem hiding this comment.
This looks like dead code given the PARTHENON_FAIL, should we just remove it (or possibly all that inherits from Interpolation and Interpolation)? To be fair, it is not really clear to me how this code is supposed to work without seeing an example.
There was a problem hiding this comment.
@Yurlungur what's your take at this (also in light of the following comment
// TODO(JMM): Is this interpolation::Do syntax reasonable? An
// alternative path would be a class called "LCInterp with all
// static functions. Then it could have an `operator()` which would
// be maybe nicer?
// TODO(JMM): Merge this w/ what Ben has done.
namespace Cent {
namespace Linear {
Did the current structure work out for you in Phoebus/Spiner or do you anticipate/suggest changes?
Shall this PR just provide the mininum with a potentially changing interface in the future or keep these abstraction given (expected) futures extensions?
There was a problem hiding this comment.
I have not needed one dimensional linear interpolation so I would remove this.
This minimal working code has worked for me in Phoebus, so I don't think there's a need to extend it or provide hooks. I would remove the dead code.
There was a problem hiding this comment.
Removed the dead code.
|
@pgrete what's the status of this? I think it would be nice to have this merged in. |
|
I should be able to make the finishing touches next week. Otherwise, if you got some spare time and would like to see this earlier, feel free to go ahead. |
|
Pulling the trigger as discussed. |
PR Summary
Add bi- and trilinear interpolation from Phoebus/Spiner.
PR Checklist