scalar_math.hh 10.3 KB
Newer Older
Philip Trettner's avatar
Philip Trettner committed
1
2
3
4
#pragma once

#include <cmath>

Philip Trettner's avatar
Philip Trettner committed
5
#include <typed-geometry/common/constants.hh>
Philip Trettner's avatar
Philip Trettner committed
6
#include <typed-geometry/detail/scalar_traits.hh>
Philip Trettner's avatar
Philip Trettner committed
7
#include <typed-geometry/detail/utility.hh>
Philip Trettner's avatar
Philip Trettner committed
8
#include <typed-geometry/types/angle.hh>
Philip Trettner's avatar
Philip Trettner committed
9
#include <typed-geometry/types/scalars/default.hh>
Philip Trettner's avatar
Philip Trettner committed
10
11

// TODO:
Philip Trettner's avatar
Philip Trettner committed
12
13
14
// - proper f8, f16
// - constexpr abs, sqrt, pow, ...
// - asm sqrt
Kersten Schuster's avatar
Kersten Schuster committed
15
// - floor, ceil, round
Philip Trettner's avatar
Philip Trettner committed
16
17
18

namespace tg
{
19
20
21
// ==================================================================
// Classification

22
23
[[nodiscard]] inline bool is_nan(f32 x) { return std::isnan(x); }
[[nodiscard]] inline bool is_nan(f64 x) { return std::isnan(x); }
24

25
26
27
28
[[nodiscard]] inline bool is_zero(i8 x) { return x == 0; }
[[nodiscard]] inline bool is_zero(i16 x) { return x == 0; }
[[nodiscard]] inline bool is_zero(i32 x) { return x == 0; }
[[nodiscard]] inline bool is_zero(i64 x) { return x == 0; }
Philip Trettner's avatar
Philip Trettner committed
29

30
31
32
33
[[nodiscard]] inline bool is_zero(u8 x) { return x == 0; }
[[nodiscard]] inline bool is_zero(u16 x) { return x == 0; }
[[nodiscard]] inline bool is_zero(u32 x) { return x == 0; }
[[nodiscard]] inline bool is_zero(u64 x) { return x == 0; }
Philip Trettner's avatar
Philip Trettner committed
34

35
36
[[nodiscard]] inline bool is_zero(f32 x) { return x == 0; }
[[nodiscard]] inline bool is_zero(f64 x) { return x == 0; }
37

38
39
[[nodiscard]] inline bool is_inf(f32 x) { return std::isinf(x); }
[[nodiscard]] inline bool is_inf(f64 x) { return std::isinf(x); }
40

41
42
[[nodiscard]] inline bool is_finite(f32 x) { return std::isfinite(x); }
[[nodiscard]] inline bool is_finite(f64 x) { return std::isfinite(x); }
43

44
45
[[nodiscard]] inline bool is_normal(f32 x) { return std::isnormal(x); }
[[nodiscard]] inline bool is_normal(f64 x) { return std::isnormal(x); }
46

47
48
[[nodiscard]] inline bool is_subnormal(f32 x) { return is_finite(x) && !is_normal(x); }
[[nodiscard]] inline bool is_subnormal(f64 x) { return is_finite(x) && !is_normal(x); }
49
50
51
52
53
54
55
56
57
58

enum class fp_class
{
    infinite,
    nan,
    normal,
    subnormal,
    zero
};

59
[[nodiscard]] inline fp_class fp_classify(f32 x)
60
61
62
63
64
65
66
67
68
69
70
71
{
    if (is_nan(x))
        return fp_class::nan;
    else if (is_inf(x))
        return fp_class::infinite;
    else if (is_zero(x))
        return fp_class::zero;
    else if (is_normal(x))
        return fp_class::normal;
    else
        return fp_class::subnormal;
}
72
[[nodiscard]] inline fp_class fp_classify(f64 x)
73
74
75
76
77
78
79
80
81
82
83
84
85
{
    if (is_nan(x))
        return fp_class::nan;
    else if (is_inf(x))
        return fp_class::infinite;
    else if (is_zero(x))
        return fp_class::zero;
    else if (is_normal(x))
        return fp_class::normal;
    else
        return fp_class::subnormal;
}

Philip Trettner's avatar
Philip Trettner committed
86
87
88
// ==================================================================
// Basics

89
// -i8 is i32 consistent with C++
90
91
92
93
[[nodiscard]] inline i32 abs(i8 v) { return v < 0 ? -v : v; }
[[nodiscard]] inline i32 abs(i16 v) { return v < 0 ? -v : v; }
[[nodiscard]] inline i32 abs(i32 v) { return v < 0 ? -v : v; }
[[nodiscard]] inline i64 abs(i64 v) { return v < 0 ? -v : v; }
Philip Trettner's avatar
Philip Trettner committed
94

95
96
97
98
[[nodiscard]] inline u8 abs(u8 v) { return v; }
[[nodiscard]] inline u16 abs(u16 v) { return v; }
[[nodiscard]] inline u32 abs(u32 v) { return v; }
[[nodiscard]] inline u64 abs(u64 v) { return v; }
Philip Trettner's avatar
Philip Trettner committed
99

100
101
102
103
[[nodiscard]] inline f8 abs(f8 v) { return v; }
[[nodiscard]] inline f16 abs(f16 v) { return std::abs(v); }
[[nodiscard]] inline f32 abs(f32 v) { return std::abs(v); }
[[nodiscard]] inline f64 abs(f64 v) { return std::abs(v); }
Philip Trettner's avatar
Philip Trettner committed
104

105
template <class T>
106
[[nodiscard]] angle_t<T> abs(angle_t<T> a)
107
{
Philip Trettner's avatar
Philip Trettner committed
108
    return radians(abs(a.radians()));
109
}
Philip Trettner's avatar
Philip Trettner committed
110

111
112
113
114
[[nodiscard]] inline f32 floor(f32 v) { return std::floor(v); }
[[nodiscard]] inline f64 floor(f64 v) { return std::floor(v); }
[[nodiscard]] inline i32 ifloor(f32 v) { return v >= 0 || f32(i32(v)) == v ? i32(v) : i32(v) - 1; }
[[nodiscard]] inline i64 ifloor(f64 v) { return v >= 0 || f64(i64(v)) == v ? i64(v) : i64(v) - 1; }
115

116
117
118
119
[[nodiscard]] inline f32 ceil(f32 v) { return std::ceil(v); }
[[nodiscard]] inline f64 ceil(f64 v) { return std::ceil(v); }
[[nodiscard]] inline i32 iceil(f32 v) { return v <= 0 || f32(i32(v)) == v ? i32(v) : i32(v) + 1; }
[[nodiscard]] inline i64 iceil(f64 v) { return v <= 0 || f64(i64(v)) == v ? i64(v) : i64(v) + 1; }
120

121
122
123
124
[[nodiscard]] inline f32 round(f32 v) { return std::round(v); }
[[nodiscard]] inline f64 round(f64 v) { return std::round(v); }
[[nodiscard]] inline i32 iround(f32 v) { return v >= 0 ? i32(v + 0.5f) : i32(v - 0.5f); }
[[nodiscard]] inline i64 iround(f64 v) { return v >= 0 ? i64(v + 0.5) : i64(v - 0.5); }
125

126
127
[[nodiscard]] inline f32 fract(f32 v) { return v - floor(v); }
[[nodiscard]] inline f64 fract(f64 v) { return v - floor(v); }
Philip Trettner's avatar
Philip Trettner committed
128

Philip Trettner's avatar
Philip Trettner committed
129
template <class A, class B>
130
[[nodiscard]] constexpr auto min(A&& a, B&& b) -> decltype(a < b ? a : b)
Kersten Schuster's avatar
Kersten Schuster committed
131
{
Philip Trettner's avatar
Philip Trettner committed
132
    return a < b ? a : b;
Kersten Schuster's avatar
Kersten Schuster committed
133
}
Philip Trettner's avatar
Philip Trettner committed
134
template <class A, class B>
135
[[nodiscard]] constexpr auto max(A&& a, B&& b) -> decltype(a < b ? b : a)
Kersten Schuster's avatar
Kersten Schuster committed
136
{
Philip Trettner's avatar
Philip Trettner committed
137
    return a < b ? b : a;
Kersten Schuster's avatar
Kersten Schuster committed
138
}
139

Philip Trettner's avatar
Philip Trettner committed
140
template <class A, class B, class C>
141
[[nodiscard]] constexpr auto clamp(A const& a, B const& min_value, C const& max_value) -> decltype(min(max(a, min_value), max_value))
Philip Trettner's avatar
Philip Trettner committed
142
143
144
145
{
    return min(max(a, min_value), max_value);
}

Philip Trettner's avatar
Philip Trettner committed
146
template <class T>
147
[[nodiscard]] constexpr T saturate(T const& a)
148
149
150
151
{
    return clamp(a, T(0), T(1));
}

Philip Trettner's avatar
Philip Trettner committed
152
template <class T, class = enable_if<is_scalar<T>>>
153
[[nodiscard]] constexpr T sign(T const& v)
Philip Trettner's avatar
Philip Trettner committed
154
155
156
157
158
159
160
161
162
163
{
    // TODO: optimize?
    // if constexpr (!is_unsigned_integer<T>)
    if (v < T(0))
        return T(-1);
    if (v > T(0))
        return T(1);
    return T(0);
}

Philip Trettner's avatar
Philip Trettner committed
164
165
166
// ==================================================================
// Powers

167
168
169
170
[[nodiscard]] inline f32 pow(f32 b, f32 e) { return std::pow(b, e); }
[[nodiscard]] inline f32 pow(f32 b, i32 e) { return f32(std::pow(b, e)); }
[[nodiscard]] inline f64 pow(f64 b, f64 e) { return std::pow(b, e); }
[[nodiscard]] inline f64 pow(f64 b, i32 e) { return std::pow(b, e); }
Philip Trettner's avatar
Philip Trettner committed
171

Philip Trettner's avatar
Philip Trettner committed
172
template <class T>
173
[[nodiscard]] constexpr auto pow2(T const& v) -> decltype(v * v)
174
175
176
{
    return v * v;
}
Philip Trettner's avatar
Philip Trettner committed
177
template <class T>
178
[[nodiscard]] constexpr auto pow3(T const& v) -> decltype(v * v * v)
179
180
181
{
    return v * v * v;
}
Philip Trettner's avatar
Philip Trettner committed
182
template <class T>
183
[[nodiscard]] constexpr auto pow4(T const& v) -> decltype((v * v) * (v * v))
184
{
Philip Trettner's avatar
Philip Trettner committed
185
    auto const v2 = v * v;
Kersten Schuster's avatar
Kersten Schuster committed
186
    return v2 * v2;
187
}
Philip Trettner's avatar
Philip Trettner committed
188
template <class T>
189
[[nodiscard]] constexpr auto pow5(T const& v) -> decltype((v * v) * (v * v) * v)
190
{
Philip Trettner's avatar
Philip Trettner committed
191
    auto const v2 = v * v;
Kersten Schuster's avatar
Kersten Schuster committed
192
    return v2 * v2 * v;
193
}
Philip Trettner's avatar
Philip Trettner committed
194

195
196
[[nodiscard]] inline f32 sqrt(f32 v) { return std::sqrt(v); }
[[nodiscard]] inline f64 sqrt(f64 v) { return std::sqrt(v); }
Philip Trettner's avatar
Philip Trettner committed
197

198
// TODO: _mm_rsqrt_ss
199
200
[[nodiscard]] inline f32 rsqrt(f32 v) { return 1 / std::sqrt(v); }
[[nodiscard]] inline f64 rsqrt(f64 v) { return 1 / std::sqrt(v); }
201

202
203
[[nodiscard]] inline f32 cbrt(f32 v) { return std::cbrt(v); }
[[nodiscard]] inline f64 cbrt(f64 v) { return std::cbrt(v); }
Philip Trettner's avatar
Philip Trettner committed
204
205
206
207

// ==================================================================
// Exponentials

208
209
[[nodiscard]] inline f32 exp(f32 v) { return std::exp(v); }
[[nodiscard]] inline f64 exp(f64 v) { return std::exp(v); }
Philip Trettner's avatar
Philip Trettner committed
210

211
212
[[nodiscard]] inline f32 exp2(f32 v) { return std::exp2(v); }
[[nodiscard]] inline f64 exp2(f64 v) { return std::exp2(v); }
Philip Trettner's avatar
Philip Trettner committed
213

214
215
[[nodiscard]] inline f32 log(f32 v) { return std::log(v); }
[[nodiscard]] inline f64 log(f64 v) { return std::log(v); }
Philip Trettner's avatar
Philip Trettner committed
216

217
218
[[nodiscard]] inline f32 log2(f32 v) { return std::log2(v); }
[[nodiscard]] inline f64 log2(f64 v) { return std::log2(v); }
Philip Trettner's avatar
Philip Trettner committed
219

220
221
[[nodiscard]] inline f32 log10(f32 v) { return std::log10(v); }
[[nodiscard]] inline f64 log10(f64 v) { return std::log10(v); }
Philip Trettner's avatar
Philip Trettner committed
222
223
224
225

// ==================================================================
// Trigonometry

226
227
228
229
230
231
[[nodiscard]] inline f32 sin(angle_t<f32> v) { return std::sin(v.radians()); }
[[nodiscard]] inline f64 sin(angle_t<f64> v) { return std::sin(v.radians()); }
[[nodiscard]] inline f32 cos(angle_t<f32> v) { return std::cos(v.radians()); }
[[nodiscard]] inline f64 cos(angle_t<f64> v) { return std::cos(v.radians()); }
[[nodiscard]] inline f32 tan(angle_t<f32> v) { return std::tan(v.radians()); }
[[nodiscard]] inline f64 tan(angle_t<f64> v) { return std::tan(v.radians()); }
Philip Trettner's avatar
Philip Trettner committed
232

233
// TODO: use SSE intrinsic to compute both directly
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
[[nodiscard]] inline pair<f32, f32> sin_cos(angle_t<f32> v) { return {sin(v), cos(v)}; }
[[nodiscard]] inline pair<f64, f64> sin_cos(angle_t<f64> v) { return {sin(v), cos(v)}; }

[[nodiscard]] inline angle_t<f32> asin(f32 v) { return radians(std::asin(v)); }
[[nodiscard]] inline angle_t<f64> asin(f64 v) { return radians(std::asin(v)); }
[[nodiscard]] inline angle_t<f32> acos(f32 v) { return radians(std::acos(v)); }
[[nodiscard]] inline angle_t<f64> acos(f64 v) { return radians(std::acos(v)); }
[[nodiscard]] inline angle_t<f32> atan(f32 v) { return radians(std::atan(v)); }
[[nodiscard]] inline angle_t<f64> atan(f64 v) { return radians(std::atan(v)); }

[[nodiscard]] inline angle_t<f32> atan2(f32 y, f32 x) { return radians(std::atan2(y, x)); }
[[nodiscard]] inline angle_t<f64> atan2(f64 y, f64 x) { return radians(std::atan2(y, x)); }

[[nodiscard]] inline f32 sinh(f32 v) { return std::sinh(v); }
[[nodiscard]] inline f64 sinh(f64 v) { return std::sinh(v); }
[[nodiscard]] inline f32 cosh(f32 v) { return std::cosh(v); }
[[nodiscard]] inline f64 cosh(f64 v) { return std::cosh(v); }
[[nodiscard]] inline f32 tanh(f32 v) { return std::tanh(v); }
[[nodiscard]] inline f64 tanh(f64 v) { return std::tanh(v); }

[[nodiscard]] inline f32 asinh(f32 v) { return std::asinh(v); }
[[nodiscard]] inline f64 asinh(f64 v) { return std::asinh(v); }
[[nodiscard]] inline f32 acosh(f32 v) { return std::acosh(v); }
[[nodiscard]] inline f64 acosh(f64 v) { return std::acosh(v); }
[[nodiscard]] inline f32 atanh(f32 v) { return std::atanh(v); }
[[nodiscard]] inline f64 atanh(f64 v) { return std::atanh(v); }
Philip Trettner's avatar
Philip Trettner committed
260

Philip Trettner's avatar
Philip Trettner committed
261
262
263
// ==================================================================
// other GLSL

264
[[nodiscard]] constexpr f32 smoothstep(f32 edge0, f32 edge1, f32 x)
Philip Trettner's avatar
Philip Trettner committed
265
266
267
268
{
    auto t = clamp((x - edge0) / (edge1 - edge0), f32(0), f32(1));
    return t * t * (f32(3) - f32(2) * t);
}
269
[[nodiscard]] constexpr f64 smoothstep(f64 edge0, f64 edge1, f64 x)
Philip Trettner's avatar
Philip Trettner committed
270
271
272
273
274
{
    auto t = clamp((x - edge0) / (edge1 - edge0), f64(0), f64(1));
    return t * t * (f64(3) - f64(2) * t);
}

275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
// ==================================================================
// other Math

// greatest common divisor
// TODO: how does it work for negative numbers?
template <class T, class = enable_if<is_integer<T>>>
T gcd(T a, T b)
{
    return b ? gcd(b, a % b) : a;
}

// least common multiple
// CAREFUL: a * b is an intermediate and might overflow
template <class T, class = enable_if<is_integer<T>>>
T lcm(T a, T b)
{
    return a * b / gcd(a, b);
}

Philip Trettner's avatar
Philip Trettner committed
294
} // namespace tg