svn commit: r319047 - in head/lib/msun: . man src

Michal Meloun mmel at FreeBSD.org
Sun May 28 06:13:41 UTC 2017


Author: mmel
Date: Sun May 28 06:13:38 2017
New Revision: 319047
URL: https://svnweb.freebsd.org/changeset/base/319047

Log:
  Implement sincos, sincosf, and sincosl.
  The primary benefit of these functions is that argument
  reduction is done once instead of twice in independent
  calls to sin() and cos().
  
  * lib/msun/Makefile:
    . Add s_sincos[fl].c to the build.
    . Add sincos.3 documentation.
    . Add appropriate MLINKS.
  
  * lib/msun/Symbol.map:
    . Expose sincos[fl] symbols in dynamic libm.so.
  
  * lib/msun/man/sincos.3:
    . Documentation for sincos[fl].
  
  * lib/msun/src/k_sincos.h:
    . Kernel for sincos() function.  This merges the individual kernels
      for sin() and cos().  The merger offered an opportunity to re-arrange
      the individual kernels for better performance.
  
  * lib/msun/src/k_sincosf.h:
     . Kernel for sincosf() function.  This merges the individual kernels
       for sinf() and cosf(). The merger offered an opportunity to re-arrange
       the individual kernels for better performance.
  
  * lib/msun/src/k_sincosl.h:
     . Kernel for sincosl() function.  This merges the individual kernels
       for sinl() and cosl(). The merger offered an opportunity to re-arrange
       the individual kernels for better performance.
  
  * lib/msun/src/math.h:
    . Add prototytpes for sincos[fl]().
  
  * lib/msun/src/math_private.h:
    . Add RETURNV macros.  This is needed to reset fpsetprec on I386
      hardware for a function with type void.
  
  * lib/msun/src/s_sincos.c:
    . Implementation of sincos() where sin() and cos() were merged into
      one routine and possibly re-arranged for better performance.
  
  * lib/msun/src/s_sincosf.c:
    . Implementation of sincosf() where sinf() and cosf() were merged into
      one routine and possibly re-arranged for better performance.
  
  * lib/msun/src/s_sincosl.c:
    . Implementation of sincosl() where sinl() and cosl() were merged into
      one routine and possibly re-arranged for better performance.
  
  PR:		215977, 218300
  Submitted by:	Steven G. Kargl <sgk at troutmask.apl.washington.edu>
  MFC after:	1 month
  Differential Revision:	https://reviews.freebsd.org/D10765

Added:
  head/lib/msun/man/sincos.3   (contents, props changed)
  head/lib/msun/src/k_sincos.h   (contents, props changed)
  head/lib/msun/src/k_sincosf.h   (contents, props changed)
  head/lib/msun/src/k_sincosl.h   (contents, props changed)
  head/lib/msun/src/s_sincos.c   (contents, props changed)
  head/lib/msun/src/s_sincosf.c   (contents, props changed)
  head/lib/msun/src/s_sincosl.c   (contents, props changed)
Modified:
  head/lib/msun/Makefile
  head/lib/msun/Symbol.map
  head/lib/msun/src/math.h
  head/lib/msun/src/math_private.h

Modified: head/lib/msun/Makefile
==============================================================================
--- head/lib/msun/Makefile	Sun May 28 05:31:18 2017	(r319046)
+++ head/lib/msun/Makefile	Sun May 28 06:13:38 2017	(r319047)
@@ -73,7 +73,8 @@ COMMON_SRCS= b_exp.c b_log.c b_tgamma.c 
 	s_nexttowardf.c s_remquo.c s_remquof.c \
 	s_rint.c s_rintf.c s_round.c s_roundf.c \
 	s_scalbln.c s_scalbn.c s_scalbnf.c s_signbit.c \
-	s_signgam.c s_significand.c s_significandf.c s_sin.c s_sinf.c \
+	s_signgam.c s_significand.c s_significandf.c s_sin.c \
+	s_sincos.c s_sincosf.c s_sinf.c \
 	s_tan.c s_tanf.c s_tanh.c s_tanhf.c s_tgammaf.c s_trunc.c s_truncf.c \
 	w_cabs.c w_cabsf.c w_drem.c w_dremf.c
 
