Template Class Mesh2DInterpolator

Class Documentation

template<size_t NF = 1>
class Mesh2DInterpolator

A class to perform interpolation on a 2D semi-tensor grid.

This class performs interpolation on a grid where the points in one direction (the y direction) are arranged in a tensor grid (i.e., in an Nx x Ny grid, there is one y value for each of the Nx rows), but not in the other direction (i.e., the values of x can be different from one row of y values to to the next). The number of quantities to be interpolated is determined by the template parameter NF.

Template Parameters:

NF – Number of quantities to inteprolate

Public Types

using Array3D = std::mdspan<double, std::dextents<size_t, 3>>
using Array2D = std::mdspan<double, std::dextents<size_t, 2>>
using Array1D = std::mdspan<double, std::dextents<size_t, 1>>

Public Functions

inline Mesh2DInterpolator(const Array2D &x, const Array1D &y, const Array3D &f, const gsl_interp_type *interpType = gsl_interp_steffen, const bool monotonic = false)

Construct a Mesh2DInterpolator.

If enforceMonotonic is set to true, the interpolation type will always be gsl_interp_linear regardless of the value of interpType.

Parameters:
  • x – A 2d array giving the x coordinates of the mesh points

  • y – A 1d array giving the y coordinates of the mesh points

  • f – A 3d array giving the value of the function to be interpolated

  • interpType – The type of interpolation to use

  • monotonic – Enforce that interpolation is monotonicity-preserving

inline Mesh2DInterpolator(const Array2D &x, const Array1D &y, const Array2D &f, const gsl_interp_type *interpType = gsl_interp_steffen, const bool monotonic = false)

Construct a Mesh2DInterpolator.

This constructor is valid only if NF == 1, and it just wraps a 2D array in a 3D view before invoking the alternative form of the constructor.

Parameters:
  • x – A 2d array giving the x coordinates of the mesh points

  • y – A 1d array giving the y coordinates of the mesh points

  • f – A 2d array giving the value of the function to be interpolated

  • interpType – The type of interpolation to use

  • monotonic – Enforce that interpolation is monotonicity-preserving

inline Mesh2DInterpolator(const Mesh2DGrid::XView &x, const Array1D &y, std::shared_ptr<const LazyValueSource<NF>> source, const gsl_interp_type *interpType = gsl_interp_steffen, const bool monotonic = false)

Construct a Mesh2DInterpolator whose mesh coordinates and function values are evaluated lazily.

Unlike the other constructors, this builds no rib or spine interpolators, and stores no function values: construction costs O(Ny) (see Mesh2DGrid’s own XView constructor), and each query instead evaluates source only at the few mesh points it needs, building a small local interpolator over just those points (see evalRib()/evalSpine()). Whatever storage x refers to, and source, must outlive this object. Query results agree with those of an eagerly-constructed Mesh2DInterpolator over the same values to within floating-point rounding. That requires a local interpolation scheme, whose value within a cell depends only on a few neighbouring points, so interpType must be gsl_interp_linear, gsl_interp_steffen, or gsl_interp_akima; any other type throws std::runtime_error.

Parameters:
  • x – The x coordinates of the mesh points (see Mesh2DGrid::XView); not copied

  • y – A 1d array giving the y coordinates of the mesh points

  • source – Supplies the function value at each mesh point on demand

  • interpType – The type of interpolation to use

  • monotonic – Enforce that interpolation is monotonicity-preserving

virtual ~Mesh2DInterpolator() = default
Mesh2DInterpolator(const Mesh2DInterpolator&) = delete
Mesh2DInterpolator(Mesh2DInterpolator&&) = default
auto operator=(const Mesh2DInterpolator&) -> Mesh2DInterpolator& = delete
auto operator=(Mesh2DInterpolator&&) noexcept -> Mesh2DInterpolator& = default
inline auto size() const

Return size of Mesh2DInterpolator.

Returns:

Size of mesh

inline auto nx() const

Return x dimension of Mesh2DInterpolator.

Returns:

Number of mesh points in x direction

inline auto ny() const

Return x dimension of Mesh2DInterpolator.

Returns:

