1 /*
2 * s_sincosf.c - single precision sine and cosine functions
3 *
4 * Copyright (c) 2009-2018, Arm Limited.
5 * SPDX-License-Identifier: MIT
6 */
7
8 /*
9 * Source: my own head, and Remez-generated polynomial approximations.
10 */
11
12 #include <fenv.h>
13 #include <math.h>
14 #include <errno.h>
15 #include "rredf.h"
16 #include "math_private.h"
17
18 #ifdef __cplusplus
19 extern "C" {
20 #endif /* __cplusplus */
21
22 #ifndef COSINE
23 #define FUNCNAME sinf
24 #define SOFTFP_FUNCNAME __softfp_sinf
25 #define DO_SIN (!(q & 1))
26 #define NEGATE_SIN ((q & 2))
27 #define NEGATE_COS ((q & 2))
28 #define TRIVIAL_RESULT(x) FLOAT_CHECKDENORM(x)
29 #define ERR_INF MATHERR_SINF_INF
30 #else
31 #define FUNCNAME cosf
32 #define SOFTFP_FUNCNAME __softfp_cosf
33 #define DO_SIN (q & 1)
34 #define NEGATE_SIN (!(q & 2))
35 #define NEGATE_COS ((q & 2))
36 #define TRIVIAL_RESULT(x) 1.0f
37 #define ERR_INF MATHERR_COSF_INF
38 #endif
39
FUNCNAME(float x)40 float FUNCNAME(float x)
41 {
42 int q;
43
44 /*
45 * Range-reduce x to the range [-pi/4,pi/4].
46 */
47 {
48 /*
49 * I enclose the call to __mathlib_rredf in braces so that
50 * the address-taken-ness of qq does not propagate
51 * throughout the rest of the function, for what that might
52 * be worth.
53 */
54 int qq;
55 x = __mathlib_rredf(x, &qq);
56 q = qq;
57 }
58 if (__builtin_expect(q < 0, 0)) { /* this signals tiny, inf, or NaN */
59 unsigned k = fai(x) << 1;
60 if (k < 0xFF000000) /* tiny */
61 return TRIVIAL_RESULT(x);
62 else if (k == 0xFF000000) /* inf */
63 return ERR_INF(x);
64 else /* NaN */
65 return FLOAT_INFNAN(x);
66 }
67
68 /*
69 * Depending on the quadrant we were in, we may have to compute
70 * a sine-like function (f(0)=0) or a cosine-like one (f(0)=1),
71 * and we may have to negate it.
72 */
73 if (DO_SIN) {
74 float x2 = x*x;
75 /*
76 * Coefficients generated by the command
77
78 ./auxiliary/remez.jl --variable=x2 --suffix=f -- '0' 'atan(BigFloat(1))^2' 2 0 'x==0 ? -1/BigFloat(6) : (x->(sin(x)-x)/x^3)(sqrt(x))' 'sqrt(x^3)'
79
80 */
81 x += x * (x2 * (
82 -1.666665066929417292436220415142244613956015227491999719404711781344783392564922e-01f+x2*(8.331978663157089651408875887703995477889496917296385733254577121461421466427772e-03f+x2*(-1.949563623766929906511886482584265500187314705496861877317774185883215997931494e-04f))
83 ));
84 if (NEGATE_SIN)
85 x = -x;
86 return x;
87 } else {
88 float x2 = x*x;
89 /*
90 * Coefficients generated by the command
91
92 ./auxiliary/remez.jl --variable=x2 --suffix=f -- '0' 'atan(BigFloat(1))^2' 2 0 'x==0 ? -1/BigFloat(6) : (x->(cos(x)-1)/x^2)(sqrt(x))' 'x'
93
94 */
95 x = 1.0f + x2*(
96 -4.999989478137016757327030935768632852012781143541026304540837816323349768666875e-01f+x2*(4.165629457842617238353362092016541041535652603456375154392942188742496860024377e-02f+x2*(-1.35978231111049428748566568960423202948250988565693107969571193763372093404347e-03f))
97 );
98 if (NEGATE_COS)
99 x = -x;
100 return x;
101 }
102 }
103
104
105 #ifdef __cplusplus
106 } /* end of extern "C" */
107 #endif /* __cplusplus */
108
109 /* end of sincosf.c */
110