nerdcruft/libcruft-util
1
0
Fork 0
Code Issues Releases Activity

maths: move remaining operations into util namespace

Browse Source
This commit is contained in:
Danny Robson 2015-11-16 11:42:20 +11:00
parent 5a26793b0a
commit b1bc54ac8c
15 changed files with 571 additions and 567 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
debug.hpp
View File

@ -82,7 +82,7 @@
DEBUG_ONLY( \
const auto __a = (A); \
const auto __b = (B); \
_CHECK_META (almost_equal (__a, __b), \
_CHECK_META (util::almost_equal (__a, __b), \
{ ; }, \
{ \
std::ostringstream __debug_os; \
@ -203,7 +203,7 @@
DEBUG_ONLY( \
const auto __a = (A); \
const auto __b = (B); \
_CHECK_META (!almost_equal (__a, __b), \
_CHECK_META (!util::almost_equal (__a, __b), \
{ ; }, \
{ \
std::ostringstream __debug_neq_os; \

6
json/schema.cpp
View File

@ -175,7 +175,7 @@ validate (json::tree::string &node,
auto maxLength = schema.find ("maxLength");
if (maxLength != schema.cend ()) {
auto cmp = maxLength->second->as_number ().native ();
if (!is_integer (cmp))
if (!util::is_integer (cmp))
throw length_error ("maxLength");
if (val.size () > cmp)
@ -186,7 +186,7 @@ validate (json::tree::string &node,
auto minLength = schema.find ("minLength");
if (minLength != schema.cend ()) {
auto cmp = minLength->second->as_number ().native ();
if (!is_integer (cmp))
if (!util::is_integer (cmp))
throw length_error ("minLength");
if (val.size () < cmp)
@ -217,7 +217,7 @@ validate (json::tree::number &node,
if (mult != schema.cend ()) {
auto div = mult->second->as_number ().native ();
if (val <= 0 || almost_equal (val, div))
if (val <= 0 || util::almost_equal (val, div))
throw json::schema_error ("multipleOf");
}

37
json/tree.cpp
View File

@ -48,24 +48,25 @@ using json::tree::null;
///////////////////////////////////////////////////////////////////////////////
template <>
bool
is_integer (const json::tree::number &node)
{
return is_integer (node.native ());
namespace util {
template <>
bool
is_integer (const json::tree::number &node)
{
return is_integer (node.native ());
}
//-----------------------------------------------------------------------------
template <>
bool
is_integer (const json::tree::node &node)
{
return node.is_number () &&
is_integer (node.as_number ());
}
}
//-----------------------------------------------------------------------------
template <>
bool
is_integer (const json::tree::node &node)
{
return node.is_number () &&
is_integer (node.as_number ());
}
///////////////////////////////////////////////////////////////////////////////
static std::vector<json::flat::item>::const_iterator
parse (std::vector<json::flat::item>::const_iterator first,
@ -305,7 +306,7 @@ size_t
json::tree::node::as_uint (void) const
{
auto val = as_number ().native ();
if (!is_integer (val))
if (!util::is_integer (val))
throw json::type_error ("cast fractional value to uint");
// TODO: use trunc_cast
@ -761,7 +762,7 @@ json::tree::number::write (std::ostream &os) const
//-----------------------------------------------------------------------------
bool
json::tree::number::operator ==(const json::tree::number &rhs) const
{ return almost_equal (rhs.m_value, m_value); }
{ return util::almost_equal (rhs.m_value, m_value); }
///////////////////////////////////////////////////////////////////////////////

91
maths.cpp
View File

@ -27,34 +27,36 @@
///////////////////////////////////////////////////////////////////////////////
template <typename T>
bool
is_pow2 (T value) {
util::is_pow2 (T value)
{
typedef typename std::enable_if<std::is_integral<T>::value, bool>::type return_type;
return (return_type)(value && !(value & (value - 1)));
}
template bool is_pow2 (uint8_t);
template bool is_pow2 (uint16_t);
template bool is_pow2 (uint32_t);
template bool is_pow2 (uint64_t);
template bool util::is_pow2 (uint8_t);
template bool util::is_pow2 (uint16_t);
template bool util::is_pow2 (uint32_t);
template bool util::is_pow2 (uint64_t);
///////////////////////////////////////////////////////////////////////////////
template <typename T>
T
log2up (T v)
util::log2up (T v)
{
return log2 ((v << 1) - 1);
}
template uint32_t log2up (uint32_t);
template uint64_t log2up (uint64_t);
template uint32_t util::log2up (uint32_t);
template uint64_t util::log2up (uint64_t);
///////////////////////////////////////////////////////////////////////////////
template <typename T>
T
log2 (T v) {
util::log2 (T v)
{
static_assert (std::is_integral<T>::value,
"log2 is only implemented for integers");
@ -65,28 +67,29 @@ log2 (T v) {
return l;
}
template uint8_t log2 (uint8_t);
template uint16_t log2 (uint16_t);
template uint32_t log2 (uint32_t);
template uint64_t log2 (uint64_t);
template uint8_t util::log2 (uint8_t);
template uint16_t util::log2 (uint16_t);
template uint32_t util::log2 (uint32_t);
template uint64_t util::log2 (uint64_t);
///////////////////////////////////////////////////////////////////////////////
template <typename T>
double
rootsquare (T a, T b)
util::rootsquare (T a, T b)
{ return sqrt (util::pow2 (a) + util::pow2 (b)); }
//-----------------------------------------------------------------------------
template double rootsquare (double, double);
template double rootsquare ( int, int);
template double util::rootsquare (double, double);
template double util::rootsquare ( int, int);
///////////////////////////////////////////////////////////////////////////////
template <typename T>
bool
is_integer (const T &value) {
util::is_integer (const T &value)
{
T integer;
return exactly_equal (std::modf (value, &integer),
static_cast<T> (0.0));
@ -94,24 +97,27 @@ is_integer (const T &value) {
//-----------------------------------------------------------------------------
template bool is_integer (const double&);
template bool is_integer (const float&);
template bool util::is_integer (const double&);
template bool util::is_integer (const float&);
///////////////////////////////////////////////////////////////////////////////
template <>
unsigned
digits (const uint32_t &v) {
return (v >= 1000000000) ? 10 :
(v >= 100000000) ? 9 :
(v >= 10000000) ? 8 :
(v >= 1000000) ? 7 :
(v >= 100000) ? 6 :
(v >= 10000) ? 5 :
(v >= 1000) ? 4 :
(v >= 100) ? 3 :
(v >= 10) ? 2 :
1;
namespace util {
template <>
unsigned
digits (const uint32_t &v)
{
return (v >= 1000000000) ? 10 :
(v >= 100000000) ? 9 :
(v >= 10000000) ? 8 :
(v >= 1000000) ? 7 :
(v >= 100000) ? 6 :
(v >= 10000) ? 5 :
(v >= 1000) ? 4 :
(v >= 100) ? 3 :
(v >= 10) ? 2 :
1;
}
}
@ -120,7 +126,8 @@ template <typename T>
std::enable_if_t<
std::is_integral<T>::value, T
>
round_pow2 (T value) {
util::round_pow2 (T value)
{
using return_type = std::enable_if_t<std::is_integral<T>::value, T>;
--value;
@ -135,15 +142,15 @@ round_pow2 (T value) {
//-----------------------------------------------------------------------------
template uint8_t round_pow2 (uint8_t);
template uint16_t round_pow2 (uint16_t);
template uint32_t round_pow2 (uint32_t);
template uint64_t round_pow2 (uint64_t);
template uint8_t util::round_pow2 (uint8_t);
template uint16_t util::round_pow2 (uint16_t);
template uint32_t util::round_pow2 (uint32_t);
template uint64_t util::round_pow2 (uint64_t);
///////////////////////////////////////////////////////////////////////////////
template const float PI<float>;
template const double PI<double>;
template const float util::PI<float>;
template const double util::PI<double>;
///////////////////////////////////////////////////////////////////////////////
@ -151,7 +158,7 @@ template const double PI<double>;
// so it's easier to instantiate early and check for broken code at library
// build time.
template float limit (float, float, float);
template float util::limit (float, float, float);
template float smoothstep (float, float, float);
template double smoothstep (double, double, double);
template float util::smoothstep (float, float, float);
template double util::smoothstep (double, double, double);

842
maths.hpp
View File

@ -44,326 +44,318 @@ namespace util {
{
return t > 0 ? t : -t;
}
}
///////////////////////////////////////////////////////////////////////////
// exponentials
///////////////////////////////////////////////////////////////////////////////
// exponentials
namespace util {
template <typename T>
constexpr T
pow2 [[gnu::const]] (T value)
{
return value * value;
}
}
///////////////////////////////////////////////////////////////////////////////
namespace util {
///////////////////////////////////////////////////////////////////////////
template <typename T>
constexpr T
pow [[gnu::const]] (T x, unsigned y);
}
//-----------------------------------------------------------------------------
template <typename T>
bool
is_pow2 (T value);
//-------------------------------------------------------------------------
template <typename T>
bool
is_pow2 (T value);
//-----------------------------------------------------------------------------
// Logarithms
template <typename T>
T
log2 (T val);
template <typename T>
T
log2up (T val);
///////////////////////////////////////////////////////////////////////////////
// Roots
template <typename T>
double
rootsquare (T a, T b);
///////////////////////////////////////////////////////////////////////////////
// Rounding
template <typename T, typename U>
inline
typename std::common_type<
std::enable_if_t<std::is_integral<T>::value,T>,
std::enable_if_t<std::is_integral<U>::value,U>
>::type
round_to (T value, U size)
{
if (value % size == 0)
return value;
return value + (size - value % size);
}
//-----------------------------------------------------------------------------
template <typename T>
std::enable_if_t<
std::is_integral<T>::value, T
>
round_pow2 (T value);
//-----------------------------------------------------------------------------
template <typename T, typename U>
constexpr std::enable_if_t<
std::is_integral<T>::value &&
std::is_integral<U>::value,
//-----------------------------------------------------------------------------
// Logarithms
template <typename T>
T
>
divup (const T a, const U b)
{
return (a + b - 1) / b;
}
log2 (T val);
///////////////////////////////////////////////////////////////////////////////
// Properties
template <typename T>
bool
is_integer (const T& value);
template <typename T>
T
log2up (T val);
//-----------------------------------------------------------------------------
template <typename T>
unsigned
digits (const T& value);
///////////////////////////////////////////////////////////////////////////////
// Roots
template <typename T>
double
rootsquare (T a, T b);
//-----------------------------------------------------------------------------
constexpr int sign (int);
constexpr float sign (float);
constexpr double sign (double);
///////////////////////////////////////////////////////////////////////////////
// Rounding
template <typename T, typename U>
inline
typename std::common_type<
std::enable_if_t<std::is_integral<T>::value,T>,
std::enable_if_t<std::is_integral<U>::value,U>
>::type
round_to (T value, U size)
{
if (value % size == 0)
return value;
return value + (size - value % size);
}
///////////////////////////////////////////////////////////////////////////////
// factorisation
template <typename T>
constexpr T
gcd (T a, T b)
{
if (a == b) return a;
if (a > b) return gcd (a - b, b);
if (b > a) return gcd (a, b - a);
unreachable ();
}
//-----------------------------------------------------------------------------
template <typename T>
std::enable_if_t<
std::is_integral<T>::value, T
>
round_pow2 (T value);
///////////////////////////////////////////////////////////////////////////////
// Comparisons
inline bool
almost_equal (const float &a, const float &b)
{
return ieee_single::almost_equal (a, b);
}
//-----------------------------------------------------------------------------
template <typename T, typename U>
constexpr std::enable_if_t<
std::is_integral<T>::value &&
std::is_integral<U>::value,
T
>
divup (const T a, const U b)
{
return (a + b - 1) / b;
}
//-----------------------------------------------------------------------------
inline bool
almost_equal (const double &a, const double &b)
{
return ieee_double::almost_equal (a, b);
}
//-----------------------------------------------------------------------------
template <typename A, typename B>
typename std::enable_if_t<
std::is_floating_point<A>::value &&
std::is_floating_point<B>::value,
///////////////////////////////////////////////////////////////////////////////
// Properties
template <typename T>
bool
>
almost_equal (const A &a, const B &b)
{
using common_t = std::common_type_t<A,B>;
return almost_equal<common_t> (static_cast<common_t> (a),
static_cast<common_t> (b));
}
is_integer (const T& value);
//-----------------------------------------------------------------------------
template <typename A, typename B>
typename std::enable_if_t<
std::is_integral<A>::value &&
std::is_integral<B>::value &&
std::is_signed<A>::value == std::is_signed<B>::value,
bool
>
almost_equal (const A &a, const B &b) {
using common_t = std::common_type_t<A,B>;
return static_cast<common_t> (a) == static_cast<common_t> (b);
}
//-----------------------------------------------------------------------------
template <typename T>
unsigned
digits (const T& value);
//-----------------------------------------------------------------------------
constexpr int sign (int);
constexpr float sign (float);
constexpr double sign (double);
//-----------------------------------------------------------------------------
template <typename Ta, typename Tb>
typename std::enable_if<
!std::is_arithmetic<Ta>::value ||
!std::is_arithmetic<Tb>::value,
bool
>::type
almost_equal (const Ta &a, const Tb &b)
{ return a == b; }
///////////////////////////////////////////////////////////////////////////////
// factorisation
template <typename T>
constexpr T
gcd (T a, T b)
{
if (a == b) return a;
if (a > b) return gcd (a - b, b);
if (b > a) return gcd (a, b - a);
unreachable ();
}
//-----------------------------------------------------------------------------
// Useful for explictly ignore equality warnings
///////////////////////////////////////////////////////////////////////////////
// Comparisons
inline bool
almost_equal (const float &a, const float &b)
{
return ieee_single::almost_equal (a, b);
}
//-----------------------------------------------------------------------------
inline bool
almost_equal (const double &a, const double &b)
{
return ieee_double::almost_equal (a, b);
}
//-----------------------------------------------------------------------------
template <typename A, typename B>
typename std::enable_if_t<
std::is_floating_point<A>::value &&
std::is_floating_point<B>::value,
bool
>
almost_equal (const A &a, const B &b)
{
using common_t = std::common_type_t<A,B>;
return almost_equal<common_t> (static_cast<common_t> (a),
static_cast<common_t> (b));
}
//-----------------------------------------------------------------------------
template <typename A, typename B>
typename std::enable_if_t<
std::is_integral<A>::value &&
std::is_integral<B>::value &&
std::is_signed<A>::value == std::is_signed<B>::value,
bool
>
almost_equal (const A &a, const B &b) {
using common_t = std::common_type_t<A,B>;
return static_cast<common_t> (a) == static_cast<common_t> (b);
}
//-----------------------------------------------------------------------------
template <typename Ta, typename Tb>
typename std::enable_if<
!std::is_arithmetic<Ta>::value ||
!std::is_arithmetic<Tb>::value,
bool
>::type
almost_equal (const Ta &a, const Tb &b)
{ return a == b; }
//-----------------------------------------------------------------------------
// Useful for explictly ignore equality warnings
#pragma GCC diagnostic push
#pragma GCC diagnostic ignored "-Wfloat-equal"
template <typename T, typename U>
bool
exactly_equal (const T &a, const U &b)
{ return a == b; }
template <typename T, typename U>
bool
exactly_equal (const T &a, const U &b)
{ return a == b; }
#pragma GCC diagnostic pop
//-----------------------------------------------------------------------------
template <typename T>
bool
almost_zero (T a)
{ return almost_equal (a, T{0}); }
//-----------------------------------------------------------------------------
template <typename T>
bool
almost_zero (T a)
{ return almost_equal (a, T{0}); }
//-----------------------------------------------------------------------------
template <typename T>
bool
exactly_zero (T a)
{ return exactly_equal (a, T{0}); }
//-----------------------------------------------------------------------------
template <typename T>
bool
exactly_zero (T a)
{ return exactly_equal (a, T{0}); }
//-----------------------------------------------------------------------------
template <typename T>
const T&
identity (const T& t)
{
return t;
}
///////////////////////////////////////////////////////////////////////////////
// angles, trig
template <typename T>
constexpr T PI = T(3.141592653589793238462643);
//-----------------------------------------------------------------------------
template <typename T>
constexpr T E = T(2.71828182845904523536028747135266250);
//-----------------------------------------------------------------------------
template <typename T>
constexpr T
to_degrees (T radians)
{
static_assert (std::is_floating_point<T>::value, "undefined for integral types");
return radians * 180 / PI<T>;
}
//-----------------------------------------------------------------------------
template <typename T>
constexpr T
to_radians (T degrees)
{
static_assert (std::is_floating_point<T>::value, "undefined for integral types");
return degrees / 180 * PI<T>;
}
//-----------------------------------------------------------------------------
//! Normalised sinc function
template <typename T>
constexpr T
sincn (T x)
{
return almost_zero (x) ? 1 : std::sin (PI<T> * x) / (PI<T> * x);
}
//-----------------------------------------------------------------------------
//! Unnormalised sinc function
template <typename T>
constexpr T
sincu (T x)
{
return almost_zero (x) ? 1 : std::sin (x) / x;
}
///////////////////////////////////////////////////////////////////////////////
// combinatorics
constexpr uintmax_t
factorial (unsigned i)
{
return i <= 1 ? 0 : i * factorial (i - 1);
}
//-----------------------------------------------------------------------------
/// stirlings approximation of factorials
constexpr uintmax_t
stirling (unsigned n)
{
return static_cast<uintmax_t> (
std::sqrt (2 * PI<float> * n) * std::pow (n / E<float>, n)
);
}
//-----------------------------------------------------------------------------
constexpr uintmax_t
combination (unsigned n, unsigned k)
{
return factorial (n) / (factorial (k) / (factorial (n - k)));
}
///////////////////////////////////////////////////////////////////////////////
// kahan summation for long floating point sequences
template <class InputIt>
typename std::iterator_traits<InputIt>::value_type
fsum (InputIt first, InputIt last)
{
using T = typename std::iterator_traits<InputIt>::value_type;
static_assert (std::is_floating_point<T>::value,
"fsum only works for floating point types");
T sum = 0;
T c = 0;
for (auto cursor = first; cursor != last; ++cursor) {
T y = *cursor - c;
T t = sum + y;
c = (t - sum) - y;
sum = t;
//-----------------------------------------------------------------------------
template <typename T>
const T&
identity (const T& t)
{
return t;
}
return sum;
}
///////////////////////////////////////////////////////////////////////////////
// angles, trig
template <typename T>
constexpr T PI = T(3.141592653589793238462643);
//-----------------------------------------------------------------------------
template <typename T>
constexpr T E = T(2.71828182845904523536028747135266250);
///////////////////////////////////////////////////////////////////////////////
/// Variadic minimum
namespace util {
//-----------------------------------------------------------------------------
template <typename T>
constexpr T
to_degrees (T radians)
{
static_assert (std::is_floating_point<T>::value, "undefined for integral types");
return radians * 180 / PI<T>;
}
//-----------------------------------------------------------------------------
template <typename T>
constexpr T
to_radians (T degrees)
{
static_assert (std::is_floating_point<T>::value, "undefined for integral types");
return degrees / 180 * PI<T>;
}
//-----------------------------------------------------------------------------
//! Normalised sinc function
template <typename T>
constexpr T
sincn (T x)
{
return almost_zero (x) ? 1 : std::sin (PI<T> * x) / (PI<T> * x);
}
//-----------------------------------------------------------------------------
//! Unnormalised sinc function
template <typename T>
constexpr T
sincu (T x)
{
return almost_zero (x) ? 1 : std::sin (x) / x;
}
///////////////////////////////////////////////////////////////////////////////
// combinatorics
constexpr uintmax_t
factorial (unsigned i)
{
return i <= 1 ? 0 : i * factorial (i - 1);
}
//-----------------------------------------------------------------------------
/// stirlings approximation of factorials
constexpr uintmax_t
stirling (unsigned n)
{
return static_cast<uintmax_t> (
std::sqrt (2 * PI<float> * n) * std::pow (n / E<float>, n)
);
}
//-----------------------------------------------------------------------------
constexpr uintmax_t
combination (unsigned n, unsigned k)
{
return factorial (n) / (factorial (k) / (factorial (n - k)));
}
///////////////////////////////////////////////////////////////////////////////
// kahan summation for long floating point sequences
template <class InputIt>
typename std::iterator_traits<InputIt>::value_type
fsum (InputIt first, InputIt last)
{
using T = typename std::iterator_traits<InputIt>::value_type;
static_assert (std::is_floating_point<T>::value,
"fsum only works for floating point types");
T sum = 0;
T c = 0;
for (auto cursor = first; cursor != last; ++cursor) {
T y = *cursor - c;
T t = sum + y;
c = (t - sum) - y;
sum = t;
}
return sum;
}
///////////////////////////////////////////////////////////////////////////
/// Variadic minimum
template <typename T>
constexpr T
min (const T a)
@ -402,164 +394,164 @@ namespace util {
{
return max (a > b ? a : b, std::forward<Args> (args)...);
}
///////////////////////////////////////////////////////////////////////////
// Limiting functions
// min/max clamping
template <typename T, typename U, typename V>
constexpr T
limit (const T val, const U lo, const V hi)
{
lo <= hi ? (void)0 : panic ();
return val > hi ? hi:
val < lo ? lo:
val;
}
//-------------------------------------------------------------------------
// clamped cubic hermite interpolation
template <typename T>
T
smoothstep (T a, T b, T x)
{
CHECK_LE(a, b);
x = limit ((x - a) / (b - a), T{0}, T{1});
return x * x * (3 - 2 * x);
}
///////////////////////////////////////////////////////////////////////////
// renormalisation of unit floating point and/or normalised integers
// int -> float
template <typename T, typename U>
constexpr
typename std::enable_if<
!std::is_floating_point<T>::value && std::is_floating_point<U>::value, U
>::type
renormalise (T t)
{
return t / static_cast<U> (std::numeric_limits<T>::max ());
}
//-------------------------------------------------------------------------
// float -> int
template <typename T, typename U>
constexpr
typename std::enable_if<
std::is_floating_point<T>::value && !std::is_floating_point<U>::value, U
>::type
renormalise (T t)
{
// Ideally std::ldexp would be involved but it complicates handing
// integers with greater precision than our floating point type. Also it
// would prohibit constexpr and involve errno.
size_t usable = std::numeric_limits<T>::digits;
size_t available = sizeof (U) * 8;
size_t shift = std::max (available, usable) - usable;
t = limit (t, 0, 1);
// construct an integer of the float's mantissa size, multiply it by our
// parameter, then shift it back into the full range of the integer type.
U in = std::numeric_limits<U>::max () >> shift;
U mid = static_cast<U> (t * in);
U out = mid << shift;
// use the top bits of the output to fill the bottom bits which through
// shifting would otherwise be zero. this gives us the full extent of the
// integer range, while varying predictably through the entire output
// space.
return out | out >> (available - shift);
}
//-------------------------------------------------------------------------
// float -> float, avoid identity conversion as we don't want to create
// ambiguous overloads
template <typename T, typename U>
constexpr
typename std::enable_if<
std::is_floating_point<T>::value &&
std::is_floating_point<U>::value &&
!std::is_same<T,U>::value, U
>::type
renormalise (T t)
{
return static_cast<U> (t);
}
//-------------------------------------------------------------------------
// hi_int -> lo_int
template <typename T, typename U>
constexpr
typename std::enable_if<
std::is_integral<T>::value &&
std::is_integral<U>::value &&
(sizeof (T) > sizeof (U)), U
>::type
renormalise (T t)
{
static_assert (sizeof (T) > sizeof (U),
"assumes right shift is sufficient");
// we have excess bits ,just shift and return
constexpr auto shift = 8 * (sizeof (T) - sizeof (U));
return t >> shift;
}
//-------------------------------------------------------------------------
// lo_int -> hi_int
template <typename T, typename U>
constexpr
typename std::enable_if<
std::is_integral<T>::value &&
std::is_integral<U>::value &&
sizeof (T) < sizeof (U), U
>::type
renormalise (T t)
{
static_assert (sizeof (T) < sizeof (U),
"assumes bit creation is required to fill space");
// we need to create bits. fill the output integer with copies of ourself.
// this is approximately correct in the general case (introducing a small
// linear positive bias), but allows us to fill the output space in the
// case of input maximum.
static_assert (sizeof (U) % sizeof (T) == 0,
"assumes integer multiple of sizes");
U out = 0;
for (size_t i = 0; i < sizeof (U) / sizeof (T); ++i)
out |= U (t) << sizeof (T) * 8 * i;
return out;
}
//-------------------------------------------------------------------------
template <typename T, typename U>
constexpr
typename std::enable_if<
std::is_same<T,U>::value, U
>::type
renormalise (T t)
{ return t; }
}
///////////////////////////////////////////////////////////////////////////////
// Limiting functions
// min/max clamping
template <typename T, typename U, typename V>
constexpr T
limit (const T val, const U lo, const V hi)
{
lo <= hi ? (void)0 : panic ();
return val > hi ? hi:
val < lo ? lo:
val;
}
//-----------------------------------------------------------------------------
// clamped cubic hermite interpolation
template <typename T>
T
smoothstep (T a, T b, T x)
{
CHECK_LE(a, b);
x = limit ((x - a) / (b - a), T{0}, T{1});
return x * x * (3 - 2 * x);
}
#include "types/string.hpp"
///////////////////////////////////////////////////////////////////////////////
// renormalisation of unit floating point and/or normalised integers
// int -> float
template <typename T, typename U>
constexpr
typename std::enable_if<
!std::is_floating_point<T>::value && std::is_floating_point<U>::value, U
>::type
renormalise (T t)
{
return t / static_cast<U> (std::numeric_limits<T>::max ());
}
//-----------------------------------------------------------------------------
// float -> int
template <typename T, typename U>
constexpr
typename std::enable_if<
std::is_floating_point<T>::value && !std::is_floating_point<U>::value, U
>::type
renormalise (T t)
{
// Ideally std::ldexp would be involved but it complicates handing
// integers with greater precision than our floating point type. Also it
// would prohibit constexpr and involve errno.
size_t usable = std::numeric_limits<T>::digits;
size_t available = sizeof (U) * 8;
size_t shift = std::max (available, usable) - usable;
t = limit (t, 0, 1);
// construct an integer of the float's mantissa size, multiply it by our
// parameter, then shift it back into the full range of the integer type.
U in = std::numeric_limits<U>::max () >> shift;
U mid = static_cast<U> (t * in);
U out = mid << shift;
// use the top bits of the output to fill the bottom bits which through
// shifting would otherwise be zero. this gives us the full extent of the
// integer range, while varying predictably through the entire output
// space.
return out | out >> (available - shift);
}
//-----------------------------------------------------------------------------
// float -> float, avoid identity conversion as we don't want to create
// ambiguous overloads
template <typename T, typename U>
constexpr
typename std::enable_if<
std::is_floating_point<T>::value &&
std::is_floating_point<U>::value &&
!std::is_same<T,U>::value, U
>::type
renormalise (T t)
{
return static_cast<U> (t);
}
//-----------------------------------------------------------------------------
// hi_int -> lo_int
template <typename T, typename U>
constexpr
typename std::enable_if<
std::is_integral<T>::value &&
std::is_integral<U>::value &&
(sizeof (T) > sizeof (U)), U
>::type
renormalise (T t)
{
static_assert (sizeof (T) > sizeof (U),
"assumes right shift is sufficient");
// we have excess bits ,just shift and return
constexpr auto shift = 8 * (sizeof (T) - sizeof (U));
return t >> shift;
}
//-----------------------------------------------------------------------------
// lo_int -> hi_int
template <typename T, typename U>
constexpr
typename std::enable_if<
std::is_integral<T>::value &&
std::is_integral<U>::value &&
sizeof (T) < sizeof (U), U
>::type
renormalise (T t)
{
static_assert (sizeof (T) < sizeof (U),
"assumes bit creation is required to fill space");
// we need to create bits. fill the output integer with copies of ourself.
// this is approximately correct in the general case (introducing a small
// linear positive bias), but allows us to fill the output space in the
// case of input maximum.
static_assert (sizeof (U) % sizeof (T) == 0,
"assumes integer multiple of sizes");
U out = 0;
for (size_t i = 0; i < sizeof (U) / sizeof (T); ++i)
out |= U (t) << sizeof (T) * 8 * i;
return out;
}
//-----------------------------------------------------------------------------
template <typename T, typename U>
constexpr
typename std::enable_if<
std::is_same<T,U>::value, U
>::type
renormalise (T t)
{ return t; }
//#include "types/string.hpp"
#include "maths.ipp"
#endif // __MATHS_HPP

6
maths.ipp
View File

@ -42,7 +42,7 @@ util::pow (T x, unsigned y)
/// problems with constexpr under clang. If you need speed then you'll probably
/// have to handcode something.
constexpr int
sign (int v)
util::sign (int v)
{
return std::signbit (v) ? -1 : 1;
}
@ -50,7 +50,7 @@ sign (int v)
//-----------------------------------------------------------------------------
constexpr float
sign (float v)
util::sign (float v)
{
return std::signbit (v) ? -1.f : 1.f;
}
@ -58,7 +58,7 @@ sign (float v)
//-----------------------------------------------------------------------------
constexpr double
sign (double v)
util::sign (double v)
{
return std::signbit (v) ? -1. : 1.f;
}

4
test/colour.cpp
View File

@ -26,11 +26,11 @@ main (int, char**)
{
// white: hue is undefined
auto white = util::rgb_to_hsv ({1,1,1});
tap.expect (exactly_equal (white.s, 0) && exactly_equal (white.v, 1), "white hsv");
tap.expect (util::exactly_zero (white.s) && util::exactly_equal (white.v, 1), "white hsv");
// black: hue is undefined
auto black = util::rgb_to_hsv ({0,0,0});
tap.expect (exactly_equal (black.s, 0) && exactly_equal (black.v, 0), "black hsv");
tap.expect (util::exactly_zero (black.s) && util::exactly_zero (black.v), "black hsv");
struct {
const char *name;

1
test/fixed.cpp
View File

@ -1,4 +1,5 @@
#include "fixed.hpp"
#include "types/string.hpp"
#include "tap.hpp"

4
test/introspection.cpp
View File

@ -46,8 +46,8 @@ int main ()
foo d_foo { 7, 42.0 };
auto f_tuple = util::as_tuple (d_foo);
tap.expect (exactly_equal (d_foo.a, std::get<0> (f_tuple)) &&
exactly_equal (d_foo.b, std::get<1> (f_tuple)),
tap.expect (util::exactly_equal (d_foo.a, std::get<0> (f_tuple)) &&
util::exactly_equal (d_foo.b, std::get<1> (f_tuple)),
"dynamic member access after conversion to tuple");
}
}

4
test/json_types.cpp
View File

@ -45,7 +45,7 @@ main (void) {
CHECK (!ref["integer"].is_object ());
CHECK (!ref["integer"].is_string ());
CHECK (
exactly_equal (
util::exactly_equal (
(unsigned)ref["integer"].as_number ().native (),
1u
)
@ -81,7 +81,7 @@ main (void) {
CHECK (!ref["double"].is_object ());
CHECK (!ref["double"].is_string ());
CHECK (
exactly_equal (
util::exactly_equal (
ref["double"].as_number ().native (),
3.14
)

100
test/maths.cpp
View File

@ -15,49 +15,49 @@ void
test_comparisons (util::TAP::logger &tap)
{
// Check pos/neg zeroes
tap.expect (almost_equal ( 0.f, 0.f), "equal float zeros +ve/+ve");
tap.expect (almost_equal ( 0.f, -0.f), "equal float zeros +ve/-ve");
tap.expect (almost_equal (-0.f, 0.f), "equal float zeros -ve/+ve");
tap.expect (almost_equal (-0.f, -0.f), "equal float zeros -ve/-ve");
tap.expect (util::almost_equal ( 0.f, 0.f), "equal float zeros +ve/+ve");
tap.expect (util::almost_equal ( 0.f, -0.f), "equal float zeros +ve/-ve");
tap.expect (util::almost_equal (-0.f, 0.f), "equal float zeros -ve/+ve");
tap.expect (util::almost_equal (-0.f, -0.f), "equal float zeros -ve/-ve");
tap.expect (almost_equal ( 0., 0.), "equal double zeroes +ve/+ve");
tap.expect (almost_equal ( 0., -0.), "equal double zeroes +ve/+ve");
tap.expect (almost_equal (-0., 0.), "equal double zeroes +ve/+ve");
tap.expect (almost_equal (-0., -0.), "equal double zeroes +ve/+ve");
tap.expect (util::almost_equal ( 0., 0.), "equal double zeroes +ve/+ve");
tap.expect (util::almost_equal ( 0., -0.), "equal double zeroes +ve/+ve");
tap.expect (util::almost_equal (-0., 0.), "equal double zeroes +ve/+ve");
tap.expect (util::almost_equal (-0., -0.), "equal double zeroes +ve/+ve");
// Check zero comparison with values near the expected cutoff
tap.expect (almost_zero (1e-45f), "almost_zero with low value");
tap.expect (!almost_zero (1e-40f), "not almost_zero with low value");
tap.expect (!exactly_zero (1e-45f), "not exactly_zero with low value");
tap.expect (util::almost_zero (1e-45f), "almost_zero with low value");
tap.expect (!util::almost_zero (1e-40f), "not almost_zero with low value");
tap.expect (!util::exactly_zero (1e-45f), "not exactly_zero with low value");
// Compare values a little away from zero
tap.expect (!almost_equal (-2.0, 0.0), "not equal floats");
tap.expect (!almost_equal (-2.f, 0.f), "not equal doubles");
tap.expect (!util::almost_equal (-2.0, 0.0), "not equal floats");
tap.expect (!util::almost_equal (-2.f, 0.f), "not equal doubles");
// Compare values at the maximum extreme
tap.expect (!almost_equal (-std::numeric_limits<float>::max (), 0.f), "not equal -max/0 float");
tap.expect (!almost_equal (-std::numeric_limits<float>::max (),
std::numeric_limits<float>::max ()),
tap.expect (!util::almost_equal (-std::numeric_limits<float>::max (), 0.f), "not equal -max/0 float");
tap.expect (!util::almost_equal (-std::numeric_limits<float>::max (),
std::numeric_limits<float>::max ()),
"not equal -max/max");
// Compare infinity values
tap.expect ( almost_equal (numeric_limits<double>::infinity (),
numeric_limits<double>::infinity ()),
tap.expect ( util::almost_equal (numeric_limits<double>::infinity (),
numeric_limits<double>::infinity ()),
"almost_equal +infinity");
tap.expect (!almost_equal (numeric_limits<double>::infinity (), 0.0),
tap.expect (!util::almost_equal (numeric_limits<double>::infinity (), 0.0),
"not almost_equal +inf/0");
// Compare NaNs
tap.expect (!almost_equal (0., numeric_limits<double>::quiet_NaN ()), "not almost_equal double 0/NaN");
tap.expect (!almost_equal (numeric_limits<double>::quiet_NaN (), 0.), "not almost_equal double NaN/0");
tap.expect (!util::almost_equal (0., numeric_limits<double>::quiet_NaN ()), "not almost_equal double 0/NaN");
tap.expect (!util::almost_equal (numeric_limits<double>::quiet_NaN (), 0.), "not almost_equal double NaN/0");
tap.expect (!almost_equal (numeric_limits<double>::quiet_NaN (),
numeric_limits<double>::quiet_NaN ()),
tap.expect (!util::almost_equal (numeric_limits<double>::quiet_NaN (),
numeric_limits<double>::quiet_NaN ()),
"not almost_equal NaN/NaN");
// Compare reasonably close values that are wrong
tap.expect (!almost_equal (1.0000f, 1.0001f), ".0001f difference inequality");
tap.expect ( almost_equal (1.0000f, 1.00001f), ".00001f difference inequality");
tap.expect (!util::almost_equal (1.0000f, 1.0001f), ".0001f difference inequality");
tap.expect ( util::almost_equal (1.0000f, 1.00001f), ".00001f difference inequality");
}
@ -66,10 +66,10 @@ test_normalisations (util::TAP::logger &tap)
{
// u8 to float
{
auto a = renormalise<uint8_t,float> (255);
auto a = util::renormalise<uint8_t,float> (255);
tap.expect_eq (a, 1.f, "normalise uint8 max");
auto b = renormalise<uint8_t,float> (0);
auto b = util::renormalise<uint8_t,float> (0);
tap.expect_eq (b, 0.f, "normalise uint8 min");
}
@ -88,8 +88,8 @@ test_normalisations (util::TAP::logger &tap)
};
for (auto i: TESTS) {
auto v = renormalise<decltype(i.a),decltype(i.b)> (i.a);
success = success && almost_equal (unsigned{v}, unsigned{i.b});
auto v = util::renormalise<decltype(i.a),decltype(i.b)> (i.a);
success = success && util::almost_equal (unsigned{v}, unsigned{i.b});
}
tap.expect (success, "float-u8 normalisation");
@ -111,17 +111,17 @@ test_normalisations (util::TAP::logger &tap)
};
for (auto t: TESTS) {
auto v = renormalise<float,uint32_t> (t.a);
success = success && almost_equal (t.b, v);
auto v = util::renormalise<float,uint32_t> (t.a);
success = success && util::almost_equal (t.b, v);
}
tap.expect (success, "float-u32 normalisation");
}
tap.expect_eq (renormalise<uint8_t,uint32_t> (0xff), 0xffffffffu, "normalise hi u8-to-u32");
tap.expect_eq (renormalise<uint8_t,uint32_t> (0x00), 0x00000000u, "normalise lo u8-to-u32");
tap.expect_eq (util::renormalise<uint8_t,uint32_t> (0xff), 0xffffffffu, "normalise hi u8-to-u32");
tap.expect_eq (util::renormalise<uint8_t,uint32_t> (0x00), 0x00000000u, "normalise lo u8-to-u32");
tap.expect_eq (renormalise<uint32_t,uint8_t> (0xffffffff), 0xffu, "normalise hi u32-to-u8");
tap.expect_eq (util::renormalise<uint32_t,uint8_t> (0xffffffff), 0xffu, "normalise hi u32-to-u8");
}
@ -143,30 +143,30 @@ main (void)
tap.expect_eq (util::pow2 (4u), 16u, "pow2");
tap.expect_eq (rootsquare (2, 2), sqrt (8), "rootsquare");
tap.expect_eq (util::rootsquare (2, 2), sqrt (8), "rootsquare");
static const double POS_ZERO = 1.0 / numeric_limits<double>::infinity ();
static const double NEG_ZERO = -1.0 / numeric_limits<double>::infinity ();
tap.expect_eq (sign (-1), -1, "sign(-1)");
tap.expect_eq (sign ( 1), 1, "sign( 1)");
tap.expect_eq (sign (POS_ZERO), 1., "sign (POS_ZERO)");
tap.expect_eq (sign (NEG_ZERO), -1., "sign (NEG_ZERO)");
tap.expect_eq (sign ( numeric_limits<double>::infinity ()), 1., "sign +inf");
tap.expect_eq (sign (-numeric_limits<double>::infinity ()), -1., "sign -inf");
tap.expect_eq (util::sign (-1), -1, "sign(-1)");
tap.expect_eq (util::sign ( 1), 1, "sign( 1)");
tap.expect_eq (util::sign (POS_ZERO), 1., "sign (POS_ZERO)");
tap.expect_eq (util::sign (NEG_ZERO), -1., "sign (NEG_ZERO)");
tap.expect_eq (util::sign ( numeric_limits<double>::infinity ()), 1., "sign +inf");
tap.expect_eq (util::sign (-numeric_limits<double>::infinity ()), -1., "sign -inf");
tap.expect_eq (to_degrees (PI< float>), 180.f, "to_degrees float");
tap.expect_eq (to_degrees (PI<double>), 180.0, "to_degrees double");
tap.expect_eq (to_radians (180.f), PI<float>, "to_radians float");
tap.expect_eq (to_radians (180.0), PI<double>, "to_radians double");
tap.expect_eq (util::to_degrees (util::PI< float>), 180.f, "to_degrees float");
tap.expect_eq (util::to_degrees (util::PI<double>), 180.0, "to_degrees double");
tap.expect_eq (util::to_radians (180.f), util::PI<float>, "to_radians float");
tap.expect_eq (util::to_radians (180.0), util::PI<double>, "to_radians double");
tap.expect_eq (log2 (8u), 3u, "log2 +ve");
tap.expect_eq (log2 (1u), 0u, "log2 zero");
tap.expect_eq (util::log2 (8u), 3u, "log2 +ve");
tap.expect_eq (util::log2 (1u), 0u, "log2 zero");
//tap.expect_eq (log2 (9u), 3, "log2up 9");
tap.expect_eq (log2up (9u), 4u, "log2up 9");
tap.expect_eq (log2up (8u), 3u, "log2up 9");
tap.expect_eq (log2up (1u), 0u, "log2up 9");
tap.expect_eq (util::log2up (9u), 4u, "log2up 9");
tap.expect_eq (util::log2up (8u), 3u, "log2up 9");
tap.expect_eq (util::log2up (1u), 0u, "log2up 9");
return tap.status ();
}

12
test/matrix.cpp
View File

@ -32,10 +32,10 @@ main (void)
auto r = m * v;
tap.expect (
almost_equal (r.x, 30.f) &&
almost_equal (r.y, 70.f) &&
almost_equal (r.z, 110.f) &&
almost_equal (r.w, 150.f),
util::almost_equal (r.x, 30.f) &&
util::almost_equal (r.y, 70.f) &&
util::almost_equal (r.z, 110.f) &&
util::almost_equal (r.w, 150.f),
"simple matrix-vector multiplication"
);
}
@ -79,9 +79,9 @@ main (void)
for (size_t r = 0; r < m.rows; ++r)
for (size_t c = 0; c < m.cols; ++c)
if (r == c)
success = success && almost_equal (m.values[r][c], 1.f);
success = success && util::almost_equal (m.values[r][c], 1.f);
else
success = success && almost_equal (m.values[r][c], 0.f);
success = success && util::almost_equal (m.values[r][c], 0.f);
tap.expect (success, "identity inversion");
}

2
test/polynomial.cpp
View File

@ -48,7 +48,7 @@ main (int, char**)
continue;
}
if (!almost_equal (i.solutions[j], s[j]))
if (!util::almost_equal (i.solutions[j], s[j]))
ok = false;
}

10
test/vector.cpp
View File

@ -28,13 +28,13 @@ test_polar (util::TAP::logger &tap)
},
{
{ 1.f, PI<float> / 2.f },
{ 1.f, util::PI<float> / 2.f },
{ 0.f, 1.f },
"unit length, rotated"
},
{
{ 1.f, 2 * PI<float> },
{ 1.f, 2 * util::PI<float> },
{ 1.f, 0.f },
"full rotation, unit length"
}
@ -53,8 +53,8 @@ test_polar (util::TAP::logger &tap)
auto in_polar = t.polar;
auto to_polar = util::cartesian_to_polar (t.cartesian);
in_polar[1] = std::fmod (in_polar[1], 2 * PI<float>);
to_polar[1] = std::fmod (to_polar[1], 2 * PI<float>);
in_polar[1] = std::fmod (in_polar[1], 2 * util::PI<float>);
to_polar[1] = std::fmod (to_polar[1], 2 * util::PI<float>);
tap.expect_eq (in_polar, to_polar, t.desc);
}
@ -83,7 +83,7 @@ test_euler (util::TAP::logger &tap)
// check that simple axis rotations look correct
for (auto i: TESTS) {
tap.expect_eq (util::to_euler (i.dir),
i.euler * PI<float>,
i.euler * util::PI<float>,
"to euler, %s", i.name);
}

15
vector.cpp
View File

@ -319,11 +319,14 @@ namespace util {
template <> vector<4,double> random (void) { util::vector<4,double> out; randomise (out.data); return out; }
}
template <>
bool
almost_equal [[gnu::pure]] (const util::vector2f &a, const util::vector2f &b)
{
bool (*comparator) (const float&, const float&) = almost_equal;
return std::equal (a.begin (), a.end (), b.begin (), comparator);
namespace util {
template <>
bool
almost_equal [[gnu::pure]] (const util::vector2f &a, const util::vector2f &b)
{
bool (*comparator) (const float&, const float&) = almost_equal;
return std::equal (a.begin (), a.end (), b.begin (), comparator);
}
}