geomc 1.0
A c++ linear algebra template library
Classes | Enumerations | Functions
Function

Tools for constructing and sampling continuous-valued objects. More...

Classes

class  Dual< T, DP >
 Class implementing the dual numbers, whose arithmetic operations perform a simultaneous calculation of the first derivative. More...
 
class  PerlinNoise< T, N >
 Real-valued smooth noise over N dimensions. More...
 
class  Raster< I, O, M, N >
 An M-dimensional grid of interpolated data which can be continuously sampled. More...
 
class  SphericalHarmonics< T, Bands >
 Class and methods for representing band-limited functions on the surface of a 3D sphere. More...
 
class  ZonalHarmonics< T, Bands >
 Represents a function on the N-sphere with axial symmetry. More...
 

Enumerations

enum class  DiscontinuityPolicy {
  NaN , Average , Left , Right ,
  Inf
}
 Policy for handling discontinuous derivatives. More...
 
enum  EdgeBehavior { EDGE_CLAMP , EDGE_PERIODIC , EDGE_MIRROR , EDGE_CONSTANT }
 Raster edge-sampling behavior. More...
 
enum  Interpolation { INTERP_NEAREST , INTERP_LINEAR , INTERP_CUBIC }
 Behavior for sampling between data points. More...
 

Functions

template<typename T , index_t N>
void legendre (SphericalHarmonics< T, N > *sh, index_t m, T x)
 
template<typename T >
legendre (index_t l, index_t m, T x)
 
template<typename T >
legendre (index_t n, T x)
 
template<typename T >
legendre_integral (index_t n, T x)
 
template<typename T >
chebyshev (index_t kind, index_t n, T x)
 
template<typename T >
spherical_harmonic_coeff (index_t l, index_t m, T cos_alt, T azi)
 
template<typename T >
diff_of_products (T a, T b, T c, T d)
 
template<typename T >
sum_of_products (T a, T b, T c, T d)
 
template<typename T >
multiply_add (T a, T b, T c)
 Fused multiply-add. More...
 
template<typename T >
clamp (const T &v, const T &lo, const T &hi)
 
template<typename S , typename T >
std::common_type< S, T >::type mix (S t, T a, T b)
 
template<typename S , typename T >
interp_linear (const S s[], const T p[], int dim)
 
template<typename S , typename T >
interp_cubic (S s, const T p[4])
 
template<typename S , typename T >
interp_cubic (const S s[], const T p[], int dim)
 
template<typename T >
bool quadratic_solve (T results[2], T a, T b, T c)
 

Detailed Description

Tools for constructing and sampling continuous-valued objects.

Enumeration Type Documentation

◆ DiscontinuityPolicy

enum class DiscontinuityPolicy
strong

Policy for handling discontinuous derivatives.

Enumerator
NaN 

Return NaN at discontinuities.

Average 

At discontinuities, return the average of the two boundary values.

Left 

At discontinuities, return the value when approaching from the left.

Right 

At discontinuities, return the value when approaching from the right.

Inf 

At discontinuities, return +/- infinity, according to whether the discontinuity is an increasing or decreasing jump.

◆ EdgeBehavior

Raster edge-sampling behavior.

Enumerator
EDGE_CLAMP 

Clip sample coordinates to the sampleable area, thus repeating edge values beyond the boundary.

EDGE_PERIODIC 

Wrap around sample coordinates to the opposite edge, thus tiling the sampled data beyond the boundary.

EDGE_MIRROR 

Mirror the sample coordinates across edges.

EDGE_CONSTANT 

Regions outside the sampleable area have a uniform, defined value (zero by default).

◆ Interpolation

Behavior for sampling between data points.

Enumerator
INTERP_NEAREST 

Return the data nearest to the sample point.

INTERP_LINEAR 

Linearly interpolate the nearest 2n data points.

INTERP_CUBIC 

Cubically interpolate the nearest 4n data points.

Function Documentation

◆ chebyshev()

T geom::chebyshev ( index_t  kind,
index_t  n,
x 
)

Evaluate a Chebyshev polynomial of order n at x.

#include <geomc/function/Basis.h>
Parameters
kind1 or 2 for the Chebyshev polynomial of the first or second kind, respectively.
nOrder of the polynomial.
xValue of x at which to evaluate.

◆ clamp()

T geom::clamp ( const T &  v,
const T &  lo,
const T &  hi 
)
inline

Clamp v between lo and hi.

Parameters
vValue to clamp.
loMinimum value of v.
hiMaximum value of v.