@@ -104,7 +105,8 @@ COMMON_SRCS+=	catrigl.c \
 	s_csqrtl.c s_erfl.c s_exp2l.c s_expl.c s_floorl.c s_fmal.c \
 	s_fmaxl.c s_fminl.c s_frexpl.c s_logbl.c s_logl.c s_nanl.c \
 	s_nextafterl.c s_nexttoward.c s_remquol.c s_rintl.c s_roundl.c \
-	s_scalbnl.c s_sinl.c s_tanhl.c s_tanl.c s_truncl.c w_cabsl.c
+	s_scalbnl.c s_sinl.c s_sincosl.c \
+	s_tanhl.c s_tanl.c s_truncl.c w_cabsl.c
 .endif
 
 # C99 complex functions
@@ -137,7 +139,8 @@ MAN=	acos.3 acosh.3 asin.3 asinh.3 atan.
 	fma.3 fmax.3 fmod.3 hypot.3 ieee.3 ieee_test.3 ilogb.3 j0.3 \
 	lgamma.3 log.3 lrint.3 lround.3 math.3 nan.3 \
 	nextafter.3 remainder.3 rint.3 \
-	round.3 scalbn.3 signbit.3 sin.3 sinh.3 sqrt.3 tan.3 tanh.3 trunc.3 \
+	round.3 scalbn.3 signbit.3 sin.3 sincos.3 \
+	sinh.3 sqrt.3 tan.3 tanh.3 trunc.3 \
 	complex.3
 
 MLINKS+=acos.3 acosf.3 acos.3 acosl.3
@@ -215,6 +218,7 @@ MLINKS+=round.3 roundf.3 round.3 roundl.
 MLINKS+=scalbn.3 scalbln.3 scalbn.3 scalblnf.3 scalbn.3 scalblnl.3
 MLINKS+=scalbn.3 scalbnf.3 scalbn.3 scalbnl.3
 MLINKS+=sin.3 sinf.3 sin.3 sinl.3
+MLINKS+=sincos.3 sincosf.3 sin.3 sincosl.3
 MLINKS+=sinh.3 sinhf.3 sinh.3 sinhl.3
 MLINKS+=sqrt.3 cbrt.3 sqrt.3 cbrtf.3 sqrt.3 cbrtl.3 sqrt.3 sqrtf.3 \
 	sqrt.3 sqrtl.3

Modified: head/lib/msun/Symbol.map
==============================================================================
--- head/lib/msun/Symbol.map	Sun May 28 05:31:18 2017	(r319046)
+++ head/lib/msun/Symbol.map	Sun May 28 06:13:38 2017	(r319047)
@@ -294,4 +294,7 @@ FBSD_1.5 {
 	casinl;
 	catanl;
 	catanhl;
+	sincos;
+	sincosf;
+	sincosl;
 };

