/src/serenity/AK/Math/Trigonometry.h
Line | Count | Source |
1 | | /* |
2 | | * Copyright (c) 2021-2026, Leon Albrecht <leon2002.la@gmail.com>. |
3 | | * Copyright (c) 2024-2026, Nico Weber <thakis@chromium.org>. |
4 | | * |
5 | | * SPDX-License-Identifier: BSD-2-Clause |
6 | | */ |
7 | | |
8 | | #pragma once |
9 | | |
10 | | #include <AK/Concepts.h> |
11 | | #include <AK/Math/Constants.h> |
12 | | #include <AK/Math/Copysign.h> |
13 | | #include <AK/Math/Sqrt.h> |
14 | | #include <math.h> |
15 | | |
16 | | #include <AK/Math/Macros.h> |
17 | | |
18 | | namespace AK { |
19 | | |
20 | | namespace Trigonometry { |
21 | | |
22 | | template<FloatingPoint T> |
23 | | constexpr T hypot(T x, T y) |
24 | 16.0M | { |
25 | | // "If x or y is ±Inf, +Inf shall be returned (even if one of x or y is NaN)." |
26 | 16.0M | if (isinf(x) || isinf(y)) |
27 | 0 | return Infinity<T>; |
28 | | // "If x or y is NaN, and the other is not ±Inf, a NaN shall be returned." |
29 | 16.0M | if (isnan(x) || isnan(y)) |
30 | 0 | return NaN<T>; |
31 | | |
32 | 16.0M | return sqrt(x * x + y * y); |
33 | 16.0M | } |
34 | | |
35 | | template<FloatingPoint T> |
36 | | constexpr T sin(T angle) |
37 | 10.7M | { |
38 | 10.7M | CONSTEXPR_STATE(sin, angle); |
39 | | |
40 | 10.7M | #if ARCH(X86_64) |
41 | 10.7M | T ret; |
42 | 10.7M | asm( |
43 | 10.7M | "fsin" |
44 | 10.7M | : "=t"(ret) |
45 | 10.7M | : "0"(angle)); |
46 | 10.7M | return ret; |
47 | | #else |
48 | | # if defined(AK_OS_SERENITY) |
49 | | T sign = 1; |
50 | | if (angle < 0) { |
51 | | angle = -angle; |
52 | | sign = -1; |
53 | | } |
54 | | |
55 | | if (angle >= 2 * Pi<T>) |
56 | | angle = fmod(angle, 2 * Pi<T>); |
57 | | |
58 | | if (angle >= Pi<T>) { |
59 | | angle = angle - Pi<T>; |
60 | | sign = -sign; |
61 | | } |
62 | | |
63 | | if (angle > Pi<T> / 2) |
64 | | angle = Pi<T> - angle; |
65 | | |
66 | | // https://github.com/samhocevar/lolremez/wiki/Tutorial-4-of-5%3A-fixing-lower-order-parameters |
67 | | auto f = [](T x) { |
68 | | if constexpr (IsSame<T, f32>) { |
69 | | // lolremez --float --degree 3 --range "1e-50:pi*pi/4" "(sin(sqrt(x))-sqrt(x))/(x*sqrt(x))" "1/(x*sqrt(x))" |
70 | | // "Estimated max error: 4.618689e-9" |
71 | | f32 u = 2.6000548e-06f; |
72 | | u = u * x + -0.00019806615f; |
73 | | u = u * x + 0.0083330171f; |
74 | | return u * x + -0.16666657f; |
75 | | } else { |
76 | | // FIXME: Could do something custom for long double. |
77 | | // lolremez --degree 6 --range "1e-50:pi*pi/4" "(sin(sqrt(x))-sqrt(x))/(x*sqrt(x))" "1/(x*sqrt(x))" |
78 | | // "Estimated max error: 1.1015766629825144e-16" |
79 | | T u = -7.3646464502210486e-13; |
80 | | u = u * x + 1.6047301196685753e-10; |
81 | | u = u * x + -2.5051851497012596e-08; |
82 | | u = u * x + 2.7557316077007725e-06; |
83 | | u = u * x + -0.00019841269820094207; |
84 | | u = u * x + 0.0083333333332628792; |
85 | | return u * x + -0.16666666666665811; |
86 | | } |
87 | | }; |
88 | | T angle_squared = angle * angle; |
89 | | return sign * (angle + angle * angle_squared * f(angle_squared)); |
90 | | # else |
91 | | CALL_BUILTIN(sin, angle); |
92 | | # endif |
93 | | #endif |
94 | 10.7M | } _ZN2AK12Trigonometry3sinITkNS_8Concepts13FloatingPointEfEET_S3_ Line | Count | Source | 37 | 10.7M | { | 38 | 10.7M | CONSTEXPR_STATE(sin, angle); | 39 | | | 40 | 10.7M | #if ARCH(X86_64) | 41 | 10.7M | T ret; | 42 | 10.7M | asm( | 43 | 10.7M | "fsin" | 44 | 10.7M | : "=t"(ret) | 45 | 10.7M | : "0"(angle)); | 46 | 10.7M | return ret; | 47 | | #else | 48 | | # if defined(AK_OS_SERENITY) | 49 | | T sign = 1; | 50 | | if (angle < 0) { | 51 | | angle = -angle; | 52 | | sign = -1; | 53 | | } | 54 | | | 55 | | if (angle >= 2 * Pi<T>) | 56 | | angle = fmod(angle, 2 * Pi<T>); | 57 | | | 58 | | if (angle >= Pi<T>) { | 59 | | angle = angle - Pi<T>; | 60 | | sign = -sign; | 61 | | } | 62 | | | 63 | | if (angle > Pi<T> / 2) | 64 | | angle = Pi<T> - angle; | 65 | | | 66 | | // https://github.com/samhocevar/lolremez/wiki/Tutorial-4-of-5%3A-fixing-lower-order-parameters | 67 | | auto f = [](T x) { | 68 | | if constexpr (IsSame<T, f32>) { | 69 | | // lolremez --float --degree 3 --range "1e-50:pi*pi/4" "(sin(sqrt(x))-sqrt(x))/(x*sqrt(x))" "1/(x*sqrt(x))" | 70 | | // "Estimated max error: 4.618689e-9" | 71 | | f32 u = 2.6000548e-06f; | 72 | | u = u * x + -0.00019806615f; | 73 | | u = u * x + 0.0083330171f; | 74 | | return u * x + -0.16666657f; | 75 | | } else { | 76 | | // FIXME: Could do something custom for long double. | 77 | | // lolremez --degree 6 --range "1e-50:pi*pi/4" "(sin(sqrt(x))-sqrt(x))/(x*sqrt(x))" "1/(x*sqrt(x))" | 78 | | // "Estimated max error: 1.1015766629825144e-16" | 79 | | T u = -7.3646464502210486e-13; | 80 | | u = u * x + 1.6047301196685753e-10; | 81 | | u = u * x + -2.5051851497012596e-08; | 82 | | u = u * x + 2.7557316077007725e-06; | 83 | | u = u * x + -0.00019841269820094207; | 84 | | u = u * x + 0.0083333333332628792; | 85 | | return u * x + -0.16666666666665811; | 86 | | } | 87 | | }; | 88 | | T angle_squared = angle * angle; | 89 | | return sign * (angle + angle * angle_squared * f(angle_squared)); | 90 | | # else | 91 | | CALL_BUILTIN(sin, angle); | 92 | | # endif | 93 | | #endif | 94 | 10.7M | } |
Unexecuted instantiation: _ZN2AK12Trigonometry3sinITkNS_8Concepts13FloatingPointEdEET_S3_ |
95 | | |
96 | | template<FloatingPoint T> |
97 | | constexpr T cos(T angle) |
98 | 5.77k | { |
99 | 5.77k | CONSTEXPR_STATE(cos, angle); |
100 | | |
101 | 5.77k | #if ARCH(X86_64) |
102 | 5.77k | T ret; |
103 | 5.77k | asm( |
104 | 5.77k | "fcos" |
105 | 5.77k | : "=t"(ret) |
106 | 5.77k | : "0"(angle)); |
107 | 5.77k | return ret; |
108 | | #else |
109 | | # if defined(AK_OS_SERENITY) |
110 | | if (angle < 0) |
111 | | angle = -angle; |
112 | | |
113 | | if (angle >= 2 * Pi<T>) |
114 | | angle = fmod(angle, 2 * Pi<T>); |
115 | | |
116 | | T sign = 1; |
117 | | if (angle >= Pi<T>) { |
118 | | angle = angle - Pi<T>; |
119 | | sign = -1; |
120 | | } |
121 | | |
122 | | if (angle > Pi<T> / 2) { |
123 | | angle = Pi<T> - angle; |
124 | | sign = -sign; |
125 | | } |
126 | | |
127 | | // https://github.com/samhocevar/lolremez/wiki/Tutorial-3-of-5:-changing-variables-for-simpler-polynomials |
128 | | // for cos(x): cos(x) - 1 is a function of x^2 terms, so we do the substitution on that page like |
129 | | // max | cos(x) - 1 - x^2 Q(x^2) | and then set y = x^2. That yields: |
130 | | // max | (cos(sqrt(y)) - 1) - y Q(y) |. Dividing through with y (instead of sqrt(y) as in the sin() case) gives us: |
131 | | auto f = [](T x) { |
132 | | if constexpr (IsSame<T, f32>) { |
133 | | // lolremez --float --degree 3 --range "1e-50:pi*pi/4" "(cos(sqrt(x)) - 1)/x" "1/x" |
134 | | // "Estimated max error: 5.2720126e-8" |
135 | | f32 u = 2.3194387e-05f; |
136 | | u = u * x + -0.0013855927f; |
137 | | u = u * x + 0.041663989f; |
138 | | return u * x + -0.49999931f; |
139 | | } else { |
140 | | // lolremez --degree 6 --range "1e-50:pi*pi/4" "(cos(sqrt(x)) - 1)/x" "1/x" |
141 | | // "Estimated max error: 2.0880269759116624e-15" |
142 | | T u = -1.101249182846601e-11; |
143 | | u = u * x + 2.0858735176345955e-09; |
144 | | u = u * x + -2.7556950755056579e-07; |
145 | | u = u * x + 2.4801583156341611e-05; |
146 | | u = u * x + -0.001388888886393419; |
147 | | u = u * x + 0.041666666665954213; |
148 | | return u * x + -0.49999999999993117; |
149 | | } |
150 | | }; |
151 | | T angle_squared = angle * angle; |
152 | | return sign * (1 + angle_squared * f(angle_squared)); |
153 | | # else |
154 | | CALL_BUILTIN(cos, angle); |
155 | | # endif |
156 | | #endif |
157 | 5.77k | } _ZN2AK12Trigonometry3cosITkNS_8Concepts13FloatingPointEfEET_S3_ Line | Count | Source | 98 | 5.77k | { | 99 | 5.77k | CONSTEXPR_STATE(cos, angle); | 100 | | | 101 | 5.77k | #if ARCH(X86_64) | 102 | 5.77k | T ret; | 103 | 5.77k | asm( | 104 | 5.77k | "fcos" | 105 | 5.77k | : "=t"(ret) | 106 | 5.77k | : "0"(angle)); | 107 | 5.77k | return ret; | 108 | | #else | 109 | | # if defined(AK_OS_SERENITY) | 110 | | if (angle < 0) | 111 | | angle = -angle; | 112 | | | 113 | | if (angle >= 2 * Pi<T>) | 114 | | angle = fmod(angle, 2 * Pi<T>); | 115 | | | 116 | | T sign = 1; | 117 | | if (angle >= Pi<T>) { | 118 | | angle = angle - Pi<T>; | 119 | | sign = -1; | 120 | | } | 121 | | | 122 | | if (angle > Pi<T> / 2) { | 123 | | angle = Pi<T> - angle; | 124 | | sign = -sign; | 125 | | } | 126 | | | 127 | | // https://github.com/samhocevar/lolremez/wiki/Tutorial-3-of-5:-changing-variables-for-simpler-polynomials | 128 | | // for cos(x): cos(x) - 1 is a function of x^2 terms, so we do the substitution on that page like | 129 | | // max | cos(x) - 1 - x^2 Q(x^2) | and then set y = x^2. That yields: | 130 | | // max | (cos(sqrt(y)) - 1) - y Q(y) |. Dividing through with y (instead of sqrt(y) as in the sin() case) gives us: | 131 | | auto f = [](T x) { | 132 | | if constexpr (IsSame<T, f32>) { | 133 | | // lolremez --float --degree 3 --range "1e-50:pi*pi/4" "(cos(sqrt(x)) - 1)/x" "1/x" | 134 | | // "Estimated max error: 5.2720126e-8" | 135 | | f32 u = 2.3194387e-05f; | 136 | | u = u * x + -0.0013855927f; | 137 | | u = u * x + 0.041663989f; | 138 | | return u * x + -0.49999931f; | 139 | | } else { | 140 | | // lolremez --degree 6 --range "1e-50:pi*pi/4" "(cos(sqrt(x)) - 1)/x" "1/x" | 141 | | // "Estimated max error: 2.0880269759116624e-15" | 142 | | T u = -1.101249182846601e-11; | 143 | | u = u * x + 2.0858735176345955e-09; | 144 | | u = u * x + -2.7556950755056579e-07; | 145 | | u = u * x + 2.4801583156341611e-05; | 146 | | u = u * x + -0.001388888886393419; | 147 | | u = u * x + 0.041666666665954213; | 148 | | return u * x + -0.49999999999993117; | 149 | | } | 150 | | }; | 151 | | T angle_squared = angle * angle; | 152 | | return sign * (1 + angle_squared * f(angle_squared)); | 153 | | # else | 154 | | CALL_BUILTIN(cos, angle); | 155 | | # endif | 156 | | #endif | 157 | 5.77k | } |
Unexecuted instantiation: _ZN2AK12Trigonometry3cosITkNS_8Concepts13FloatingPointEdEET_S3_ |
158 | | |
159 | | template<FloatingPoint T> |
160 | | constexpr void sincos(T angle, T& sin_val, T& cos_val) |
161 | 25.4M | { |
162 | 25.4M | if (is_constant_evaluated()) { |
163 | 0 | sin_val = sin(angle); |
164 | 0 | cos_val = cos(angle); |
165 | 0 | return; |
166 | 0 | } |
167 | 25.4M | #if ARCH(X86_64) |
168 | 25.4M | asm( |
169 | 25.4M | "fsincos" |
170 | 25.4M | : "=t"(cos_val), "=u"(sin_val) |
171 | 25.4M | : "0"(angle)); |
172 | | #else |
173 | | sin_val = sin(angle); |
174 | | cos_val = cos(angle); |
175 | | #endif |
176 | 25.4M | } _ZN2AK12Trigonometry6sincosITkNS_8Concepts13FloatingPointEfEEvT_RS3_S4_ Line | Count | Source | 161 | 24.0M | { | 162 | 24.0M | if (is_constant_evaluated()) { | 163 | 0 | sin_val = sin(angle); | 164 | 0 | cos_val = cos(angle); | 165 | 0 | return; | 166 | 0 | } | 167 | 24.0M | #if ARCH(X86_64) | 168 | 24.0M | asm( | 169 | 24.0M | "fsincos" | 170 | 24.0M | : "=t"(cos_val), "=u"(sin_val) | 171 | 24.0M | : "0"(angle)); | 172 | | #else | 173 | | sin_val = sin(angle); | 174 | | cos_val = cos(angle); | 175 | | #endif | 176 | 24.0M | } |
_ZN2AK12Trigonometry6sincosITkNS_8Concepts13FloatingPointEdEEvT_RS3_S4_ Line | Count | Source | 161 | 1.38M | { | 162 | 1.38M | if (is_constant_evaluated()) { | 163 | 0 | sin_val = sin(angle); | 164 | 0 | cos_val = cos(angle); | 165 | 0 | return; | 166 | 0 | } | 167 | 1.38M | #if ARCH(X86_64) | 168 | 1.38M | asm( | 169 | 1.38M | "fsincos" | 170 | 1.38M | : "=t"(cos_val), "=u"(sin_val) | 171 | 1.38M | : "0"(angle)); | 172 | | #else | 173 | | sin_val = sin(angle); | 174 | | cos_val = cos(angle); | 175 | | #endif | 176 | 1.38M | } |
|
177 | | |
178 | | template<FloatingPoint T> |
179 | | constexpr T tan(T angle) |
180 | 10.9M | { |
181 | 10.9M | CONSTEXPR_STATE(tan, angle); |
182 | | |
183 | 10.9M | #if ARCH(X86_64) |
184 | 10.9M | T ret, one; |
185 | 10.9M | asm( |
186 | 10.9M | "fptan" |
187 | 10.9M | : "=t"(one), "=u"(ret) |
188 | 10.9M | : "0"(angle)); |
189 | | |
190 | 10.9M | return ret; |
191 | | #else |
192 | | # if defined(AK_OS_SERENITY) |
193 | | return sin(angle) / cos(angle); |
194 | | # else |
195 | | CALL_BUILTIN(tan, angle); |
196 | | # endif |
197 | | #endif |
198 | 10.9M | } |
199 | | |
200 | | template<FloatingPoint T> |
201 | | constexpr T asin(T x) |
202 | | { |
203 | | CONSTEXPR_STATE(asin, x); |
204 | | |
205 | | T sign = 1; |
206 | | |
207 | | // asin(x) = -asin(-x) |
208 | | // Use this to make x always >= 0. |
209 | | if (x < 0) { |
210 | | x = -x; |
211 | | sign = -1; |
212 | | } |
213 | | |
214 | | if (x > 1) |
215 | | return NaN<T>; |
216 | | |
217 | | T offset = 0; |
218 | | if (x > (T)0.5) { |
219 | | // asin(x) = pi/2 - 2 * asin(sqrt((1 - x) / 2)) |
220 | | // If 0.5 < x <= 1, then sqrt((1 - x) / 2 ) < 0.5, |
221 | | // and the recursion will go down the `<= 0.5` branch. |
222 | | offset = sign * Pi<T> / 2; |
223 | | x = sqrt((1 - x) / 2); |
224 | | sign = sign * -2; |
225 | | } |
226 | | |
227 | | // https://github.com/samhocevar/lolremez/wiki/Tutorial-4-of-5%3A-fixing-lower-order-parameters |
228 | | auto f = [](T x) { |
229 | | if constexpr (IsSame<T, f32>) { |
230 | | // lolremez --float --degree 5 --range "1e-50:1/4" "(asin(sqrt(x))-sqrt(x))/(x*sqrt(x))" "1/(x*sqrt(x))" |
231 | | // "Estimated max error: 8.5377476e-11" |
232 | | f32 u = 0.039088517f; |
233 | | u = u * x + 0.013210787f; |
234 | | u = u * x + 0.032170232f; |
235 | | u = u * x + 0.044467181f; |
236 | | u = u * x + 0.075008146f; |
237 | | return u * x + 0.16666654f; |
238 | | } else { |
239 | | // FIXME: Could do something custom for long double. |
240 | | // lolremez --degree 11 --range "1e-50:1/4" "(asin(sqrt(x))-sqrt(x))/(x*sqrt(x))" "1/(x*sqrt(x))" |
241 | | // "Estimated max error: 4.3580651842210834e-18" |
242 | | T u = 0.033660618737294104; |
243 | | u = u * x + -0.018930050045116921; |
244 | | u = u * x + 0.021367266503207524; |
245 | | u = u * x + 0.005806729168702944; |
246 | | u = u * x + 0.012349571686866897; |
247 | | u = u * x + 0.013855247852028872; |
248 | | u = u * x + 0.017363058003675668; |
249 | | u = u * x + 0.022371508713932262; |
250 | | u = u * x + 0.03038197103721944; |
251 | | u = u * x + 0.044642856488345831; |
252 | | u = u * x + 0.075000000008444367; |
253 | | return u * x + 0.16666666666662552; |
254 | | } |
255 | | }; |
256 | | T x2 = x * x; |
257 | | return offset + sign * (x + x * x2 * f(x2)); |
258 | | } |
259 | | |
260 | | template<FloatingPoint T> |
261 | | constexpr T acos(T value) |
262 | | { |
263 | | CONSTEXPR_STATE(acos, value); |
264 | | |
265 | | // FIXME: I am naive |
266 | | return static_cast<T>(0.5) * Pi<T> - asin<T>(value); |
267 | | } |
268 | | |
269 | | template<FloatingPoint T> |
270 | | constexpr T atan(T value) |
271 | | { |
272 | | CONSTEXPR_STATE(atan, value); |
273 | | |
274 | | #if ARCH(X86_64) |
275 | | T ret; |
276 | | asm( |
277 | | "fld1\n" |
278 | | "fpatan\n" |
279 | | : "=t"(ret) |
280 | | : "0"(value)); |
281 | | return ret; |
282 | | #else |
283 | | # if defined(AK_OS_SERENITY) |
284 | | if (value < 0) |
285 | | return -atan(-value); |
286 | | |
287 | | if (value > 1) |
288 | | return Pi<T> / 2 - atan(1 / value); |
289 | | |
290 | | // atan is an odd function a leading factor of 1, so this uses the same substitutions as sin(). |
291 | | // See there for a description. |
292 | | auto f = [](T x) { |
293 | | if constexpr (IsSame<T, f32>) { |
294 | | // lolremez --float --degree 7 --range "1e-50:1" "(atan(sqrt(x))-sqrt(x))/(x*sqrt(x))" "1/(x*sqrt(x))" |
295 | | // "Estimated max error: 7.351572e-9" |
296 | | float u = 0.0026222446f; |
297 | | u = u * x + -0.015132537f; |
298 | | u = u * x + 0.041121863f; |
299 | | u = u * x + -0.073667064f; |
300 | | u = u * x + 0.10573932f; |
301 | | u = u * x + -0.14185975f; |
302 | | u = u * x + 0.19990396f; |
303 | | return u * x + -0.33332986f; |
304 | | } else { |
305 | | // FIXME: Could do something custom for long double. |
306 | | // lolremez --degree 18 --range "1e-50:1" "(atan(sqrt(x))-sqrt(x))/(x*sqrt(x))" "1/(x*sqrt(x))" |
307 | | // "Estimated max error: 1.1404872427303659e-17" |
308 | | T u = -1.7655716820720429e-05; |
309 | | u = u * x + 0.00019779001389463795; |
310 | | u = u * x + -0.001051241688667142; |
311 | | u = u * x + 0.0035473387064041432; |
312 | | u = u * x + -0.008607083731273335; |
313 | | u = u * x + 0.016195140256854403; |
314 | | u = u * x + -0.025032045474239515; |
315 | | u = u * x + 0.033452143515758453; |
316 | | u = u * x + -0.040558834967395346; |
317 | | u = u * x + 0.046569903058368577; |
318 | | u = u * x + -0.052332028769674181; |
319 | | u = u * x + 0.05875670126323182; |
320 | | u = u * x + -0.066655258619483626; |
321 | | u = u * x + 0.076921627408511856; |
322 | | u = u * x + -0.090908958935242423; |
323 | | u = u * x + 0.11111110296047642; |
324 | | u = u * x + -0.14285714254349127; |
325 | | u = u * x + 0.1999999999935175; |
326 | | return u * x + -0.33333333333328197; |
327 | | } |
328 | | }; |
329 | | T squared = value * value; |
330 | | return value + value * squared * f(squared); |
331 | | # endif |
332 | | CALL_BUILTIN(atan, value); |
333 | | #endif |
334 | | } |
335 | | |
336 | | template<FloatingPoint T> |
337 | | constexpr T atan2(T y, T x) |
338 | 2.76M | { |
339 | 2.76M | CONSTEXPR_STATE(atan2, y, x); |
340 | | |
341 | 2.76M | #if ARCH(X86_64) |
342 | 2.76M | T ret; |
343 | 2.76M | asm("fpatan" |
344 | 2.76M | : "=t"(ret) |
345 | 2.76M | : "0"(x), "u"(y) |
346 | 2.76M | : "st(1)"); |
347 | 2.76M | return ret; |
348 | | #else |
349 | | # if defined(AK_OS_SERENITY) |
350 | | if (__builtin_isnan(y)) |
351 | | return y; |
352 | | if (__builtin_isnan(x)) |
353 | | return x; |
354 | | |
355 | | // SPECIAL VALUES |
356 | | // atan2(±0, -0) returns ±pi. |
357 | | if (y == 0 && x == 0 && signbit(x)) |
358 | | return copysign(Pi<T>, y); |
359 | | |
360 | | // atan2(±0, +0) returns ±0. |
361 | | if (y == 0 && x == 0 && !signbit(x)) |
362 | | return y; |
363 | | |
364 | | // atan2(±0, x) returns ±pi for x < 0. |
365 | | if (y == 0 && x < 0) |
366 | | return copysign(Pi<T>, y); |
367 | | |
368 | | // atan2(±0, x) returns ±0 for x > 0. |
369 | | if (y == 0 && x > 0) |
370 | | return y; |
371 | | |
372 | | // atan2(y, ±0) returns +pi/2 for y > 0. |
373 | | if (y > 0 && x == 0) |
374 | | return Pi<T> / 2; |
375 | | |
376 | | // atan2(y, ±0) returns -pi/2 for y < 0. |
377 | | if (y < 0 && x == 0) |
378 | | return -Pi<T> / 2; |
379 | | |
380 | | // atan2(±y, -infinity) returns ±pi for finite y > 0. |
381 | | if (!__builtin_isinf(y) && y > 0 && __builtin_isinf(x) && signbit(x)) |
382 | | return copysign(Pi<T>, y); |
383 | | |
384 | | // atan2(±y, +infinity) returns ±0 for finite y > 0. |
385 | | if (!__builtin_isinf(y) && y > 0 && __builtin_isinf(x) && !signbit(x)) |
386 | | return copysign(static_cast<T>(0), y); |
387 | | |
388 | | // atan2(±infinity, x) returns ±pi/2 for finite x. |
389 | | if (__builtin_isinf(y) && !__builtin_isinf(x)) |
390 | | return copysign(Pi<T> / 2, y); |
391 | | |
392 | | // atan2(±infinity, -infinity) returns ±3*pi/4. |
393 | | if (__builtin_isinf(y) && __builtin_isinf(x) && signbit(x)) |
394 | | return copysign(3 * Pi<T> / 4, y); |
395 | | |
396 | | // atan2(±infinity, +infinity) returns ±pi/4. |
397 | | if (__builtin_isinf(y) && __builtin_isinf(x) && !signbit(x)) |
398 | | return copysign(Pi<T> / 4, y); |
399 | | |
400 | | // Check quadrant, going counterclockwise. |
401 | | if (y > 0 && x > 0) |
402 | | return atan(y / x); |
403 | | if (y > 0 && x < 0) |
404 | | return atan(y / x) + Pi<T>; |
405 | | if (y < 0 && x < 0) |
406 | | return atan(y / x) - Pi<T>; |
407 | | // y < 0 && x > 0 |
408 | | return atan(y / x); |
409 | | # else |
410 | | CALL_BUILTIN(atan2, y, x); |
411 | | # endif |
412 | | #endif |
413 | 2.76M | } _ZN2AK12Trigonometry5atan2ITkNS_8Concepts13FloatingPointEfEET_S3_S3_ Line | Count | Source | 338 | 3.05k | { | 339 | 3.05k | CONSTEXPR_STATE(atan2, y, x); | 340 | | | 341 | 3.05k | #if ARCH(X86_64) | 342 | 3.05k | T ret; | 343 | 3.05k | asm("fpatan" | 344 | 3.05k | : "=t"(ret) | 345 | 3.05k | : "0"(x), "u"(y) | 346 | 3.05k | : "st(1)"); | 347 | 3.05k | return ret; | 348 | | #else | 349 | | # if defined(AK_OS_SERENITY) | 350 | | if (__builtin_isnan(y)) | 351 | | return y; | 352 | | if (__builtin_isnan(x)) | 353 | | return x; | 354 | | | 355 | | // SPECIAL VALUES | 356 | | // atan2(±0, -0) returns ±pi. | 357 | | if (y == 0 && x == 0 && signbit(x)) | 358 | | return copysign(Pi<T>, y); | 359 | | | 360 | | // atan2(±0, +0) returns ±0. | 361 | | if (y == 0 && x == 0 && !signbit(x)) | 362 | | return y; | 363 | | | 364 | | // atan2(±0, x) returns ±pi for x < 0. | 365 | | if (y == 0 && x < 0) | 366 | | return copysign(Pi<T>, y); | 367 | | | 368 | | // atan2(±0, x) returns ±0 for x > 0. | 369 | | if (y == 0 && x > 0) | 370 | | return y; | 371 | | | 372 | | // atan2(y, ±0) returns +pi/2 for y > 0. | 373 | | if (y > 0 && x == 0) | 374 | | return Pi<T> / 2; | 375 | | | 376 | | // atan2(y, ±0) returns -pi/2 for y < 0. | 377 | | if (y < 0 && x == 0) | 378 | | return -Pi<T> / 2; | 379 | | | 380 | | // atan2(±y, -infinity) returns ±pi for finite y > 0. | 381 | | if (!__builtin_isinf(y) && y > 0 && __builtin_isinf(x) && signbit(x)) | 382 | | return copysign(Pi<T>, y); | 383 | | | 384 | | // atan2(±y, +infinity) returns ±0 for finite y > 0. | 385 | | if (!__builtin_isinf(y) && y > 0 && __builtin_isinf(x) && !signbit(x)) | 386 | | return copysign(static_cast<T>(0), y); | 387 | | | 388 | | // atan2(±infinity, x) returns ±pi/2 for finite x. | 389 | | if (__builtin_isinf(y) && !__builtin_isinf(x)) | 390 | | return copysign(Pi<T> / 2, y); | 391 | | | 392 | | // atan2(±infinity, -infinity) returns ±3*pi/4. | 393 | | if (__builtin_isinf(y) && __builtin_isinf(x) && signbit(x)) | 394 | | return copysign(3 * Pi<T> / 4, y); | 395 | | | 396 | | // atan2(±infinity, +infinity) returns ±pi/4. | 397 | | if (__builtin_isinf(y) && __builtin_isinf(x) && !signbit(x)) | 398 | | return copysign(Pi<T> / 4, y); | 399 | | | 400 | | // Check quadrant, going counterclockwise. | 401 | | if (y > 0 && x > 0) | 402 | | return atan(y / x); | 403 | | if (y > 0 && x < 0) | 404 | | return atan(y / x) + Pi<T>; | 405 | | if (y < 0 && x < 0) | 406 | | return atan(y / x) - Pi<T>; | 407 | | // y < 0 && x > 0 | 408 | | return atan(y / x); | 409 | | # else | 410 | | CALL_BUILTIN(atan2, y, x); | 411 | | # endif | 412 | | #endif | 413 | 3.05k | } |
_ZN2AK12Trigonometry5atan2ITkNS_8Concepts13FloatingPointEdEET_S3_S3_ Line | Count | Source | 338 | 2.75M | { | 339 | 2.75M | CONSTEXPR_STATE(atan2, y, x); | 340 | | | 341 | 2.75M | #if ARCH(X86_64) | 342 | 2.75M | T ret; | 343 | 2.75M | asm("fpatan" | 344 | 2.75M | : "=t"(ret) | 345 | 2.75M | : "0"(x), "u"(y) | 346 | 2.75M | : "st(1)"); | 347 | 2.75M | return ret; | 348 | | #else | 349 | | # if defined(AK_OS_SERENITY) | 350 | | if (__builtin_isnan(y)) | 351 | | return y; | 352 | | if (__builtin_isnan(x)) | 353 | | return x; | 354 | | | 355 | | // SPECIAL VALUES | 356 | | // atan2(±0, -0) returns ±pi. | 357 | | if (y == 0 && x == 0 && signbit(x)) | 358 | | return copysign(Pi<T>, y); | 359 | | | 360 | | // atan2(±0, +0) returns ±0. | 361 | | if (y == 0 && x == 0 && !signbit(x)) | 362 | | return y; | 363 | | | 364 | | // atan2(±0, x) returns ±pi for x < 0. | 365 | | if (y == 0 && x < 0) | 366 | | return copysign(Pi<T>, y); | 367 | | | 368 | | // atan2(±0, x) returns ±0 for x > 0. | 369 | | if (y == 0 && x > 0) | 370 | | return y; | 371 | | | 372 | | // atan2(y, ±0) returns +pi/2 for y > 0. | 373 | | if (y > 0 && x == 0) | 374 | | return Pi<T> / 2; | 375 | | | 376 | | // atan2(y, ±0) returns -pi/2 for y < 0. | 377 | | if (y < 0 && x == 0) | 378 | | return -Pi<T> / 2; | 379 | | | 380 | | // atan2(±y, -infinity) returns ±pi for finite y > 0. | 381 | | if (!__builtin_isinf(y) && y > 0 && __builtin_isinf(x) && signbit(x)) | 382 | | return copysign(Pi<T>, y); | 383 | | | 384 | | // atan2(±y, +infinity) returns ±0 for finite y > 0. | 385 | | if (!__builtin_isinf(y) && y > 0 && __builtin_isinf(x) && !signbit(x)) | 386 | | return copysign(static_cast<T>(0), y); | 387 | | | 388 | | // atan2(±infinity, x) returns ±pi/2 for finite x. | 389 | | if (__builtin_isinf(y) && !__builtin_isinf(x)) | 390 | | return copysign(Pi<T> / 2, y); | 391 | | | 392 | | // atan2(±infinity, -infinity) returns ±3*pi/4. | 393 | | if (__builtin_isinf(y) && __builtin_isinf(x) && signbit(x)) | 394 | | return copysign(3 * Pi<T> / 4, y); | 395 | | | 396 | | // atan2(±infinity, +infinity) returns ±pi/4. | 397 | | if (__builtin_isinf(y) && __builtin_isinf(x) && !signbit(x)) | 398 | | return copysign(Pi<T> / 4, y); | 399 | | | 400 | | // Check quadrant, going counterclockwise. | 401 | | if (y > 0 && x > 0) | 402 | | return atan(y / x); | 403 | | if (y > 0 && x < 0) | 404 | | return atan(y / x) + Pi<T>; | 405 | | if (y < 0 && x < 0) | 406 | | return atan(y / x) - Pi<T>; | 407 | | // y < 0 && x > 0 | 408 | | return atan(y / x); | 409 | | # else | 410 | | CALL_BUILTIN(atan2, y, x); | 411 | | # endif | 412 | | #endif | 413 | 2.75M | } |
|
414 | | |
415 | | } |
416 | | |
417 | | #include <AK/Math/UndefMacros.h> |
418 | | |
419 | | using Trigonometry::acos; |
420 | | using Trigonometry::asin; |
421 | | using Trigonometry::atan; |
422 | | using Trigonometry::atan2; |
423 | | using Trigonometry::cos; |
424 | | using Trigonometry::hypot; |
425 | | using Trigonometry::sin; |
426 | | using Trigonometry::sincos; |
427 | | using Trigonometry::tan; |
428 | | |
429 | | } |