Gridding for Interpolation

nSpiner performs interpolation on Cartesian-product grids. There are two lower-level objects:

• RegularGrid1D

• PiecewiseGrid1D

• NonUniformGrid1D

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:

using RegularGrid1D = Spiner::RegularGrid1D<double>;
using PiecewiseGrid1D = Spiner::PiecewiseGrid1D<double>;
using NonUniformGrid1D = Spiner::NonUniformGrid1D<double>;

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:

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

T RegularGrid1D::x(const int i) const;

returns a “physical” position on the grid given an index i.

The function

int RegularGrid1D::index(const T x) const;

returns the index on the grid of a “physical” value x.

The function

T RegularGrid1D::min() const;

returns the minimum value on the independent variable grid.

The function

T RegularGrid1D::max() const;

returns the maximum value on the independent variable grid.

The 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<double, NonUniformGrid1D> 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:

NonUniformGrid1D grid({-1.0, -0.5, 0.25, 2.0});

Construction always allocates host coordinates. Use getOnDevice() when a device copy is needed:

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):

std::vector<double> 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:

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.

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.,

// Maximum number of "pieces" in a grid = 10
using PiecewiseGrid1D = Spiner::PiecewiseGrid1D<double, 10>;

Constructiong a PiecewiseGrid1D

A PiecewiseGrid1D is constructed from either a std::vector or a std::initializer_list of RegularGrid1D s. For example:

// 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<double> g1(0, 0.25, 3);
Spiner::RegularGrid1D<double> g2(0.25, 0.75, 11);
Spiner::RegularGrid1D<double> 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<double> 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

T PiecewiseGrid1D::x(const int i) const;

returns a “physical” position on the grid given an index i.

The function

int PiecewiseGrid1D::index(const T x) const;

returns the index on the grid of a “physical” value x.

The function

T PiecewiseGrid1D::min() const;

returns the minimum value on the independent variable grid.

The function

T PiecewiseGrid1D::max() const;

returns the maximum value on the independent variable grid.

The 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:

std::size_t Grid::serializedSizeInBytes() const;

std::size_t Grid::dynamicMemorySizeInBytes() const;

std::size_t Grid::dumpDynamicMemory(char *dst) const;

std::size_t Grid::serialize(char *dst) const;

std::size_t Grid::deSerialize(char *src);

std::size_t Grid::setPointer(char *src);

void Grid::copy(const Grid &other);

Grid Grid::getOnDevice() const;

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 the DataBox serialization documentation for details.

Generative AI was used to assist with modifications to this page.