Added: head/lib/msun/man/sincos.3
==============================================================================
--- /dev/null	00:00:00 1970	(empty, because file is newly added)
+++ head/lib/msun/man/sincos.3	Sun May 28 06:13:38 2017	(r319047)
@@ -0,0 +1,82 @@
+.\" Copyright (c) 2011 Steven G. Kargl.
+.\"
+.\" Redistribution and use in source and binary forms, with or without
+.\" modification, are permitted provided that the following conditions
+.\" are met:
+.\" 1. Redistributions of source code must retain the above copyright
+.\"    notice, this list of conditions and the following disclaimer.
+.\" 2. Redistributions in binary form must reproduce the above copyright
+.\"    notice, this list of conditions and the following disclaimer in the
+.\"    documentation and/or other materials provided with the distribution.
+.\"
+.\" THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
+.\" ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+.\" IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
+.\" ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
+.\" FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+.\" DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
+.\" OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
+.\" HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+.\" LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
+.\" OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+.\" SUCH DAMAGE.
+.\"
+.\" $FreeBSD$
+.\"
+.Dd March 12, 2011
+.Dt SINCOS 3
+.Os
+.Sh NAME
+.Nm sincos ,
+.Nm sincosf ,
+.Nm sincosl
+.Nd sine and cosine functions
+.Sh LIBRARY
+.Lb libm
+.Sh SYNOPSIS
+.In math.h
+.Ft void
+.Fn sincos "double x" "double *s" "double *c"
+.Ft void
+.Fn sincosf "float x" "float *s" "float *c"
+.Ft void
+.Fn sincosl "long double x" "long double *s" "long double *c"
+.Sh DESCRIPTION
+The
+.Fn sincos ,
+.Fn sincosf ,
+and
+.Fn sincosl
+functions compute the sine and cosine of
+.Fa x .
+Using these functions allows argument reduction to occur only
+once instead of twice with individual invocations of 
+.Fn sin
+and 
+.Fn cos .
+Like 
+.Fn sin
+and 
+.Fn cos ,
+a large magnitude argument may yield a result with little
+or no significance.
+.Sh RETURN VALUES
+Upon returning from 
+.Fn sincos ,
+.Fn sincosf ,
+and
+.Fn sincosl ,
+the memory pointed to by 
+.Ar "*s" 
+and
+.Ar "*c" 
+are assigned the values of sine and cosine, respectively.
+.Sh SEE ALSO
+.Xr cos 3 ,
+.Xr sin 3 ,
+.Sh HISTORY
+These functions were added to 
+.Fx 9.0
+to aid in writing various complex function contained in 
+.St -isoC-99 .
+

Added: head/lib/msun/src/k_sincos.h
==============================================================================
--- /dev/null	00:00:00 1970	(empty, because file is newly added)
+++ head/lib/msun/src/k_sincos.h	Sun May 28 06:13:38 2017	(r319047)
@@ -0,0 +1,52 @@
+/*-
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice 
+ * is preserved.
+ * ====================================================
+ * 
+ * k_sin.c and k_cos.c merged by Steven G. Kargl.
+ */
+
+#include <sys/cdefs.h>
+__FBSDID("$FreeBSD$");
+
+static const double
+S1  = -1.66666666666666324348e-01, /* 0xBFC55555, 0x55555549 */
+S2  =  8.33333333332248946124e-03, /* 0x3F811111, 0x1110F8A6 */
+S3  = -1.98412698298579493134e-04, /* 0xBF2A01A0, 0x19C161D5 */
+S4  =  2.75573137070700676789e-06, /* 0x3EC71DE3, 0x57B1FE7D */
+S5  = -2.50507602534068634195e-08, /* 0xBE5AE5E6, 0x8A2B9CEB */
+S6  =  1.58969099521155010221e-10; /* 0x3DE5D93A, 0x5ACFD57C */
+
+static const double
+C1  =  4.16666666666666019037e-02, /* 0x3FA55555, 0x5555554C */
+C2  = -1.38888888888741095749e-03, /* 0xBF56C16C, 0x16C15177 */
+C3  =  2.48015872894767294178e-05, /* 0x3EFA01A0, 0x19CB1590 */
+C4  = -2.75573143513906633035e-07, /* 0xBE927E4F, 0x809C52AD */
+C5  =  2.08757232129817482790e-09, /* 0x3E21EE9E, 0xBDB4B1C4 */
+C6  = -1.13596475577881948265e-11; /* 0xBDA8FAE9, 0xBE8838D4 */
+
+static inline void
+__kernel_sincos(double x, double y, int iy, double *sn, double *cs)
+{
+	double hz, r, v, w, z;
+
+	z = x * x;
+	w = z * z;
+	r = S2 + z * (S3 + z * S4) + z * w * (S5 + z * S6);
+	v = z * x;
+
+	if (iy == 0)
+		*sn = x + v * (S1 + z * r);
+	else
+		*sn = x - ((z * (y / 2 - v * r) - y) - v * S1);
+
+	r = z * (C1 + z * (C2 + z * C3)) + w * w * (C4 + z * (C5 + z * C6));
+	hz = z / 2;
+	w = 1 - hz;
+	*cs = w + (((1 - w) - hz) + (z * r - x * y));
+}

