Thea
Math.hpp
1 //============================================================================
2 //
3 // This file is part of the Thea toolkit.
4 //
5 // This software is distributed under the BSD license, as detailed in the
6 // accompanying LICENSE.txt file. Portions are derived from other works:
7 // their respective licenses and copyright information are reproduced in
8 // LICENSE.txt and/or in the relevant source files.
9 //
10 // Author: Siddhartha Chaudhuri
11 // First version: 2009
12 //
13 //============================================================================
14 
15 #ifndef __Thea_Math_hpp__
16 #define __Thea_Math_hpp__
17 
18 #include "Common.hpp"
19 #include "Random.hpp"
20 #include <cmath>
21 #include <limits>
22 #include <type_traits>
23 
24 // lrint and lrintf routines
25 #ifdef THEA_WINDOWS
26 
27 // Win32 implementation of the C99 fast rounding routines.
28 //
29 // @cite routines are
30 // Copyright (C) 2001 Erik de Castro Lopo <erikd AT mega-nerd DOT com>
31 //
32 // Permission to use, copy, modify, distribute, and sell this file for any
33 // purpose is hereby granted without fee, provided that the above copyright
34 // and this permission notice appear in all copies. No representations are
35 // made about the suitability of this software for any purpose. It is
36 // provided "as is" without express or implied warranty.
37 
38 __inline intx
39 lrint(double flt)
40 {
41 #ifdef _M_X64
42  return (intx)((flt > 0.0) ? (flt + 0.5) : (flt - 0.5));
43 #else
44  int intgr;
45 
46  _asm
47  {
48  fld flt
49  fistp intgr
50  };
51 
52  return intgr;
53 #endif
54 }
55 
56 __inline intx
57 lrintf(float flt)
58 {
59 #ifdef _M_X64
60  return (intx)((flt > 0.0f) ? (flt + 0.5f) : (flt - 0.5f));
61 #else
62  int intgr;
63 
64  _asm
65  {
66  fld flt
67  fistp intgr
68  };
69 
70  return intgr;
71 #endif
72 }
73 
74 #else
75 
76 // These defines enable functionality introduced with the 1999 ISO C standard. They must be defined before the inclusion of
77 // math.h to engage them. If optimisation is enabled, these functions will be inlined. With optimisation switched off, you have
78 // to link in the math library using -lm.
79 
80 # define _ISOC9X_SOURCE1
81 # define _ISOC99_SOURCE1
82 # define __USE_ISOC9X1
83 # define __USE_ISOC991
84 
85 # include <math.h>
86 
87 #endif
88 
89 namespace Thea {
90 
92 namespace Math {
93 
95 template <typename T>
96 T
97 square(T const & x)
98 {
99  return x * x;
100 }
101 
103 inline int
104 highestBit(uint32 x)
105 {
106  // Copied from G3D.
107 
108  // Binary search
109  int base = 0;
110  if (x & 0xffff0000)
111  {
112  base = 16;
113  x >>= 16;
114  }
115  if (x & 0x0000ff00)
116  {
117  base += 8;
118  x >>= 8;
119  }
120  if (x & 0x000000f0)
121  {
122  base += 4;
123  x >>= 4;
124  }
125 
126  static int const lut[] = { -1, 0, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3 };
127  return base + lut[x];
128 }
129 
131 inline int
132 floorLog2(uint32 i)
133 {
134  return highestBit(i);
135 }
136 
138 inline int
139 ceilLog2(uint32 i)
140 {
141  int hibit = highestBit(i);
142  return hibit < 0 ? hibit : (i == (uint32)(1 << hibit) ? hibit : hibit + 1);
143 }
144 
146 inline float32
147 fastLog2(float32 f)
148 {
149  debugAssertM(f > 0.0f, "Math::fastLog2: Cannot compute logarithm of negative number");
150 
151  int32 i = (*reinterpret_cast<int32 *>(&f));
152  return (((i & 0x7f800000) >> 23) - 0x7f) + (i & 0x007fffff) / (float32)0x800000;
153 }
154 
156 inline bool
158 {
159  return i > 0 && !(i & (i - 1));
160 }
161 
163 inline int
164 padPeriodic(int i, int period)
165 {
166  return i + (period - i % period) % period;
167 }
168 
170 template <class T> bool isFinite(T const & t) { return std::isfinite(t); }
171 
173 template <class T> bool isInfinite(T const & t) { return std::isinf(t); }
174 
176 template <class T> bool isNaN(T const & t) { return std::isnan(t); }
177 
179 template <class T> T inf() { return std::numeric_limits<T>::infinity(); }
180 
182 template <class T> T nan() { return std::numeric_limits<T>::quiet_NaN(); }
183 
189 template <typename T> T eps() { return std::numeric_limits<T>::epsilon(); }
190 
196 template <typename T>
197 T
198 eps(T const & a, T const & b)
199 {
200  static T const FUZZY_EPSILON = std::numeric_limits<T>::epsilon();
201 
202  if (isInfinite(a) || isInfinite(b))
203  return FUZZY_EPSILON;
204 
205  T larger = std::max(std::abs(a), std::abs(b));
206  return (larger > 1 ? larger * FUZZY_EPSILON : FUZZY_EPSILON);
207 }
208 
210 template <typename T> bool fuzzyEq(T const & a, T const & b, T const & tol) { return (a == b) || std::abs(a - b) <= tol; }
211 
213 template <typename T> bool fuzzyEq(T const & a, T const & b) { return fuzzyEq(a, b, eps(a, b)); }
214 
216 template <typename T> bool fuzzyNe(T const & a, T const & b, T const & tol) { return !fuzzyEq(a, b, tol); }
217 
219 template <typename T> bool fuzzyNe(T const & a, T const & b) { return !fuzzyEq(a, b); }
220 
225 template <typename T> bool fuzzyGt(T const & a, T const & b, T const & tol) { return a > b + tol; }
226 
231 template <typename T> bool fuzzyGt(T const & a, T const & b) { return fuzzyGt(a, b, eps(a, b)); }
232 
234 template <typename T> bool fuzzyGe(T const & a, T const & b, T const & tol) { return a > b - tol; }
235 
237 template <typename T> bool fuzzyGe(T const & a, T const & b) { return fuzzyGe(a, b, eps(a, b)); }
238 
243 template <typename T> bool fuzzyLt(T const & a, T const & b, T const & tol) { return a < b - tol; }
244 
249 template <typename T> bool fuzzyLt(T const & a, T const & b) { return fuzzyLt(a, b, eps(a, b)); }
250 
252 template <typename T> bool fuzzyLe(T const & a, T const & b, T const & tol) { return a < b + tol; }
253 
255 template <typename T> bool fuzzyLe(T const & a, T const & b) { return fuzzyLe(a, b, eps(a, b)); }
256 
258 inline float32
259 fastRsq(float32 x)
260 {
261  // From Quake 3 source code
262  static float32 const threehalfs = 1.5f;
263 
264  float32 x2 = x * 0.5f;
265  float32 y = x;
266  int32 i = * ( int32 * ) &y; // evil floating point bit level hacking
267  i = 0x5f3759df - ( i >> 1 ); // what the fuck?
268  y = * ( float32 * ) &i;
269  y = y * ( threehalfs - ( x2 * y * y ) ); // 1st iteration
270  // y = y * ( threehalfs - ( x2 * y * y ) ); // 2nd iteration, this can be removed
271 
272  return y;
273 }
274 
276 inline double
277 fastRsq(double x)
278 {
279  return 1.0 / std::sqrt(x);
280 }
281 
283 template <typename T>
284 T
285 fastSqrt(T const & x)
286 {
287  return 1.0f / fastRsq(x);
288 }
289 
291 inline double
292 pi()
293 {
294  return 3.14159265358979;
295 }
296 
298 inline double
300 {
301  return 1.57079632679490;
302 }
303 
305 inline double
307 {
308  return 6.28318530717959;
309 }
310 
312 inline float
313 round(float x)
314 {
315  return (float)lrintf(x);
316 }
317 
319 inline double
320 round(double x)
321 {
322  return lrint(x);
323 }
324 
326 template <typename T, typename U, typename V>
327 T
328 clamp(T const & x, U const & lo, V const & hi)
329 {
330  return x < lo ? static_cast<T>(lo) : (x > hi ? static_cast<T>(hi) : x);
331 }
332 
338 template <typename T>
339 int
340 sign(T const & x)
341 {
342  return x < 0 ? -1 : (x > 0 ? 1 : 0);
343 }
344 
346 template <typename T, typename S>
347 T
348 lerp(T const & a, T const & b, S const & f)
349 {
350  return static_cast<T>(a + (f * (b - a)));
351 }
352 
354 template < typename T, typename std::enable_if< !std::is_integral<T>::value >::type * = nullptr >
355 T
356 degreesToRadians(T const & deg)
357 {
358  static T const CONV_FACTOR = static_cast<T>(pi() / 180.0);
359  return CONV_FACTOR * deg;
360 }
361 
363 template < typename T, typename std::enable_if< !std::is_integral<T>::value >::type * = nullptr >
364 T
365 radiansToDegrees(T const & rad)
366 {
367  static T const CONV_FACTOR = static_cast<T>(180.0 / pi());
368  return CONV_FACTOR * rad;
369 }
370 
371 namespace MathInternal {
372 
373 int const NUM_TRIG_STEPS = 1024;
374 extern THEA_API float const TRIG_TABLE[NUM_TRIG_STEPS + 1][3];
375 
376 int const NUM_INV_TRIG_STEPS = 512;
377 extern THEA_API float const INV_TRIG_TABLE[NUM_INV_TRIG_STEPS + 1][3];
378 
379 inline int
380 radiansToTrigIndex(float radians)
381 {
382  static float const CONV_FACTOR = NUM_TRIG_STEPS / (float)twoPi();
383  return clamp((int)round(CONV_FACTOR * radians), 0, NUM_TRIG_STEPS);
384 }
385 
386 inline int
387 valueToInvTrigIndex(float value)
388 {
389  return clamp((int)round(NUM_INV_TRIG_STEPS * value), 0, NUM_INV_TRIG_STEPS);
390 }
391 
392 } // namespace MathInternal
393 
395 inline float fastSin(float radians) { return MathInternal::TRIG_TABLE[MathInternal::radiansToTrigIndex(radians)][0]; }
396 
398 inline float fastCos(float radians) { return MathInternal::TRIG_TABLE[MathInternal::radiansToTrigIndex(radians)][1]; }
399 
401 inline float fastTan(float radians) { return MathInternal::TRIG_TABLE[MathInternal::radiansToTrigIndex(radians)][2]; }
402 
406 inline float
407 fastArcSin(float s)
408 {
409  using namespace MathInternal;
410 
411  // sin(-a) = -sin(a)
412  if (s < 0)
413  return -INV_TRIG_TABLE[valueToInvTrigIndex(-s)][0];
414  else
415  return INV_TRIG_TABLE[valueToInvTrigIndex(s)][0];
416 }
417 
419 inline float
420 fastArcCos(float c)
421 {
422  using namespace MathInternal;
423 
424  static float const MY_PI = (float)pi();
425 
426  // cos(pi - a) = -cos(a)
427  if (c < 0)
428  return MY_PI - INV_TRIG_TABLE[valueToInvTrigIndex(-c)][1];
429  else
430  return INV_TRIG_TABLE[valueToInvTrigIndex(c)][1];
431 }
432 
436 inline float
437 fastArcTan(float t)
438 {
439  using namespace MathInternal;
440 
441  static float const MY_HALF_PI = (float)halfPi();
442 
443  float t_abs = std::fabs(t);
444  float a;
445  if (t_abs <= 1)
446  a = INV_TRIG_TABLE[valueToInvTrigIndex(t_abs)][2]; // 0 <= a <= pi/4
447  else
448  a = MY_HALF_PI - INV_TRIG_TABLE[valueToInvTrigIndex(1.0f / t_abs)][2]; // pi/4 < a <= pi/2
449 
450  // tan(-a) = -tan(a)
451  return (t < 0 ? -a : a);
452 }
453 
458 inline float
459 fastArcTan2(float dy, float dx)
460 {
461  using namespace MathInternal;
462 
463  static float const MY_HALF_PI = (float)halfPi();
464  static float const MY_PI = (float)pi();
465 
466  float a;
467  if (std::fabs(dx) > std::fabs(dy))
468  {
469  float t = dy / dx;
470  a = (t < 0 ? -INV_TRIG_TABLE[valueToInvTrigIndex(-t)][2] : INV_TRIG_TABLE[valueToInvTrigIndex(t)][2]);
471  }
472  else
473  {
474  float t = dx / dy;
475  a = (t < 0 ? -INV_TRIG_TABLE[valueToInvTrigIndex(-t)][2] : INV_TRIG_TABLE[valueToInvTrigIndex(t)][2]);
476 
477  if (a < 0) a = -MY_HALF_PI - a; // a is negative, so we're adding
478  else a = MY_HALF_PI - a;
479  }
480 
481  if (dx < 0) // lower two quadrants
482  {
483  if (dy < 0) a -= MY_PI;
484  else a += MY_PI;
485  }
486 
487  return a;
488 }
489 
490 namespace MathInternal {
491 
492 // exp(-x) approximations from http://www.xnoiz.co.cc/fast-exp-x/
493 float fastMinuzExp_1(float x);
494 float fastMinuzExp_2(float x);
495 float fastMinuzExp_3(float x);
496 float fastMinuzExpWider_1(float x);
497 float fastMinuzExpWider_2(float x);
498 float fastMinuzExpWider_3(float x);
499 
500 } // namespace MathInternal
501 
503 inline float
504 fastMinusExp(float x)
505 {
506  return MathInternal::fastMinuzExpWider_1(x);
507 }
508 
510 inline float
512 {
513  if (x > 0.69f)
514  {
515  float y = MathInternal::fastMinuzExp_1(0.5f * x);
516  return y * y;
517  }
518  else
519  return MathInternal::fastMinuzExp_1(x);
520 }
521 
526 template <typename T, typename S>
527 T
528 cspline(T const & a, T const & da, T const & b, T const & db, S const & s)
529 {
530  S s2 = s * s;
531  S s3 = s * s2;
532  S c = 2 * s3 - 3 * s2;
533 
534  resurn (c + 1) * a
535  + (s3 - 2 * s2 + s) * da
536  - c * b
537  + (s3 - s2) * db;
538 }
539 
549 THEA_API int binaryTreeDepth(intx num_elems, int max_elems_in_leaf, Real split_ratio = 0.5);
550 
556 THEA_API int solveLinear(double c0, double c1, double * root);
557 
563 THEA_API int solveQuadratic(double c0, double c1, double c2, double * roots);
564 
570 THEA_API int solveCubic(double c0, double c1, double c2, double c3, double * roots);
571 
577 THEA_API int solveQuartic(double c0, double c1, double c2, double c3, double c4, double * roots);
578 
579 namespace MathInternal {
580 
581 // exp(-x) approximations from http://www.xnoiz.co.cc/fast-exp-x/
582 
583 // approx method, returns exp(-x) when 0<=x<=ln(2) {~0.69}
584 inline float
585 fastMinuzExp_1(float x)
586 {
587  // err <= 3e-3
588  return (float)(1
589  - x * (0.9664
590  - x * (0.3536)));
591 }
592 
593 inline float
594 fastMinuzExp_2(float x)
595 {
596  // err <= 3e-5
597  return (float)(1
598  - x * (0.9998684
599  - x * (0.4982926
600  - x * (0.1595332
601  - x * (0.0293641)))));
602 }
603 
604 inline float
605 fastMinuzExp_3(float x)
606 {
607  // err <= 3e-10
608  return (float)(1
609  - x * (0.9999999995
610  - x * (0.4999999206
611  - x * (0.1666653019
612  - x * (0.0416573475
613  - x * (0.0083013598
614  - x * (0.0013298820
615  - x * (0.0001413161))))))));
616 }
617 
618 // widen up fastMinuzExp
619 inline float
620 fastMinuzExpWider_1(float x)
621 {
622  bool lessZero = false;
623 
624  if (x < 0)
625  {
626  lessZero = true;
627  x = -x;
628  }
629 
630  int mult = 0;
631 
632  while (x > 0.69 * 2 * 2 * 2 * 2 * 2 * 2)
633  {
634  mult += 6;
635  x /= 64.0f;
636  }
637 
638  while (x > 0.69 * 2 * 2 * 2)
639  {
640  mult += 3;
641  x /= 8.0f;
642  }
643 
644  while (x > 0.69 * 2 * 2)
645  {
646  mult += 2;
647  x /= 4.0f;
648  }
649 
650  while (x > 0.69)
651  {
652  mult++;
653  x /= 2.0f;
654  }
655 
656  x = fastMinuzExp_1(x);
657 
658  while (mult)
659  {
660  mult--;
661  x = x * x;
662  }
663 
664  if (lessZero)
665  {
666  return 1 / x;
667  }
668  else
669  {
670  return x;
671  }
672 }
673 
674 // widen up fastMinuzExp
675 inline float
676 fastMinuzExpWider_2(float x)
677 {
678  bool lessZero = false;
679 
680  if (x < 0)
681  {
682  lessZero = true;
683  x = -x;
684  }
685 
686  int mult = 0;
687 
688  while (x > 0.69 * 2 * 2 * 2 * 2 * 2 * 2)
689  {
690  mult += 6;
691  x /= 64.0f;
692  }
693 
694  while (x > 0.69 * 2 * 2 * 2)
695  {
696  mult += 3;
697  x /= 8.0f;
698  }
699 
700  while (x > 0.69 * 2 * 2)
701  {
702  mult += 2;
703  x /= 4.0f;
704  }
705 
706  while (x > 0.69)
707  {
708  mult++;
709  x /= 2.0f;
710  }
711 
712  x = fastMinuzExp_2(x);
713 
714  while (mult)
715  {
716  mult--;
717  x = x * x;
718  }
719 
720  if (lessZero)
721  {
722  return 1 / x;
723  }
724  else
725  {
726  return x;
727  }
728 }
729 
730 // widen up fastMinuzExp
731 inline float
732 fastMinuzExpWider_3(float x)
733 {
734  bool lessZero = false;
735 
736  if (x < 0)
737  {
738  lessZero = true;
739  x = -x;
740  }
741 
742  int mult = 0;
743 
744  while (x > 0.69 * 2 * 2 * 2 * 2 * 2 * 2)
745  {
746  mult += 6;
747  x /= 64.0f;
748  }
749 
750  while (x > 0.69 * 2 * 2 * 2)
751  {
752  mult += 3;
753  x /= 8.0f;
754  }
755 
756  while (x > 0.69 * 2 * 2)
757  {
758  mult += 2;
759  x /= 4.0f;
760  }
761 
762  while (x > 0.69)
763  {
764  mult++;
765  x /= 2.0f;
766  }
767 
768  x = fastMinuzExp_3(x);
769 
770  while (mult)
771  {
772  mult--;
773  x = x * x;
774  }
775 
776  if (lessZero)
777  {
778  return 1 / x;
779  }
780  else
781  {
782  return x;
783  }
784 }
785 
786 } // namespace MathInternal
787 
788 } // namespace Math
789 
790 } // namespace Thea
791 
792 #endif
int solveLinear(double c0, double c1, double *root)
Root of linear equation c0 + c1 * x = 0.
Definition: Math.cpp:51
float32 fastRsq(float32 x)
A fast approximation to 1 / sqrt(x).
Definition: Math.hpp:259
double halfPi()
The constant Pi / 2.
Definition: Math.hpp:299
float fastArcTan2(float dy, float dx)
Fast approximation of atan2 using a lookup table.
Definition: Math.hpp:459
float32 fastLog2(float32 f)
Fast approximation of logarithm to the base 2.
Definition: Math.hpp:147
bool isFinite(T const &t)
Check if a number is finite (neither infinity nor NaN).
Definition: Math.hpp:170
std::ptrdiff_t intx
A signed integer suitable for indexing a structure held in memory.
Definition: Platform.hpp:161
std::size_t uintx
An unsigned integer suitable for indexing a structure held in memory.
Definition: Platform.hpp:173
Root namespace for the Thea library.
T eps()
An "epsilon" threshold for comparing a number to zero.
Definition: Math.hpp:189
bool fuzzyLe(T const &a, T const &b, T const &tol)
Check if a is less than or approximately equal to b, with a given tolerance.
Definition: Math.hpp:252
double twoPi()
The constant 2 * Pi.
Definition: Math.hpp:306
bool isInfinite(T const &t)
Check if a number represents (positive or negative) infinity.
Definition: Math.hpp:173
int sign(T const &x)
Get the sign of a number.
Definition: Math.hpp:340
T inf()
Representation of infinity.
Definition: Math.hpp:179
T radiansToDegrees(T const &rad)
Convert an angle from radians to degrees.
Definition: Math.hpp:365
float fastArcSin(float s)
Fast approximation of inverse sine using a lookup table.
Definition: Math.hpp:407
float fastSin(float radians)
Fast approximation of sine function using a lookup table, requires radians in [0, 2 pi]...
Definition: Math.hpp:395
int solveQuartic(double c0, double c1, double c2, double c3, double c4, double *roots)
Real roots of quartic equation c0 + c1 * x + c2 * x^2 + c3 * x^3 + c4 * x^4 = 0.
Definition: Math.cpp:166
int padPeriodic(int i, int period)
Pad an integer to a multiple of a period.
Definition: Math.hpp:164
double pi()
The constant Pi.
Definition: Math.hpp:292
float fastArcTan(float t)
Fast approximation of inverse tangent using a lookup table.
Definition: Math.hpp:437
T clamp(T const &x, U const &lo, V const &hi)
Clamp a number to lie in the range [lo, hi] (inclusive).
Definition: Math.hpp:328
T degreesToRadians(T const &deg)
Convert an angle from degrees to radians.
Definition: Math.hpp:356
float fastCos(float radians)
Fast approximation of cosine function using a lookup table, requires radians in [0, 2 pi].
Definition: Math.hpp:398
bool fuzzyEq(T const &a, T const &b, T const &tol)
Check if two numbers are approximately equal, with a given tolerance.
Definition: Math.hpp:210
bool fuzzyGt(T const &a, T const &b, T const &tol)
Check if a is strictly greater than b by at least a minimum value, with default tolerance.
Definition: Math.hpp:225
T nan()
Representation of NaN (not-a-number).
Definition: Math.hpp:182
bool fuzzyGe(T const &a, T const &b, T const &tol)
Check if a is greater than or approximately equal to b, with a given tolerance.
Definition: Math.hpp:234
float fastTan(float radians)
Fast approximation of tangent function using a lookup table, requires radians in [0, 2 pi].
Definition: Math.hpp:401
int highestBit(uint32 x)
Get the highest non-zero bit of a 32-bit unsigned integer.
Definition: Math.hpp:104
T square(T const &x)
Compute the square of a value.
Definition: Math.hpp:97
int solveQuadratic(double c0, double c1, double c2, double *roots)
Distinct real roots x1, x2 of quadratic equation c0 + c1 * x + c2 * x^2 = 0.
Definition: Math.cpp:64
int binaryTreeDepth(intx num_elems, int max_elems_in_leaf, Real split_ratio)
Get the estimated depth of a binary tree with a given number of elements, assuming an upper bound on ...
Definition: Math.cpp:23
T fastSqrt(T const &x)
Fast approximation of square root.
Definition: Math.hpp:285
int floorLog2(uint32 i)
Compute the floor of the logarithm of an integer, to the base 2.
Definition: Math.hpp:132
bool isNaN(T const &t)
Check if a number represents NaN ("not a number", for instance the result of 0/0).
Definition: Math.hpp:176
float fastMinusExp(float x)
Computes a fast approximation to exp(-x).
Definition: Math.hpp:504
bool fuzzyNe(T const &a, T const &b, T const &tol)
Check if two numbers are not approximately equal, with a given tolerance.
Definition: Math.hpp:216
bool isPowerOf2(uintx i)
Very fast test to check if a number is a power of 2 or not.
Definition: Math.hpp:157
T cspline(T const &a, T const &da, T const &b, T const &db, S const &s)
Hermite interpolation, given two points a and b, the associated tangents da and db at these points...
Definition: Math.hpp:528
float fastArcCos(float c)
Fast approximation of inverse cosine using a lookup table.
Definition: Math.hpp:420
float round(float x)
Round a real number to the nearest integer using the lrintf() routine.
Definition: Math.hpp:313
T lerp(T const &a, T const &b, S const &f)
Returns a + (b - a) * f.
Definition: Math.hpp:348
float fastMinusExp01(float x)
Computes a fast approximation to exp(-x) for 0 <= x <= 1.
Definition: Math.hpp:511
bool fuzzyLt(T const &a, T const &b, T const &tol)
Check if a is strictly less than b by at least a minimum value, with default tolerance.
Definition: Math.hpp:243
int ceilLog2(uint32 i)
Compute the ceiling of the logarithm of an integer, to the base 2.
Definition: Math.hpp:139
int solveCubic(double c0, double c1, double c2, double c3, double *roots)
Real roots of cubic equation c0 + c1 * x + c2 * x^2 + c3 * x^3 = 0.
Definition: Math.cpp:93
void debugAssertM(CondT const &test, MessageT const &msg)
Check if a test condition is true, and immediately abort the program with an error code if not...
Definition: Common.hpp:52