1/*-
2 * Copyright (c) 2017,2023 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
27#include <sys/cdefs.h>
28
29/*
30 * See ../src/s_sinpi.c for implementation details.
31 */
32
33#define	INLINE_KERNEL_SINDF
34#define	INLINE_KERNEL_COSDF
35
36#include "namespace.h"
37__weak_alias(sinpif, _sinpif)
38
39#include "math.h"
40#include "math_private.h"
41#include "k_cosdf.c"
42#include "k_sindf.c"
43
44#define	__kernel_cospif(x)	(__kernel_cosdf(M_PI * (x)))
45#define	__kernel_sinpif(x)	(__kernel_sindf(M_PI * (x)))
46
47static const float
48pi_hi =  3.14160156e+00F,	/* 0x40491000 */
49pi_lo = -8.90890988e-06F;	/* 0xb715777a */
50
51static volatile const float vzero = 0;
52
53float
54sinpif(float x)
55{
56	float ax, hi, lo, s;
57	uint32_t hx, ix, j0;
58
59	GET_FLOAT_WORD(hx, x);
60	ix = hx & 0x7fffffff;
61	SET_FLOAT_WORD(ax, ix);
62
63	if (ix < 0x3f800000) {			/* |x| < 1 */
64		if (ix < 0x3e800000) {		/* |x| < 0.25 */
65	 		if (ix < 0x38800000) {	/* |x| < 0x1p-14 */
66				if (x == 0)
67					return (x);
68				SET_FLOAT_WORD(hi, hx & 0xffff0000);
69				hi *= 0x1p23F;
70				lo = x * 0x1p23F - hi;
71				s = (pi_lo + pi_hi) * lo + pi_lo * hi +
72				    pi_hi * hi;
73				return (s * 0x1p-23F);
74			}
75
76			s = __kernel_sinpif(ax);
77			return ((hx & 0x80000000) ? -s : s);
78		}
79
80		if (ix < 0x3f000000)		/* |x| < 0.5 */
81			s = __kernel_cospif(0.5F - ax);
82		else if (ix < 0x3f400000)	/* |x| < 0.75 */
83			s = __kernel_cospif(ax - 0.5F);
84		else
85			s = __kernel_sinpif(1 - ax);
86		return ((hx & 0x80000000) ? -s : s);
87	}
88
89	if (ix < 0x4b000000) {		/* 1 <= |x| < 0x1p23 */
90		FFLOORF(x, j0, ix);	/* Integer part of ax. */
91		ax -= x;
92		GET_FLOAT_WORD(ix, ax);
93
94		if (ix == 0)
95			s = 0;
96		else {
97			if (ix < 0x3f000000) {		/* |x| < 0.5 */
98				if (ix < 0x3e800000)	/* |x| < 0.25 */
99					s = __kernel_sinpif(ax);
100				else
101					s = __kernel_cospif(0.5F - ax);
102			} else {
103				if (ix < 0x3f400000)	/* |x| < 0.75 */
104					s = __kernel_cospif(ax - 0.5F);
105				else
106					s = __kernel_sinpif(1 - ax);
107			}
108
109			j0 = (uint32_t)x;
110			s = (j0 & 1) ? -s : s;
111		}
112		return ((hx & 0x80000000) ? -s : s);
113	}
114
115	/* x = +-inf or nan. */
116	if (ix >= 0x7f800000)
117		return (vzero / vzero);
118
119	/*
120	 * |x| >= 0x1p23 is always an integer, so return +-0.
121	 */
122	return (copysignf(0, x));
123}
124