sysdeps/ieee754/flt-32/sincosf_poly.h

   1 /* Used by sinf, cosf and sincosf functions.
   2    Copyright (C) 2018-2024 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 /* The constants and polynomials for sine and cosine.  */
  20 typedef struct
  21 {
  22   double sign[4];               /* Sign of sine in quadrants 0..3.  */
  23   double hpi_inv;               /* 2 / PI ( * 2^24 if !TOINT_INTRINSICS).  */
  24   double hpi;                   /* PI / 2.  */
  25   double c0, c1, c2, c3, c4;    /* Cosine polynomial.  */
  26   double s1, s2, s3;            /* Sine polynomial.  */
  27 } sincos_t;
  28
  29 /* Compute the sine and cosine of inputs X and X2 (X squared), using the
  30    polynomial P and store the results in SINP and COSP.  N is the quadrant,
  31    if odd the cosine and sine polynomials are swapped.  */
  32 static inline void
  33 sincosf_poly (double x, double x2, const sincos_t *p, int n, float *sinp,
  34               float *cosp)
  35 {
  36   double x3, x4, x5, x6, s, c, c1, c2, s1;
  37
  38   x4 = x2 * x2;
  39   x3 = x2 * x;
  40   c2 = p->c3 + x2 * p->c4;
  41   s1 = p->s2 + x2 * p->s3;
  42
  43   /* Swap sin/cos result based on quadrant.  */
  44   float *tmp = (n & 1 ? cosp : sinp);
  45   cosp = (n & 1 ? sinp : cosp);
  46   sinp = tmp;
  47
  48   c1 = p->c0 + x2 * p->c1;
  49   x5 = x3 * x2;
  50   x6 = x4 * x2;
  51
  52   s = x + x3 * p->s1;
  53   c = c1 + x4 * p->c2;
  54
  55   *sinp = s + x5 * s1;
  56   *cosp = c + x6 * c2;
  57 }
  58
  59 /* Return the sine of inputs X and X2 (X squared) using the polynomial P.
  60    N is the quadrant, and if odd the cosine polynomial is used.  */
  61 static inline float
  62 sinf_poly (double x, double x2, const sincos_t *p, int n)
  63 {
  64   double x3, x4, x6, x7, s, c, c1, c2, s1;
  65
  66   if ((n & 1) == 0)
  67     {
  68       x3 = x * x2;
  69       s1 = p->s2 + x2 * p->s3;
  70
  71       x7 = x3 * x2;
  72       s = x + x3 * p->s1;
  73
  74       return s + x7 * s1;
  75     }
  76   else
  77     {
  78       x4 = x2 * x2;
  79       c2 = p->c3 + x2 * p->c4;
  80       c1 = p->c0 + x2 * p->c1;
  81
  82       x6 = x4 * x2;
  83       c = c1 + x4 * p->c2;
  84
  85       return c + x6 * c2;
  86     }
  87 }