1 /** 2 Math Functions 3 4 Important: 5 When possible compiler level intrinsics will be used, but in some cases 6 this is not possible; in those cases we opt to focus on readability over performance. 7 8 As such these math functions may not be the most performant. 9 10 Copyright: Copyright The D Language Foundation 2000 - 2011. 11 Copyright: 12 Copyright © 2023-2025, Kitsunebi Games 13 Copyright © 2023-2025, Inochi2D Project 14 15 License: $(HTTP www.boost.org/LICENSE_1_0.txt, Boost License 1.0). 16 Authors: 17 $(HTTP digitalmars.com, Walter Bright), 18 Don Clugston, 19 Luna Nielsen 20 21 Macros: 22 SUB = $1<sub>$2</sub> 23 PI = π 24 SQRT = √ 25 HALF = ½ 26 */ 27 module nulib.math; 28 public import nulib.math.intrinsics; 29 public import nulib.math.floating; 30 31 // Values obtained from Wolfram Alpha. 116 bits ought to be enough for anybody. 32 // Wolfram Alpha LLC. 2011. Wolfram|Alpha. http://www.wolframalpha.com/input/?i=e+in+base+16 (access July 6, 2011). 33 enum real E = 0x1.5bf0a8b1457695355fb8ac404e7a8p+1L; /** e = 2.718281... */ 34 enum real LOG2T = 0x1.a934f0979a3715fc9257edfe9b5fbp+1L; /** $(SUB log, 2)10 = 3.321928... */ 35 enum real LOG2E = 0x1.71547652b82fe1777d0ffda0d23a8p+0L; /** $(SUB log, 2)e = 1.442695... */ 36 enum real LOG2 = 0x1.34413509f79fef311f12b35816f92p-2L; /** $(SUB log, 10)2 = 0.301029... */ 37 enum real LOG10E = 0x1.bcb7b1526e50e32a6ab7555f5a67cp-2L; /** $(SUB log, 10)e = 0.434294... */ 38 enum real LN2 = 0x1.62e42fefa39ef35793c7673007e5fp-1L; /** ln 2 = 0.693147... */ 39 enum real LN10 = 0x1.26bb1bbb5551582dd4adac5705a61p+1L; /** ln 10 = 2.302585... */ 40 enum real PI = 0x1.921fb54442d18469898cc51701b84p+1L; /** π = 3.141592... */ 41 enum real PI_2 = PI/2; /** $(PI) / 2 = 1.570796... */ 42 enum real PI_4 = PI/4; /** $(PI) / 4 = 0.785398... */ 43 enum real M_1_PI = 0x1.45f306dc9c882a53f84eafa3ea69cp-2L; /** 1 / $(PI) = 0.318309... */ 44 enum real M_2_PI = 2*M_1_PI; /** 2 / $(PI) = 0.636619... */ 45 enum real M_2_SQRTPI = 0x1.20dd750429b6d11ae3a914fed7fd8p+0L; /** 2 / $(SQRT)$(PI) = 1.128379... */ 46 enum real SQRT2 = 0x1.6a09e667f3bcc908b2fb1366ea958p+0L; /** $(SQRT)2 = 1.414213... */ 47 enum real SQRT1_2 = SQRT2/2; /** $(SQRT)$(HALF) = 0.707106... */ 48 49 @safe pure: 50 51 /** 52 Returns the smaller of the 2 given scalar values. 53 54 Params: 55 rhs = value 56 lhs = value 57 58 Returns: 59 The smallest of the 2 given values. 60 */ 61 T min(T)(T lhs, T rhs) { 62 return lhs < rhs ? lhs : rhs; 63 } 64 65 /** 66 Returns the larger of the 2 given scalar values. 67 68 Params: 69 rhs = value 70 lhs = value 71 72 Returns: 73 The largest of the 2 given values. 74 */ 75 T max(T)(T lhs, T rhs) { 76 return lhs > rhs ? lhs : rhs; 77 } 78 79 /** 80 Clamps scalar value into the given range. 81 82 Params: 83 value = The value to clamp, 84 min_ = The minimum value 85 max_ = The maximum value. 86 87 Returns: 88 $(D value) clamped between $(D min_) and $(D max_), 89 equivalent of $(D min(max(value, min_), max_)) 90 */ 91 T clamp(T)(T value, T min_, T max_) { 92 return min(max(value, min_), max_); 93 } 94 95 /** 96 Modulates the given value, preserving sign bit. 97 98 Params: 99 value = The value to modulate. 100 delta = The modulation delta. 101 102 Returns: 103 The modulated value. 104 */ 105 pragma(inline, true) 106 auto mod(T)(T value, T delta) { 107 return copysign(abs(value) % abs(delta), value); 108 } 109 110 /** 111 Linearly interpolates between $(D a) and $(D b) 112 113 Params: 114 a = The first value to interpolate 115 b = The second value to interpolate 116 t = The interpolation step from 0..1 117 118 Returns: 119 The interpolated value between $(D a) and $(D b) 120 */ 121 T lerp(T)(T a, T b, float t) { 122 return a * (1 - t) + b * t; 123 } 124 125 /** 126 Quadilaterally interpolates between $(D p0) and $(D p2), 127 with $(D p1) as a control point. 128 129 Params: 130 p0 = The first value to interpolate 131 p1 = The control value for the curve. 132 p2 = The second value to interpolate 133 t = The interpolation step from 0..1 134 135 Returns: 136 The interpolated value between $(D p0) and $(D p2) 137 */ 138 T quad(T)(T p0, T p1, T p2, float t) { 139 float tm = 1.0 - t; 140 float a = tm * tm; 141 float b = 2.0 * tm * t; 142 float c = t * t; 143 144 return a * p0 + b * p1 + c * p2; 145 } 146 147 /** 148 Interpolates between $(D p0) and $(D p3), using a cubic 149 spline with $(D p1) and $(D p2) as control points. 150 151 Params: 152 p0 = The first value to interpolate 153 p1 = The first control value for the curve. 154 p2 = The second control value for the curve. 155 p3 = The second value to interpolate 156 t = The interpolation step from 0..1 157 158 Returns: 159 The interpolated value between $(D p0) and $(D p3) 160 */ 161 T cubic(T)(T p0, T p1, T p2, T p3, float t) { 162 T a = -0.5 * p0 + 1.5 * p1 - 1.5 * p2 + 0.5 * p3; 163 T b = p0 - 2.5 * p1 + 2 * p2 - 0.5 * p3; 164 T c = -0.5 * p0 + 0.5 * p2; 165 T d = p1; 166 167 return a * (t ^^ 3) + b * (t ^^ 2) + c * t + d; 168 } 169 170 /** 171 Gets whether the value's sign bit is set. 172 173 Params: 174 x = The value to check 175 176 Returns: 177 $(D true) if the value is signed (positive), 178 $(D false) otherwise. 179 */ 180 bool signbit(T)(T x) @trusted @nogc nothrow pure if (__traits(isScalar, T)) { 181 static if (__traits(isFloating, T)) { 182 183 double tmp = cast(double)x; 184 return 0 > (*cast(long*)&tmp); 185 } else { 186 187 return 0 > x; 188 } 189 } 190 191 /** 192 Copies the sign-bit from one value to another. 193 194 Params: 195 to = The value to copy to 196 from = The value to copy from 197 198 Returns: 199 The value of $(D to) with the sign bit flipped 200 to match $(D from). 201 */ 202 T copysign(T)(T to, T from) @safe @nogc nothrow pure if (__traits(isScalar, T)) { 203 return signbit(to) == signbit(from) ? to : -to; 204 }