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 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