corecel/grid/FindInterp.hh

0001 //----------------------------------*-C++-*----------------------------------//
0002 // Copyright 2021-2024 UT-Battelle, LLC, and other Celeritas developers.
0003 // See the top-level COPYRIGHT file for details.
0004 // SPDX-License-Identifier: (Apache-2.0 OR MIT)
0005 //---------------------------------------------------------------------------//
0006 //! \file corecel/grid/FindInterp.hh
0007 //---------------------------------------------------------------------------//
0008 #pragma once
0009
0010 #include "corecel/Assert.hh"
0011 #include "corecel/Macros.hh"
0012 #include "corecel/Types.hh"
0013
0014 namespace celeritas
0015 {
0016 //---------------------------------------------------------------------------//
0017 /*!
0018  * Result of finding a point on a grid for interpolating.
0019  *
0020  * The resulting index will be in [0, grid.size() - 1)
0021  * and the fraction will be in [0, 1).
0022  */
0023 template<class T>
0024 struct FindInterp
0025 {
0026     size_type index{};  //!< Lower index into the grid
0027     T fraction{};  //!< Fraction of the value between its neighbors
0028 };
0029
0030 //---------------------------------------------------------------------------//
0031 /*!
0032  * Find the index of the value and its fraction between neighboring points.
0033  *
0034  * The grid class should have a floating point value and must have methods \c
0035  * find, \c front, \c back, and \c operator[] .
0036  *
0037  * The value must be bounded by the grid and less than the final value. The
0038  * result will always have an index such that its neighbor to the right is a
0039  * valid point on the grid, and the fraction between neghbors may be zero (in
0040  * the case where the value is exactly on a grid point) but is always less than
0041  * one.
0042  */
0043 template<class Grid>
0044 inline CELER_FUNCTION FindInterp<typename Grid::value_type>
0045 find_interp(Grid const& grid, typename Grid::value_type value)
0046 {
0047     CELER_EXPECT(value >= grid.front() && value < grid.back());
0048
0049     FindInterp<typename Grid::value_type> result;
0050     result.index = grid.find(value);
0051     CELER_ASSERT(result.index + 1 < grid.size());
0052     auto const lower_val = grid[result.index];
0053     auto const upper_val = grid[result.index + 1];
0054     result.fraction = (value - lower_val) / (upper_val - lower_val);
0055
0056     return result;
0057 }
0058
0059 //---------------------------------------------------------------------------//
0060 }  // namespace celeritas