Added: head/lib/msun/src/k_sincosf.h
==============================================================================
--- /dev/null	00:00:00 1970	(empty, because file is newly added)
+++ head/lib/msun/src/k_sincosf.h	Sun May 28 06:13:38 2017	(r319047)
@@ -0,0 +1,43 @@
+/*-
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ *
+ * k_sinf.c and k_cosf.c merged by Steven G. Kargl.
+ */
+
+#include <sys/cdefs.h>
+__FBSDID("$FreeBSD$");
+
+/* |sin(x)/x - s(x)| < 2**-37.5 (~[-4.89e-12, 4.824e-12]). */
+static const double
+S1 = -0x15555554cbac77.0p-55,	/* -0.166666666416265235595 */
+S2 =  0x111110896efbb2.0p-59,	/*  0.0083333293858894631756 */
+S3 = -0x1a00f9e2cae774.0p-65,	/* -0.000198393348360966317347 */
+S4 =  0x16cd878c3b46a7.0p-71;	/*  0.0000027183114939898219064 */
+
+/* |cos(x) - c(x)| < 2**-34.1 (~[-5.37e-11, 5.295e-11]). */
+static const double
+C0  = -0x1ffffffd0c5e81.0p-54,	/* -0.499999997251031003120 */
+C1  =  0x155553e1053a42.0p-57,	/*  0.0416666233237390631894 */
+C2  = -0x16c087e80f1e27.0p-62,	/* -0.00138867637746099294692 */
+C3  =  0x199342e0ee5069.0p-68;	/*  0.0000243904487962774090654 */
+
+static inline void
+__kernel_sincosdf(double x, float *sn, float *cs)
+{
+	double r, s, w, z;
+
+	z = x * x;
+	w = z * z;
+	r = S3 + z * S4;
+	s = z * x;
+	*sn = (x + s * (S1 + z * S2)) + s * w * r;
+	r = C2 + z * C3;
+	*cs = ((1 + z * C0) + w * C1) + (w * z) * r;
+}