Number of mesh points in x direction

inline auto convex() const

Return if the mesh is convex.

Returns:

True if the mesh is convex

inline auto x() const -> const auto&

Return a const reference to the x array.

Returns:

A const reference to the x array

inline auto y() const -> const auto&

Return a const reference to the y array.

Returns:

A const reference to the y array

inline auto s() const -> const auto&

Return a const reference to the s array.

Returns:

A const reference to the s array

inline auto xMin() const

Minimum x value in mesh.

Returns:

Minimum x value in mesh

inline auto xMin(const double y) const

Find minimum x value at given y.

Returns:

Minimum x value at given y

inline auto xMax() const

Maximum x value in mesh.

Returns:

Maximum x value in mesh

inline auto xMax(const double y) const

Find maximum x value at given y.

Returns:

Maximum x value at given y

inline auto yMin() const

Minimum y value in mesh.

Returns:

Minimum y value in mesh

inline auto yMax() const

Maximum y value in mesh.

Returns:

Maximum y value in mesh

inline auto yLim(double x) const

Find minimum and maximum y position at a given x in the mesh.

The return value for this function is a vectors of pairs; each pair contains the minimum and maximum y values where a line segment at a given x value intersects the mesh. For a convex mesh there is at most one minimum and one maximum, so this vector will be of length 0 (if the mesh is missed entirely) or 1, but for a non-convex mesh a line at fixed x can intersect the edge an arbitrary number of times, so the vector can be longer.

Parameters:

x – x position

Returns:

Minimum and maximum y values

inline auto yEdgeSlope(double x, bool leftEdge) const

Find y and dy/dx at every point where a line of constant x intersects one mesh edge.

Thin wrapper around Mesh2DGrid::yEdgeSlope() – see its own comment.

Parameters:
  • x – x position

  • leftEdge – True to find intersections with the mesh’s left edge (column 0), false for its right edge (column nx()-1)

Returns:

One (y, dy/dx) pair per intersection with the requested edge, in ascending y order

inline auto contains(const double x, const double y) const

Check whether a point is contained in the mesh.

Parameters:
  • x – x position

  • y – y position

Returns:

True if the point is within the mesh

inline auto interpConstX(double x) const -> std::vector<std::unique_ptr<Interpolator1D<NF>>>

Function to create a 1D interpolator at constant x.

This function returns a vector of Interpolator1D objeects rather than a single one because, for a non-convex mesh, the intersection of a line of constant x with the mesh may consist of multiple disconnected sgments. In such cases, the vector will contain one Interpolator1D for each segment. If the value of x provided is outside the mesh entirely, so there are no intersections, the returned vector will be empty. A segment consisting of a single point (e.g. a line of constant x that is tangent to the mesh, touching it at only one mass) carries no interpolable information and is likewise dropped, exactly as if that portion of the mesh had been missed.

Parameters:

x – x coordinate

Returns:

A vector of Interpolator1D’s that interpolates in y at the given x

inline auto interpConstY(double y) const -> std::unique_ptr<Interpolator1D<NF>>

Function to create a 1D interpolator at constant y.

If the value of y provided is outside the mesh entirely (i.e., below yMin or above yMax), this routine returns a null pointer.

Parameters:

y – y coordinate

Returns:

An Interpolator1D that interpolates in x and the given y

inline auto operator()(double x, double y, bool linear = false) const

Interpolate to a single point.

Return type is a double if NF = 1, and a std::array otherwise. This routine builds a local Interpolator1D from just enough points around (x, y) to support the requested interpolation type, then evaluates it at y; unlike interpConstX, it does not need to handle multiple disconnected segments, since xIntersectN only ever searches out from a single starting point. If xIntersectN cannot find enough points to support the requested type but still finds at least 2, this falls back to linear interpolation instead, mirroring interpConstX’s handling of undersized segments. If it finds only a single point, (x, y) must be an isolated tangent point of the mesh, so there is nothing to interpolate and that point’s own value is returned directly.

Parameters:
  • x – x coordinate

  • y – y coordinate

  • linear – If true, force linear interpolation instead of this object’s default interpolation type

Returns:

The interpolated value(s)