Home | History | Annotate | Line # | Download | only in src 
      1 /*-
      2  * Copyright (c) 2007, 2010-2013 Steven G. Kargl
      3  * All rights reserved.
      4  *
      5  * Redistribution and use in source and binary forms, with or without
      6  * modification, are permitted provided that the following conditions
      7  * are met:
      8  * 1. Redistributions of source code must retain the above copyright
      9  *    notice unmodified, this list of conditions, and the following
     10  *    disclaimer.
     11  * 2. Redistributions in binary form must reproduce the above copyright
     12  *    notice, this list of conditions and the following disclaimer in the
     13  *    documentation and/or other materials provided with the distribution.
     14  *
     15  * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
     16  * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
     17  * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
     18  * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT,
     19  * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
     20  * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
     21  * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
     22  * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
     23  * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
     24  * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
     25  *
     26  * s_sinl.c and s_cosl.c merged by Steven G. Kargl.
     27  */
     28 
     29 #include <sys/cdefs.h>
     30 #if 0
     31 __FBSDID("$FreeBSD: head/lib/msun/src/s_sincosl.c 319047 2017-05-28 06:13:38Z mmel $");
     32 #endif
     33 #if defined(LIBM_SCCS) && !defined(lint)
     34 __RCSID("$NetBSD: s_sincosl.c,v 1.4 2024/04/03 18:53:42 christos Exp $");
     35 #endif
     36 
     37 #include "namespace.h"
     38 #include <float.h>
     39 
     40 #include "math.h"
     41 #include "math_private.h"
     42 
     43 #ifdef __weak_alias
     44 __weak_alias(sincosl,_sincosl)
     45 #endif
     46 
     47 #ifdef __HAVE_LONG_DOUBLE
     48 
     49 #include "k_sincosl.h"
     50 
     51 #if LDBL_MANT_DIG == 64
     52 #include "../ld80/e_rem_pio2l.h"
     53 #elif LDBL_MANT_DIG == 113
     54 #include "../ld128/e_rem_pio2l.h"
     55 #else
     56 #error "Unsupported long double format"
     57 #endif
     58 
     59 void
     60 sincosl(long double x, long double *sn, long double *cs)
     61 {
     62 	union ieee_ext_u z;
     63 	int e0;
     64 	long double y[2];
     65 
     66 	z.extu_ld = x;
     67 	z.extu_sign = 0;
     68 
     69 	ENTERV();
     70 
     71 	/* Optimize the case where x is already within range. */
     72 	if (z.extu_ld < M_PI_4) {
     73 		/*
     74 		 * If x = +-0 or x is a subnormal number, then sin(x) = x and
     75 		 * cos(x) = 1.
     76 		 */
     77 		if (z.extu_exp == 0) {
     78 			*sn = x;
     79 			*cs = 1;
     80 		} else
     81 			__kernel_sincosl(x, 0, 0, sn, cs);
     82 		RETURNV();
     83 	}
     84 
     85 	/* If x = NaN or Inf, then sin(x) and cos(x) are NaN. */
     86 	if (z.extu_exp == 32767) {
     87 		*sn = x - x;
     88 		*cs = x - x;
     89 		RETURNV();
     90 	}
     91 
     92 	/* Range reduction. */
     93 	e0 = __ieee754_rem_pio2l(x, y);
     94 
     95 	switch (e0 & 3) {
     96 	case 0:
     97 		__kernel_sincosl(y[0], y[1], 1, sn, cs);
     98 		break;
     99 	case 1:
    100 		__kernel_sincosl(y[0], y[1], 1, cs, sn);
    101 		*cs = -*cs;
    102 		break;
    103 	case 2:
    104 		__kernel_sincosl(y[0], y[1], 1, sn, cs);
    105 		*sn = -*sn;
    106 		*cs = -*cs;
    107 		break;
    108 	default:
    109 		__kernel_sincosl(y[0], y[1], 1, cs, sn);
    110 		*sn = -*sn;
    111 	}
    112 
    113 	RETURNV();
    114 }
    115 
    116 #else
    117 
    118 void
    119 sincosl(long double x, long double *sn, long double *cs)
    120 {
    121 	double snn, css;
    122 	sincos(x, &snn, &css);
    123 	*sn = snn;
    124 	*cs = css;
    125 }
    126 
    127 #endif
    128