sysdeps/ieee754/flt-32/s_cosf.c

   1 /* Compute cosine of argument.
   2    Copyright (C) 2018 Free Software Foundation, Inc.
   3    This file is part of the GNU C Library.
   4
   5    The GNU C Library is free software; you can redistribute it and/or
   6    modify it under the terms of the GNU Lesser General Public
   7    License as published by the Free Software Foundation; either
   8    version 2.1 of the License, or (at your option) any later version.
   9
  10    The GNU C Library is distributed in the hope that it will be useful,
  11    but WITHOUT ANY WARRANTY; without even the implied warranty of
  12    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13    Lesser General Public License for more details.
  14
  15    You should have received a copy of the GNU Lesser General Public
  16    License along with the GNU C Library; if not, see
  17    <http://www.gnu.org/licenses/>.  */
  18
  19 #include <stdint.h>
  20 #include <math.h>
  21 #include <math-barriers.h>
  22 #include <libm-alias-float.h>
  23 #include "math_config.h"
  24 #include "s_sincosf.h"
  25
  26 #ifndef COSF
  27 # define COSF_FUNC __cosf
  28 #else
  29 # define COSF_FUNC COSF
  30 #endif
  31
  32 /* Fast cosf implementation.  Worst-case ULP is 0.5607, maximum relative
  33    error is 0.5303 * 2^-23.  A single-step range reduction is used for
  34    small values.  Large inputs have their range reduced using fast integer
  35    arithmetic.
  36 */
  37 float
  38 COSF_FUNC (float y)
  39 {
  40   double x = y;
  41   double s;
  42   int n;
  43   const sincos_t *p = &__sincosf_table[0];
  44
  45   if (abstop12 (y) < abstop12 (pio4))
  46     {
  47       double x2 = x * x;
  48
  49       if (__glibc_unlikely (abstop12 (y) < abstop12 (0x1p-12f)))
  50         return 1.0f;
  51
  52       return sinf_poly (x, x2, p, 1);
  53     }
  54   else if (__glibc_likely (abstop12 (y) < abstop12 (120.0f)))
  55     {
  56       x = reduce_fast (x, p, &n);
  57
  58       /* Setup the signs for sin and cos.  */
  59       s = p->sign[n & 3];
  60
  61       if (n & 2)
  62         p = &__sincosf_table[1];
  63
  64       return sinf_poly (x * s, x * x, p, n ^ 1);
  65     }
  66   else if (abstop12 (y) < abstop12 (INFINITY))
  67     {
  68       uint32_t xi = asuint (y);
  69       int sign = xi >> 31;
  70
  71       x = reduce_large (xi, &n);
  72
  73       /* Setup signs for sin and cos - include original sign.  */
  74       s = p->sign[(n + sign) & 3];
  75
  76       if ((n + sign) & 2)
  77         p = &__sincosf_table[1];
  78
  79       return sinf_poly (x * s, x * x, p, n ^ 1);
  80     }
  81   else
  82     return __math_invalidf (y);
  83 }
  84
  85 #ifndef COSF
  86 libm_alias_float (__cos, cos)
  87 #endif