Initial commit

test_audio_wasm/math.c

@@ -0,0 +1,59 @@
+#include "math.h"
+
+__attribute__((always_inline))
+static float fminf(float a, float b) {
+    return a < b ? a : b;
+}
+
+/* Approach 1: polynomial; input in radians. */
+static float sin_poly_radians(float x)
+{
+    x *= 0.159155f;          /* approximately 1 / (2*pi) */
+    x -= __builtin_floorf(x);          /* reduce to [0, 1) */
+    float xx = x * x;
+    float y = -6.87897f;
+    y = y * xx + 33.7755f;
+    y = y * xx - 72.5257f;
+    y = y * xx + 80.5874f;
+    y = y * xx - 41.2408f;
+    y = y * xx + 6.28077f;
+    return x * y;
+}
+
+// Approach 2: rational approximation; input in turns (1.0 = 2*pi radians)
+// Calculate sin(tau x) (where x is in the range [0, 1]) using Bhaskara I's approximation
+// https://en.wikipedia.org/wiki/Bhaskara_I%27s_sine_approximation_formula
+// Accurate to about 0.2%
+static float sin_rational_tau(float x)
+{
+    float frac = x - __builtin_floorf(x);
+    float sign = 0.5f - frac;
+    float y = 2.0f * fminf(frac, 1.0f - frac);
+    float q = y * (1.0f - y);
+    float approximation = (16.0f * q) / (5.0f - 4.0f * q);
+    return __builtin_copysignf(approximation, sign);
+}
+
+/* Approach 3: folded polynomial; input in turns (1.0 = 2*pi radians). */
+static float sin_folded_poly_tau(float x)
+{
+    union { uint32_t i; float f; } coeffs[6] = {
+        { 0x4261c73fu }, { 0x4113e0d9u }, { 0xc22aaa4bu },
+        { 0x3da6b273u }, { 0x40c90140u }, { 0x36d0957eu }
+    };
+
+    float frac = x - __builtin_floorf(x);
+    float x1 = __builtin_fabsf(frac - 0.5f);
+    float t = fminf(0.5f - x1, x1); /* fold to [0, 0.25] */
+
+    float y = 0.0f;
+    for (int i = 0; i < 6; ++i) {
+        y = y * t + coeffs[i].f;
+    }
+
+    return __builtin_copysignf(y, 0.5f - frac);
+}
+
+float sintau(float x) {
+    return sin_rational_tau(x);
+}