◆ diff_of_products()

T diff_of_products ( a,
b,
c,
d 
)
inline

A high-precision method for computing (a * b) - (c * d). In cases where the two products are large and close in value, the result can be wildly wrong due to catastrophic cancellation when using the "naive" method. This function will return a result within 1.5 ULP of the exact answer.

A high-precision method is available for float, double, long double, and Dual numbers composed of those types. Other types will fall back to the "naive" method.

◆ interp_cubic() [1/2]

T geom::interp_cubic ( const S  s[],
const T  p[],
int  dim 
)

Cubic interpolation across dim input dimensions.

Parameters
sAn array of dim values; a point on (0, 1)dim representing the position in the center grid cell at which to interpolate.
pAn array of the 4dim values adjacent to s.
dimNumber of input dimensions.
Returns
Interpolated value.

◆ interp_cubic() [2/2]

T geom::interp_cubic ( s,
const T  p[4] 
)

Cubic interpolation in 1 dimension.

Parameters
sInterpolation parameter between 0 and 1, representing the point between p[1] and p[2] at which to interpolate.
pAn array of 4 values surrounding s.
Returns
Interpolated value.

◆ interp_linear()

T geom::interp_linear ( const S  s[],
const T  p[],
int  dim 
)

Linearly interpolate across dim input dimensions.

Parameters
sAn array of dim values; a point on (0, 1)dim at which to interpolate.
pAn array of the 2dim values adjacent to s.
dimNumber of input dimensions.
Returns
Interpolated value.

◆ legendre() [1/2]

T geom::legendre ( index_t  l,
index_t  m,
x 
)

Evaluate an associated Legendre polynomial.

#include <geomc/function/Basis.h>
Parameters
lBand index.
mSub-band.
xValue of x at which to evaluate the polynomial.

◆ legendre() [2/2]

T geom::legendre ( index_t  n,
x 
)

Evaluate a Legendre polynomal of order n at x. Equivalent to the associated legendre polynomial with m = 0.

#include <geomc/function/Basis.h>
Parameters
nOrder of polynomial.
xValue of x at which to evaluate the polynomial.

◆ legendre_integral()

T geom::legendre_integral ( index_t  n,
x 
)

Evaluate the integral of the Legendre polynomal of order n between -1 and x.

#include <geomc/function/Basis.h>
Parameters
nOrder of polynomial.
xValue of x at which to evaluate the integral.

◆ mix()

std::common_type< S, T >::type geom::mix ( t,
a,
b 
)
inline

Linearly interpolate between a and b.

Parameters
tInterpolation parameter, between 0 and 1.
aValue at t = 0.
bValue at t = 1.
Returns
Interpolated value.

◆ multiply_add()

T multiply_add ( a,
b,
c 
)
inline

Fused multiply-add.

Compute a * b + c using a fused multiply add operation (if available), which may be both more performant and more precise than multiple instructions.

Overloaded such that the calculation falls back to an ordinary composite multiply-add if a fused instruction is not available for T. (Compare to std::fma() which will promote to floating point, possibly losing precision).

◆ quadratic_solve()

bool geom::quadratic_solve ( results[2],
a,
b,
c 
)

Solve the quadratic equation ax2+ bx + c = 0. Formulated to preserve precision and avoid any divisions by zero.

In the case where there is exactly one real solution, the two roots will be set equal to each other.

Parameters
[out]resultsRoots of the specified quadratic equation.
aCoefficient of x2.
bCoefficient of x.
cConstant.
Returns
true if real roots exist and the coefficients are not all zero; false otherwise (in which case the contents of results will be unaltered).

◆ spherical_harmonic_coeff()

T geom::spherical_harmonic_coeff ( index_t  l,
index_t  m,
cos_alt,
azi 
)

Evaluate a spherical harmonic coefficient.

#include <geomc/function/SphericalHarmonics.h>
Parameters
lBand index.
mSub-band.
cos_altCosine of angle from polar axis.
aziAzimuthal (longitude) angle.

◆ sum_of_products()

T sum_of_products ( a,
b,
c,
d 
)
inline

A high-precision method for computing (a * b) + (c * d). In cases where the two products are large and similar in magnitude but opposite in sign, the result can be wildly wrong due to catastrophic cancellation when using the "naive" method. This function will return a result within 1.5 ULP of the exact answer.

A high-precision method is available for float, double, long double, and Dual numbers composed of those types. Other types will fall back to the "naive" method.