1  /* Compute cosine of argument.
       2     Copyright (C) 2018-2023 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     <https://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 SECTION
      27  # define SECTION
      28  #endif
      29  
      30  #ifndef COSF
      31  # define COSF_FUNC __cosf
      32  #else
      33  # define COSF_FUNC COSF
      34  #endif
      35  
      36  /* Fast cosf implementation.  Worst-case ULP is 0.5607, maximum relative
      37     error is 0.5303 * 2^-23.  A single-step range reduction is used for
      38     small values.  Large inputs have their range reduced using fast integer
      39     arithmetic.
      40  */
      41  float
      42  SECTION
      43  COSF_FUNC (float y)
      44  {
      45    double x = y;
      46    double s;
      47    int n;
      48    const sincos_t *p = &__sincosf_table[0];
      49  
      50    if (abstop12 (y) < abstop12 (pio4))
      51      {
      52        double x2 = x * x;
      53  
      54        if (__glibc_unlikely (abstop12 (y) < abstop12 (0x1p-12f)))
      55  	return 1.0f;
      56  
      57        return sinf_poly (x, x2, p, 1);
      58      }
      59    else if (__glibc_likely (abstop12 (y) < abstop12 (120.0f)))
      60      {
      61        x = reduce_fast (x, p, &n);
      62  
      63        /* Setup the signs for sin and cos.  */
      64        s = p->sign[n & 3];
      65  
      66        if (n & 2)
      67  	p = &__sincosf_table[1];
      68  
      69        return sinf_poly (x * s, x * x, p, n ^ 1);
      70      }
      71    else if (abstop12 (y) < abstop12 (INFINITY))
      72      {
      73        uint32_t xi = asuint (y);
      74        int sign = xi >> 31;
      75  
      76        x = reduce_large (xi, &n);
      77  
      78        /* Setup signs for sin and cos - include original sign.  */
      79        s = p->sign[(n + sign) & 3];
      80  
      81        if ((n + sign) & 2)
      82  	p = &__sincosf_table[1];
      83  
      84        return sinf_poly (x * s, x * x, p, n ^ 1);
      85      }
      86    else
      87      return __math_invalidf (y);
      88  }
      89  
      90  #ifndef COSF
      91  libm_alias_float (__cos, cos)
      92  #endif