freertos: release the generic version source code
freertos runs on the second core (small one) of the CPU
This commit is contained in:
@ -0,0 +1,153 @@
|
||||
/*
|
||||
* Header for sinf, cosf and sincosf.
|
||||
*
|
||||
* Copyright (c) 2018, Arm Limited.
|
||||
* SPDX-License-Identifier: MIT
|
||||
*/
|
||||
|
||||
#include <stdint.h>
|
||||
#include <math.h>
|
||||
#include "math_config.h"
|
||||
|
||||
/* 2PI * 2^-64. */
|
||||
static const double pi63 = 0x1.921FB54442D18p-62;
|
||||
/* PI / 4. */
|
||||
static const double pio4 = 0x1.921FB54442D18p-1;
|
||||
|
||||
/* The constants and polynomials for sine and cosine. */
|
||||
typedef struct
|
||||
{
|
||||
double sign[4]; /* Sign of sine in quadrants 0..3. */
|
||||
double hpi_inv; /* 2 / PI ( * 2^24 if !TOINT_INTRINSICS). */
|
||||
double hpi; /* PI / 2. */
|
||||
double c0, c1, c2, c3, c4; /* Cosine polynomial. */
|
||||
double s1, s2, s3; /* Sine polynomial. */
|
||||
} sincos_t;
|
||||
|
||||
/* Polynomial data (the cosine polynomial is negated in the 2nd entry). */
|
||||
extern const sincos_t __sincosf_table[2] HIDDEN;
|
||||
|
||||
/* Table with 4/PI to 192 bit precision. */
|
||||
extern const uint32_t __inv_pio4[] HIDDEN;
|
||||
|
||||
/* Top 12 bits of the float representation with the sign bit cleared. */
|
||||
static inline uint32_t
|
||||
abstop12 (float x)
|
||||
{
|
||||
return (asuint (x) >> 20) & 0x7ff;
|
||||
}
|
||||
|
||||
/* Compute the sine and cosine of inputs X and X2 (X squared), using the
|
||||
polynomial P and store the results in SINP and COSP. N is the quadrant,
|
||||
if odd the cosine and sine polynomials are swapped. */
|
||||
static inline void
|
||||
sincosf_poly (double x, double x2, const sincos_t *p, int n, float *sinp,
|
||||
float *cosp)
|
||||
{
|
||||
double x3, x4, x5, x6, s, c, c1, c2, s1;
|
||||
|
||||
x4 = x2 * x2;
|
||||
x3 = x2 * x;
|
||||
c2 = p->c3 + x2 * p->c4;
|
||||
s1 = p->s2 + x2 * p->s3;
|
||||
|
||||
/* Swap sin/cos result based on quadrant. */
|
||||
float *tmp = (n & 1 ? cosp : sinp);
|
||||
cosp = (n & 1 ? sinp : cosp);
|
||||
sinp = tmp;
|
||||
|
||||
c1 = p->c0 + x2 * p->c1;
|
||||
x5 = x3 * x2;
|
||||
x6 = x4 * x2;
|
||||
|
||||
s = x + x3 * p->s1;
|
||||
c = c1 + x4 * p->c2;
|
||||
|
||||
*sinp = s + x5 * s1;
|
||||
*cosp = c + x6 * c2;
|
||||
}
|
||||
|
||||
/* Return the sine of inputs X and X2 (X squared) using the polynomial P.
|
||||
N is the quadrant, and if odd the cosine polynomial is used. */
|
||||
static inline float
|
||||
sinf_poly (double x, double x2, const sincos_t *p, int n)
|
||||
{
|
||||
double x3, x4, x6, x7, s, c, c1, c2, s1;
|
||||
|
||||
if ((n & 1) == 0)
|
||||
{
|
||||
x3 = x * x2;
|
||||
s1 = p->s2 + x2 * p->s3;
|
||||
|
||||
x7 = x3 * x2;
|
||||
s = x + x3 * p->s1;
|
||||
|
||||
return s + x7 * s1;
|
||||
}
|
||||
else
|
||||
{
|
||||
x4 = x2 * x2;
|
||||
c2 = p->c3 + x2 * p->c4;
|
||||
c1 = p->c0 + x2 * p->c1;
|
||||
|
||||
x6 = x4 * x2;
|
||||
c = c1 + x4 * p->c2;
|
||||
|
||||
return c + x6 * c2;
|
||||
}
|
||||
}
|
||||
|
||||
/* Fast range reduction using single multiply-subtract. Return the modulo of
|
||||
X as a value between -PI/4 and PI/4 and store the quadrant in NP.
|
||||
The values for PI/2 and 2/PI are accessed via P. Since PI/2 as a double
|
||||
is accurate to 55 bits and the worst-case cancellation happens at 6 * PI/4,
|
||||
the result is accurate for |X| <= 120.0. */
|
||||
static inline double
|
||||
reduce_fast (double x, const sincos_t *p, int *np)
|
||||
{
|
||||
double r;
|
||||
#if TOINT_INTRINSICS
|
||||
/* Use fast round and lround instructions when available. */
|
||||
r = x * p->hpi_inv;
|
||||
*np = converttoint (r);
|
||||
return x - roundtoint (r) * p->hpi;
|
||||
#else
|
||||
/* Use scaled float to int conversion with explicit rounding.
|
||||
hpi_inv is prescaled by 2^24 so the quadrant ends up in bits 24..31.
|
||||
This avoids inaccuracies introduced by truncating negative values. */
|
||||
r = x * p->hpi_inv;
|
||||
int n = ((int32_t)r + 0x800000) >> 24;
|
||||
*np = n;
|
||||
return x - n * p->hpi;
|
||||
#endif
|
||||
}
|
||||
|
||||
/* Reduce the range of XI to a multiple of PI/2 using fast integer arithmetic.
|
||||
XI is a reinterpreted float and must be >= 2.0f (the sign bit is ignored).
|
||||
Return the modulo between -PI/4 and PI/4 and store the quadrant in NP.
|
||||
Reduction uses a table of 4/PI with 192 bits of precision. A 32x96->128 bit
|
||||
multiply computes the exact 2.62-bit fixed-point modulo. Since the result
|
||||
can have at most 29 leading zeros after the binary point, the double
|
||||
precision result is accurate to 33 bits. */
|
||||
static inline double
|
||||
reduce_large (uint32_t xi, int *np)
|
||||
{
|
||||
const uint32_t *arr = &__inv_pio4[(xi >> 26) & 15];
|
||||
int shift = (xi >> 23) & 7;
|
||||
uint64_t n, res0, res1, res2;
|
||||
|
||||
xi = (xi & 0xffffff) | 0x800000;
|
||||
xi <<= shift;
|
||||
|
||||
res0 = xi * arr[0];
|
||||
res1 = (uint64_t)xi * arr[4];
|
||||
res2 = (uint64_t)xi * arr[8];
|
||||
res0 = (res2 >> 32) | (res0 << 32);
|
||||
res0 += res1;
|
||||
|
||||
n = (res0 + (1ULL << 61)) >> 62;
|
||||
res0 -= n << 62;
|
||||
double x = (int64_t)res0;
|
||||
*np = n;
|
||||
return x * pi63;
|
||||
}
|
||||
Reference in New Issue
Block a user