Added: head/lib/msun/src/k_sincosl.h
==============================================================================
--- /dev/null	00:00:00 1970	(empty, because file is newly added)
+++ head/lib/msun/src/k_sincosl.h	Sun May 28 06:13:38 2017	(r319047)
@@ -0,0 +1,134 @@
+/*-
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ * Copyright (c) 2008 Steven G. Kargl, David Schultz, Bruce D. Evans.
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice 
+ * is preserved.
+ * ====================================================
+ *
+ * k_sinl.c and k_cosl.c merged by Steven G. Kargl
+ */
+
+#include <sys/cdefs.h>
+__FBSDID("$FreeBSD$");
+
+#if LDBL_MANT_DIG == 64		/* ld80 version of k_sincosl.c. */
+
+#if defined(__amd64__) || defined(__i386__)
+/* Long double constants are slow on these arches, and broken on i386. */
+static const volatile double
+C1hi = 0.041666666666666664,		/*  0x15555555555555.0p-57 */
+C1lo = 2.2598839032744733e-18,		/*  0x14d80000000000.0p-111 */
+S1hi = -0.16666666666666666,		/* -0x15555555555555.0p-55 */
+S1lo = -9.2563760475949941e-18;		/* -0x15580000000000.0p-109 */
+#define	S1	((long double)S1hi + S1lo)
+#define	C1	((long double)C1hi + C1lo)
+#else
+static const long double
+C1 =  0.0416666666666666666136L;	/*  0xaaaaaaaaaaaaaa9b.0p-68 */
+S1 = -0.166666666666666666671L,		/* -0xaaaaaaaaaaaaaaab.0p-66 */
+#endif
+
+static const double
+C2 = -0.0013888888888888874,		/* -0x16c16c16c16c10.0p-62 */
+C3 =  0.000024801587301571716,		/*  0x1a01a01a018e22.0p-68 */
+C4 = -0.00000027557319215507120,	/* -0x127e4fb7602f22.0p-74 */
+C5 =  0.0000000020876754400407278,	/*  0x11eed8caaeccf1.0p-81 */
+C6 = -1.1470297442401303e-11,		/* -0x19393412bd1529.0p-89 */
+C7 =  4.7383039476436467e-14,		/*  0x1aac9d9af5c43e.0p-97 */
+S2 =  0.0083333333333333332,		/*  0x11111111111111.0p-59 */
+S3 = -0.00019841269841269427,		/* -0x1a01a01a019f81.0p-65 */
+S4 =  0.0000027557319223597490,		/*  0x171de3a55560f7.0p-71 */
+S5 = -0.000000025052108218074604,	/* -0x1ae64564f16cad.0p-78 */
+S6 =  1.6059006598854211e-10,		/*  0x161242b90243b5.0p-85 */
+S7 = -7.6429779983024564e-13,		/* -0x1ae42ebd1b2e00.0p-93 */
+S8 =  2.6174587166648325e-15;		/*  0x179372ea0b3f64.0p-101 */
+
+static inline void
+__kernel_sincosl(long double x, long double y, int iy, long double *sn,
+    long double *cs)
+{
+	long double hz, r, v, w, z;
+
+	z = x * x;
+	v = z * x;
+	/*
+	 * XXX Replace Horner scheme with an algorithm suitable for CPUs
+	 * with more complex pipelines.
+	 */
+	r = S2 + z * (S3 + z * (S4 + z * (S5 + z * (S6 + z * (S7 + z * S8)))));
+
+	if (iy == 0)
+		*sn = x + v * (S1 + z * r);
+	else
+		*sn = x - ((z * (y / 2 - v * r) - y) - v * S1);
+
+	hz = z / 2;
+	w = 1 - hz;
+	r = z * (C1 + z * (C2 + z * (C3 + z * (C4 + z * (C5 + z * (C6 +
+	    z * C7))))));
+	*cs = w + (((1 - w) - hz) + (z * r - x * y));
+}
+
+#elif LDBL_MANT_DIG == 113	/* ld128 version of k_sincosl.c. */
+
+static const long double
+C1 =  0.04166666666666666666666666666666658424671L,
+C2 = -0.001388888888888888888888888888863490893732L,
+C3 =  0.00002480158730158730158730158600795304914210L,
+C4 = -0.2755731922398589065255474947078934284324e-6L,
+C5 =  0.2087675698786809897659225313136400793948e-8L,
+C6 = -0.1147074559772972315817149986812031204775e-10L,
+C7 =  0.4779477332386808976875457937252120293400e-13L,
+S1 = -0.16666666666666666666666666666666666606732416116558L,
+S2 =  0.0083333333333333333333333333333331135404851288270047L,
+S3 = -0.00019841269841269841269841269839935785325638310428717L,
+S4 =  0.27557319223985890652557316053039946268333231205686e-5L,
+S5 = -0.25052108385441718775048214826384312253862930064745e-7L,
+S6 =  0.16059043836821614596571832194524392581082444805729e-9L,
+S7 = -0.76471637318198151807063387954939213287488216303768e-12L,
+S8 =  0.28114572543451292625024967174638477283187397621303e-14L;
+
+static const double
+C8  = -0.1561920696721507929516718307820958119868e-15,
+C9  =  0.4110317413744594971475941557607804508039e-18,
+C10 = -0.8896592467191938803288521958313920156409e-21,
+C11 =  0.1601061435794535138244346256065192782581e-23,
+S9  = -0.82206352458348947812512122163446202498005154296863e-17,
+S10 =  0.19572940011906109418080609928334380560135358385256e-19,
+S11 = -0.38680813379701966970673724299207480965452616911420e-22,
+S12 =  0.64038150078671872796678569586315881020659912139412e-25;
+
+static inline void
+__kernel_sincosl(long double x, long double y, int iy, long double *sn, 
+    long double *cs)
+{
+	long double hz, r, v, w, z;
+
+	z = x * x;
+	v = z * x;
+	/*
+	 * XXX Replace Horner scheme with an algorithm suitable for CPUs
+	 * with more complex pipelines.
+	 */
+	r = S2 + z * (S3 + z * (S4 + z * (S5 + z * (S6 + z * (S7 + z * (S8 +
+	    z * (S9 + z * (S10 + z * (S11 + z * S12)))))))));
+
+	if (iy == 0)
+		*sn = x + v * (S1 + z * r);
+	else
+		*cs = x - ((z * (y / 2 - v * r) - y) - v * S1);
+
+	hz = z / 2;
+	w = 1 - hz;
+	r = z * (C1 + z * (C2 + z * (C3 + z * (C4 + z * (C5 + z * (C6 + 
+	    z * (C7 + z * (C8 + z * (C9 + z * (C10 + z * C11))))))))));
+
+	*cs =  w + (((1 - w) - hz) + (z * r - x * y));
+}
+#else
+#error "Unsupported long double format"
+#endif

