1  /* Compute sine and cosine of argument.
       2     Copyright (C) 1997-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 <errno.h>
      20  #include <math.h>
      21  
      22  #include <math_private.h>
      23  #include <libm-alias-ldouble.h>
      24  
      25  
      26  void
      27  __sincosl (long double x, long double *sinx, long double *cosx)
      28  {
      29    int32_t se, i0, i1 __attribute__ ((unused));
      30  
      31    /* High word of x. */
      32    GET_LDOUBLE_WORDS (se, i0, i1, x);
      33  
      34    /* |x| ~< pi/4 */
      35    se &= 0x7fff;
      36    if (se < 0x3ffe || (se == 0x3ffe && i0 <= 0xc90fdaa2))
      37      {
      38        *sinx = __kernel_sinl (x, 0.0, 0);
      39        *cosx = __kernel_cosl (x, 0.0);
      40      }
      41    else if (se == 0x7fff)
      42      {
      43        /* sin(Inf or NaN) is NaN */
      44        *sinx = *cosx = x - x;
      45        if (isinf (x))
      46  	__set_errno (EDOM);
      47      }
      48    else
      49      {
      50        /* Argument reduction needed.  */
      51        long double y[2];
      52        int n;
      53  
      54        n = __ieee754_rem_pio2l (x, y);
      55        switch (n & 3)
      56  	{
      57  	case 0:
      58  	  *sinx = __kernel_sinl (y[0], y[1], 1);
      59  	  *cosx = __kernel_cosl (y[0], y[1]);
      60  	  break;
      61  	case 1:
      62  	  *sinx = __kernel_cosl (y[0], y[1]);
      63  	  *cosx = -__kernel_sinl (y[0], y[1], 1);
      64  	  break;
      65  	case 2:
      66  	  *sinx = -__kernel_sinl (y[0], y[1], 1);
      67  	  *cosx = -__kernel_cosl (y[0], y[1]);
      68  	  break;
      69  	default:
      70  	  *sinx = -__kernel_cosl (y[0], y[1]);
      71  	  *cosx = __kernel_sinl (y[0], y[1], 1);
      72  	  break;
      73  	}
      74      }
      75  }
      76  libm_alias_ldouble (__sincos, sincos)