.. _interpolation: Gridding for Interpolation =========================== nSpiner performs interpolation on Cartesian-product grids. There are two lower-level objects: * ``RegularGrid1D`` * ``PiecewiseGrid1D`` * ``NonUniformGrid1D`` * ``FastNonUniformGrid1D`` These objects contain the metadata required for interpolation operations and have a few useful userspace functions, which are described here. Like ``DataBox``, these grid objects are templated on underlying data type, the default type being a ``Real`` as provided by ``ports-of-call``. You may wish to specialize to a specific type with a type alias such as: .. code-block:: cpp using RegularGrid1D = Spiner::RegularGrid1D; using PiecewiseGrid1D = Spiner::PiecewiseGrid1D; using NonUniformGrid1D = Spiner::NonUniformGrid1D; using FastNonUniformGrid1D = Spiner::FastNonUniformGrid1D; .. note:: In the function signature below we refer to ``T`` and ``Real`` as the underlying arithmetic data type. When constructing a ``DataBox``, you may wish to specify which interpolation object you are using. It is a template parameter. ``RegularGrid1D`` ------------------ We begin by discussing ``RegularGrid1D``, as the ``PiecewiseGrid1D`` object is built on top of it. Construction ^^^^^^^^^^^^^ A ``RegularGrid1D`` requires three values to specify an interpolation grid: the minimum value of the independent variable, the maximum value of the independent variable, and the number of points on the grid. These are passed into the constructor: .. cpp:function:: RegularGrid1D::RegularGrid1D(T min, T max, size_t N); Default constructors and copy constructors are also provided. Mapping an index to a real number and vice-versa ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ The function .. cpp:function:: T RegularGrid1D::x(const int i) const; returns a "physical" position on the grid given an index ``i``. The function .. cpp:function:: int RegularGrid1D::index(const T x) const; returns the index on the grid of a "physical" value ``x``. The function .. cpp:function:: T RegularGrid1D::min() const; returns the minimum value on the independent variable grid. The function .. cpp:function:: T RegularGrid1D::max() const; returns the maximum value on the independent variable grid. The function .. cpp:function:: int RegularGrid1D::nPoints() const; returns the number of points in the independent variable grid. ``NonUniformGrid1D`` ---------------------- ``NonUniformGrid1D`` stores one explicitly supplied coordinate for each grid point. It is appropriate when spacing cannot be represented efficiently as regular or piecewise-regular intervals. Each grid has independent coordinate storage, so a ``DataBox`` can use a different non-uniform coordinate sequence for every axis. Construction ^^^^^^^^^^^^^ Constructing from a ``std::vector`` or initializer list makes an owning host copy of the points: .. code-block:: cpp NonUniformGrid1D grid({-1.0, -0.5, 0.25, 2.0}); Construction always allocates host coordinates. Use ``getOnDevice()`` when a device copy is needed: .. code-block:: cpp NonUniformGrid1D host_grid({-1.0, -0.5, 0.25, 2.0}); NonUniformGrid1D device_grid = host_grid.getOnDevice(); The coordinates must contain at least two finite, strictly increasing values. The pointer constructor borrows caller-owned memory instead and validates those requirements in the supplied execution space (the host by default): .. code-block:: cpp std::vector points = {-1.0, -0.5, 0.25, 2.0}; NonUniformGrid1D view(points.data(), points.size()); For device-resident points, pass the execution space that can access them to the pointer constructor. The supplied execution space must match the memory space of the borrowed points; using an inaccessible pointer can cause a segmentation fault. Borrowed points must remain alive and unchanged for the grid's lifetime. Like ``DataBox``, ordinary grid copies are shallow reference-style copies; finalize an owned grid exactly once. Every grid type also provides an explicit ``copy(other)`` method. For ``NonUniformGrid1D``, it is a deep copy that: .. cpp:function:: void NonUniformGrid1D::copy(const NonUniformGrid1D& other) allocates independent host-owned coordinate storage and copies the coordinates from a host-resident ``other`` grid. It must not be used with a device-resident source; use ``getOnDevice()`` to make a deep device copy. For ``RegularGrid1D`` and ``PiecewiseGrid1D``, which do not own non-trivial dynamic memory, ``copy(other)`` is correspondingly a trivial metadata copy. Mapping and interpolation ^^^^^^^^^^^^^^^^^^^^^^^^^ ``x(i)`` returns the stored coordinate. ``index(x)`` finds its lower bracketing point with a portable binary search, and ``weights(x, ix, w)`` uses the spacing of that local interval. Lookup is ``O(log N)``. Values below or above the coordinate range use the first or last interval respectively, just as ``RegularGrid1D`` does, so interpolation extrapolates linearly. ``FastNonUniformGrid1D`` ------------------------ ``FastNonUniformGrid1D`` owns the same explicit physical coordinates as ``NonUniformGrid1D`` and can add an integer lookup table that makes interval lookup ``O(1)``. The table is uniform after transforming coordinates with the first-order Ports-of-Call NQT ``asinh`` function. This provides signed logarithmic spacing far from zero and linear spacing near zero. The ``Settings`` struct supplies construction settings. A positive ``scale`` sets the transition length in the physical coordinate's units. By default it is negative, which infers the scale as the smallest coordinate magnitude at least ``PortsOfCall::Robust::SMALL()``. The resolved positive scale is available through ``scale()`` and is stored in HDF5: .. code-block:: cpp FastNonUniformGrid1D grid(points); By default, the lookup table may contain at most 32 times as many entries as the physical grid has points. If an exact table would exceed that limit, the grid transparently uses the wrapped binary search. Pass a ``Settings`` object to select a scale, policy, or different limit: .. code-block:: cpp using Settings = FastNonUniformGrid1D::Settings; using Policy = FastNonUniformGrid1D::Policy; FastNonUniformGrid1D grid( points, Settings{ .scale = scale, .policy = Policy::RequireFast, .max_lookup_ratio = 16}); ``Automatic`` permits fallback, ``RequireFast`` fails construction if the table cannot fit, and ``ForceBinary`` skips the table. The selected mode can be inspected with ``usesFastLookup()`` and ``lookupSize()``. A host-owned grid can be reconfigured later, including after HDF5 loading: .. code-block:: cpp auto settings = grid.settings(); settings.scale = -1.0; // infer again from the physical coordinates settings.policy = FastNonUniformGrid1D::Policy::ForceBinary; grid.reconfigureLookup(settings); Reconfiguration follows the same explicit ownership convention as ``finalize()``: do not reconfigure an owner while shallow aliases depend on its lookup table. Reconfigure host storage before making device copies. Ordinary copies remain shallow. ``copy()``, ``getOnDevice()``, binary serialization, and HDF5 otherwise follow the ``NonUniformGrid1D`` lifecycle. HDF5 stores the physical coordinates and resolved construction settings, then rebuilds the derived lookup table when loading. Coordinate access is read-only because changing a coordinate would invalidate the table. The ``PiecewiseGrid1D`` ------------------------ A ``PiecewiseGrid1D`` is a non-intersecting, contiguous, ordered collection ``RegularGrid1D`` s. It can be used to construct grids with non-uniform spacing, so long as the grid spacing is piecewise constant. The maximum number of ``RegularGrid1D``s that can be used to construct a ``PiecewiseGrid1D`` is a compile-time parameter (default is 5). You can specify a different value with, e.g., .. code-block:: cpp // Maximum number of "pieces" in a grid = 10 using PiecewiseGrid1D = Spiner::PiecewiseGrid1D; Constructiong a ``PiecewiseGrid1D`` ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ A ``PiecewiseGrid1D`` is constructed from either a ``std::vector`` or a ``std::initializer_list`` of ``RegularGrid1D`` s. For example: .. code-block:: cpp // Initialize the regular grids // Note that the start and end points match // for each consecutive pair of grids. // g1 ends when g2 starts, etc. Spiner::RegularGrid1D g1(0, 0.25, 3); Spiner::RegularGrid1D g2(0.25, 0.75, 11); Spiner::RegularGrid1D g3(0.75, 1, 7); // Build the piecewise grid. The double bracket notation // is an "initalizer list" and is very convenient, // as it is a C++ language feature. Spiner::PiecewiseGrid1D h = {{g1, g2, g3}}; Default constructors and copy constructors are also provided. Index Mapping with ``PiecewiseGrid1D`` ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ A ``PiecewiseGrid1D`` has all the same functionality as ``RegularGrid1D``, but it automatically uses the relevant underlying grid spacing. The function .. cpp:function:: T PiecewiseGrid1D::x(const int i) const; returns a "physical" position on the grid given an index ``i``. The function .. cpp:function:: int PiecewiseGrid1D::index(const T x) const; returns the index on the grid of a "physical" value ``x``. The function .. cpp:function:: T PiecewiseGrid1D::min() const; returns the minimum value on the independent variable grid. The function .. cpp:function:: T PiecewiseGrid1D::max() const; returns the maximum value on the independent variable grid. The function .. cpp:function:: int PiecewiseGrid1D::nPoints() const; returns the number of points in the independent variable grid. Developer functionality ------------------------ All grid types implement the resource-management interface used by ``DataBox``: .. cpp:function:: std::size_t Grid::serializedSizeInBytes() const; .. cpp:function:: std::size_t Grid::dynamicMemorySizeInBytes() const; .. cpp:function:: std::size_t Grid::dumpDynamicMemory(char *dst) const; .. cpp:function:: std::size_t Grid::serialize(char *dst) const; .. cpp:function:: std::size_t Grid::deSerialize(char *src); .. cpp:function:: std::size_t Grid::setPointer(char *src); .. cpp:function:: void Grid::copy(const Grid& other); .. cpp:function:: Grid Grid::getOnDevice() const; .. cpp:function:: void Grid::finalize(); These operations are currently trivial for ``RegularGrid1D`` because it contains only inline data. ``PiecewiseGrid1D`` applies them recursively to its component grids. ``NonUniformGrid1D`` uses them to manage its coordinate array. This interface allows grid types to own dynamically allocated host or device data without requiring grid-specific resource handling in ``DataBox``. ``serialize`` writes the inline grid object followed by its dynamic memory. ``dumpDynamicMemory`` writes only the latter, allowing a grid embedded in a ``DataBox`` to avoid serializing its inline bytes twice. The binary serialization methods have the same transient, build-dependent compatibility limitations as ``DataBox`` serialization. See :ref:`the DataBox serialization documentation ` for details. Generative AI was used to assist with modifications to this page.