Modified: head/lib/msun/src/math.h
==============================================================================
--- head/lib/msun/src/math.h	Sun May 28 05:31:18 2017	(r319046)
+++ head/lib/msun/src/math.h	Sun May 28 06:13:38 2017	(r319047)
@@ -500,6 +500,9 @@ long double	truncl(long double);
 
 #if __BSD_VISIBLE
 long double	lgammal_r(long double, int *);
+void		sincos(double, double *, double *);
+void		sincosf(float, float *, float *);
+void		sincosl(long double, long double *, long double *);
 #endif
 
 __END_DECLS

Modified: head/lib/msun/src/math_private.h
==============================================================================
--- head/lib/msun/src/math_private.h	Sun May 28 05:31:18 2017	(r319046)
+++ head/lib/msun/src/math_private.h	Sun May 28 06:13:38 2017	(r319047)
@@ -306,9 +306,21 @@ do {								\
 		fpsetprec(__oprec);		\
 	RETURNF(__retval);			\
 } while (0)
+#define	ENTERV()				\
+	fp_prec_t __oprec;			\
+						\
+	if ((__oprec = fpgetprec()) != FP_PE)	\
+		fpsetprec(FP_PE)
+#define	RETURNV() do {				\
+	if (__oprec != FP_PE)			\
+		fpsetprec(__oprec);		\
+	return;			\
+} while (0)
 #else
-#define	ENTERI(x)
+#define	ENTERI()
 #define	RETURNI(x)	RETURNF(x)
+#define	ENTERV()
+#define	RETURNV()	return
 #endif
 
 /* Default return statement if hack*_t() is not used. */

Added: head/lib/msun/src/s_sincos.c
==============================================================================
--- /dev/null	00:00:00 1970	(empty, because file is newly added)
+++ head/lib/msun/src/s_sincos.c	Sun May 28 06:13:38 2017	(r319047)
@@ -0,0 +1,80 @@
+/*-
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ * 
+ * s_sin.c and s_cos.c merged by Steven G. Kargl.  Descriptions of the
+ * algorithms are contained in the original files.
+ */
+
+#include <sys/cdefs.h>
+__FBSDID("$FreeBSD$");
+
+#include <float.h>
+
+#include "math.h"
+#define INLINE_REM_PIO2
+#include "math_private.h"
+#include "e_rem_pio2.c"
+#include "k_sincos.h"
+
+void
+sincos(double x, double *sn, double *cs)
+{
+	double y[2];
+	int32_t n, ix;
+
+	/* High word of x. */
+	GET_HIGH_WORD(ix, x);
+
+	/* |x| ~< pi/4 */
+	ix &= 0x7fffffff;
+	if (ix <= 0x3fe921fb) {
+		if (ix < 0x3e400000) {		/* |x| < 2**-27 */
+			if ((int)x == 0) {	/* Generate inexact. */
+				*sn = x;
+				*cs = 1;
+				return;
+			}
+		}
+		__kernel_sincos(x, 0, 0, sn, cs);
+		return;
+	}
+
+	/* If x = Inf or NaN, then sin(x) = NaN and cos(x) = NaN. */
+	if (ix >= 0x7ff00000) {
+		*sn = x - x;
+		*cs = x - x;
+		return;
+	}
+
+	/* Argument reduction. */
+	n = __ieee754_rem_pio2(x, y);
+
+	switch(n & 3) {
+	case 0:
+		__kernel_sincos(y[0], y[1], 1, sn, cs);
+		break;
+	case 1:
+		__kernel_sincos(y[0], y[1], 1, cs, sn);
+		*cs = -*cs;
+		break;
+	case 2:
+		__kernel_sincos(y[0], y[1], 1, sn, cs);
+		*sn = -*sn;
+		*cs = -*cs;
+		break;
+	default:
+		__kernel_sincos(y[0], y[1], 1, cs, sn);
+		*sn = -*sn;
+	}
+}
+
+#if (LDBL_MANT_DIG == 53)
+__weak_reference(sincos, sincosl);
+#endif

