geomc 1.0
A c++ linear algebra template library
|
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 > | |
T | legendre (index_t l, index_t m, T x) |
template<typename T > | |
T | legendre (index_t n, T x) |
template<typename T > | |
T | legendre_integral (index_t n, T x) |
template<typename T > | |
T | chebyshev (index_t kind, index_t n, T x) |
template<typename T > | |
T | spherical_harmonic_coeff (index_t l, index_t m, T cos_alt, T azi) |
template<typename T > | |
T | diff_of_products (T a, T b, T c, T d) |
template<typename T > | |
T | sum_of_products (T a, T b, T c, T d) |
template<typename T > | |
T | multiply_add (T a, T b, T c) |
Fused multiply-add. More... | |
template<typename T > | |
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 > | |
T | interp_linear (const S s[], const T p[], int dim) |
template<typename S , typename T > | |
T | interp_cubic (S s, const T p[4]) |
template<typename S , typename T > | |
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) |
Tools for constructing and sampling continuous-valued objects.
|
strong |
Policy for handling discontinuous derivatives.
enum EdgeBehavior |
Raster edge-sampling behavior.
enum Interpolation |
T geom::chebyshev | ( | index_t | kind, |
index_t | n, | ||
T | x | ||
) |
Evaluate a Chebyshev polynomial of order n
at x
.
kind | 1 or 2 for the Chebyshev polynomial of the first or second kind, respectively. |
n | Order of the polynomial. |
x | Value of x at which to evaluate. |
|
inline |
Clamp v
between lo
and hi
.
v | Value to clamp. |
lo | Minimum value of v . |
hi | Maximum value of v . |
|
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.
T geom::interp_cubic | ( | const S | s[], |
const T | p[], | ||
int | dim | ||
) |
Cubic interpolation across dim
input dimensions.
s | An array of dim values; a point on (0, 1) dim representing the position in the center grid cell at which to interpolate. |
p | An array of the 4dim values adjacent to s . |
dim | Number of input dimensions. |
T geom::interp_cubic | ( | S | s, |
const T | p[4] | ||
) |
Cubic interpolation in 1 dimension.
s | Interpolation parameter between 0 and 1, representing the point between p[1] and p[2] at which to interpolate. |
p | An array of 4 values surrounding s . |
T geom::interp_linear | ( | const S | s[], |
const T | p[], | ||
int | dim | ||
) |
Linearly interpolate across dim
input dimensions.
s | An array of dim values; a point on (0, 1) dim at which to interpolate. |
p | An array of the 2dim values adjacent to s . |
dim | Number of input dimensions. |
T geom::legendre | ( | index_t | l, |
index_t | m, | ||
T | x | ||
) |
Evaluate an associated Legendre polynomial.
l | Band index. |
m | Sub-band. |
x | Value of x at which to evaluate the polynomial. |
T geom::legendre | ( | index_t | n, |
T | x | ||
) |
Evaluate a Legendre polynomal of order n
at x
. Equivalent to the associated legendre polynomial with m = 0.
n | Order of polynomial. |
x | Value of x at which to evaluate the polynomial. |
T geom::legendre_integral | ( | index_t | n, |
T | x | ||
) |
Evaluate the integral of the Legendre polynomal of order n
between -1 and x.
n | Order of polynomial. |
x | Value of x at which to evaluate the integral. |
|
inline |
Linearly interpolate between a
and b
.
t | Interpolation parameter, between 0 and 1. |
a | Value at t = 0 . |
b | Value at t = 1 . |
|
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).
bool geom::quadratic_solve | ( | T | results[2], |
T | a, | ||
T | b, | ||
T | c | ||
) |
Solve the quadratic equation ax
2+ 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.
[out] | results | Roots of the specified quadratic equation. |
a | Coefficient of x 2. | |
b | Coefficient of x . | |
c | Constant. |
true
if real roots exist and the coefficients are not all zero; false
otherwise (in which case the contents of results
will be unaltered). T geom::spherical_harmonic_coeff | ( | index_t | l, |
index_t | m, | ||
T | cos_alt, | ||
T | azi | ||
) |
Evaluate a spherical harmonic coefficient.
l | Band index. |
m | Sub-band. |
cos_alt | Cosine of angle from polar axis. |
azi | Azimuthal (longitude) angle. |
|
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.