1  /* Used by sinf, cosf and sincosf functions.
       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  /* 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  }