Added: head/lib/msun/src/s_sincosf.c
==============================================================================
--- /dev/null	00:00:00 1970	(empty, because file is newly added)
+++ head/lib/msun/src/s_sincosf.c	Sun May 28 06:13:38 2017	(r319047)
@@ -0,0 +1,126 @@
+/*-
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+/* s_sincosf.c -- float version of s_sincos.c.
+ * Conversion to float by Ian Lance Taylor, Cygnus Support, ian at cygnus.com.
+ * Optimized by Bruce D. Evans.
+ * Merged s_sinf.c and s_cosf.c by Steven G. Kargl.
+ */
+
+#include <sys/cdefs.h>
+__FBSDID("$FreeBSD$");
+
+#include <float.h>
+
+#include "math.h"
+#define INLINE_REM_PIO2F
+#include "math_private.h"
+#include "e_rem_pio2f.c"
+#include "k_sincosf.h"
+
+/* Small multiples of pi/2 rounded to double precision. */
+static const double
+p1pio2 = 1*M_PI_2,			/* 0x3FF921FB, 0x54442D18 */
+p2pio2 = 2*M_PI_2,			/* 0x400921FB, 0x54442D18 */
+p3pio2 = 3*M_PI_2,			/* 0x4012D97C, 0x7F3321D2 */
+p4pio2 = 4*M_PI_2;			/* 0x401921FB, 0x54442D18 */
+
+void
+sincosf(float x, float *sn, float *cs)
+{
+	float c, s;
+	double y;
+	int32_t n, hx, ix;
+
+	GET_FLOAT_WORD(hx, x);
+	ix = hx & 0x7fffffff;
+
+	if (ix <= 0x3f490fda) {		/* |x| ~<= pi/4 */
+		if (ix < 0x39800000) {	/* |x| < 2**-12 */
+			if ((int)x == 0) {
+				*sn = x;	/* x with inexact if x != 0 */
+				*cs = 1;
+				return;
+			}
+		}
+		__kernel_sincosdf(x, sn, cs);
+		return;
+	}
+
+	if (ix <= 0x407b53d1) {		/* |x| ~<= 5*pi/4 */
+		if (ix <= 0x4016cbe3) {	/* |x| ~<= 3pi/4 */
+			if (hx > 0) {
+				__kernel_sincosdf(x - p1pio2, cs, sn);
+				*cs = -*cs;
+			} else {
+				__kernel_sincosdf(x + p1pio2, cs, sn);
+				*sn = -*sn;
+			}
+		} else {
+			if (hx > 0)
+				__kernel_sincosdf(x - p2pio2, sn, cs);
+			else
+				__kernel_sincosdf(x + p2pio2, sn, cs);
+			*sn = -*sn;
+			*cs = -*cs;
+		}
+		return;
+	}
+
+	if (ix <= 0x40e231d5) {		/* |x| ~<= 9*pi/4 */
+		if (ix <= 0x40afeddf) {	/* |x| ~<= 7*pi/4 */
+			if (hx > 0) {
+				__kernel_sincosdf(x - p3pio2, cs, sn);
+				*sn = -*sn;
+			} else {
+				__kernel_sincosdf(x + p3pio2, cs, sn);
+				*cs = -*cs;
+			}
+		} else {
+			if (hx > 0)
+				__kernel_sincosdf(x - p4pio2, sn, cs);
+			else
+				__kernel_sincosdf(x + p4pio2, sn, cs);
+		}
+		return;
+	}
+
+	/* If x = Inf or NaN, then sin(x) = NaN and cos(x) = NaN. */
+	if (ix >= 0x7f800000) {
+		*sn = x - x;
+		*cs = x - x;
+		return;
+	}
+
+	/* Argument reduction. */
+	n = __ieee754_rem_pio2f(x, &y);
+	__kernel_sincosdf(y, &s, &c);
+
+	switch(n & 3) {
+	case 0:
+		*sn = s;
+		*cs = c;
+		break;
+	case 1:
+		*sn = c;
+		*cs = -s;
+		break;
+	case 2:
+		*sn = -s;
+		*cs = -c;
+		break;
+	default:
+		*sn = -c;
+		*cs = s;
+	}
+}
+
+

Added: head/lib/msun/src/s_sincosl.c
==============================================================================
--- /dev/null	00:00:00 1970	(empty, because file is newly added)
+++ head/lib/msun/src/s_sincosl.c	Sun May 28 06:13:38 2017	(r319047)
@@ -0,0 +1,105 @@
+/*-
+ * Copyright (c) 2007, 2010-2013 Steven G. Kargl
+ * All rights reserved.
+ *
+ * Redistribution and use in source and binary forms, with or without
+ * modification, are permitted provided that the following conditions
+ * are met:
+ * 1. Redistributions of source code must retain the above copyright
+ *    notice unmodified, this list of conditions, and the following
+ *    disclaimer.
+ * 2. Redistributions in binary form must reproduce the above copyright
+ *    notice, this list of conditions and the following disclaimer in the
+ *    documentation and/or other materials provided with the distribution.
+ *
+ * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
+ * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
+ * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
+ * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT,
+ * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
+ * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
+ * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
+ * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
+ * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
+ * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ *
+ * s_sinl.c and s_cosl.c merged by Steven G. Kargl.
+ */
+
+#include <sys/cdefs.h>
+__FBSDID("$FreeBSD$");
+
+#include <float.h>
+#ifdef __i386__
+#include <ieeefp.h>
+#endif
+
+#include "math.h"
+#include "math_private.h"
+#include "k_sincosl.h"
+
+#if LDBL_MANT_DIG == 64
+#include "../ld80/e_rem_pio2l.h"
+#elif LDBL_MANT_DIG == 113
+#include "../ld128/e_rem_pio2l.h"
+#else
+#error "Unsupported long double format"
+#endif
+
+void
+sincosl(long double x, long double *sn, long double *cs)
+{
+	union IEEEl2bits z;
+	int e0, sgn;
+	long double y[2];
+
+	z.e = x;
+	sgn = z.bits.sign;
+	z.bits.sign = 0;
+
+	ENTERV();
+
+	/* Optimize the case where x is already within range. */
+	if (z.e < M_PI_4) {
+		/*
+		 * If x = +-0 or x is a subnormal number, then sin(x) = x and
+		 * cos(x) = 1.
+		 */
+		if (z.bits.exp == 0) {
+			*sn = x;
+			*cs = 1;
+		} else
+			__kernel_sincosl(x, 0, 0, sn, cs);
+		RETURNV();
+	}
+
+	/* If x = NaN or Inf, then sin(x) and cos(x) are NaN. */
+	if (z.bits.exp == 32767) {
+		*sn = x - x;
+		*cs = x - x;
+		RETURNV();
+	}
+
+	/* Range reduction. */
+	e0 = __ieee754_rem_pio2l(x, y);
+
+	switch (e0 & 3) {
+	case 0:
+		__kernel_sincosl(y[0], y[1], 1, sn, cs);
+		break;
+	case 1:
+		__kernel_sincosl(y[0], y[1], 1, cs, sn);
+		*cs = -*cs;
+		break;
+	case 2:
+		__kernel_sincosl(y[0], y[1], 1, sn, cs);
+		*sn = -*sn;
+		*cs = -*cs;
+		break;
+	default:
+		__kernel_sincosl(y[0], y[1], 1, cs, sn);
+		*sn = -*sn;
+	}
+
+	RETURNV();
+}


More information about the svn-src-all mailing list