1251241Sdas/*-
2251241Sdas * Copyright (c) 2005-2013 David Schultz <das@FreeBSD.org>
3251241Sdas * All rights reserved.
4251241Sdas *
5251241Sdas * Redistribution and use in source and binary forms, with or without
6251241Sdas * modification, are permitted provided that the following conditions
7251241Sdas * are met:
8251241Sdas * 1. Redistributions of source code must retain the above copyright
9251241Sdas *    notice, this list of conditions and the following disclaimer.
10251241Sdas * 2. Redistributions in binary form must reproduce the above copyright
11251241Sdas *    notice, this list of conditions and the following disclaimer in the
12251241Sdas *    documentation and/or other materials provided with the distribution.
13251241Sdas *
14251241Sdas * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
15251241Sdas * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
16251241Sdas * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
17251241Sdas * ARE DISCLAIMED.  IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
18251241Sdas * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
19251241Sdas * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
20251241Sdas * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
21251241Sdas * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
22251241Sdas * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
23251241Sdas * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
24251241Sdas * SUCH DAMAGE.
25251241Sdas *
26251241Sdas * $FreeBSD: stable/11/lib/msun/tests/test-utils.h 315121 2017-03-12 04:52:09Z ngie $
27251241Sdas */
28251241Sdas
29251241Sdas#ifndef	_TEST_UTILS_H_
30251241Sdas#define	_TEST_UTILS_H_
31251241Sdas
32251241Sdas#include <complex.h>
33251241Sdas#include <fenv.h>
34251241Sdas
35251241Sdas/*
36251241Sdas * Implementations are permitted to define additional exception flags
37251241Sdas * not specified in the standard, so it is not necessarily true that
38251241Sdas * FE_ALL_EXCEPT == ALL_STD_EXCEPT.
39251241Sdas */
40251241Sdas#define	ALL_STD_EXCEPT	(FE_DIVBYZERO | FE_INEXACT | FE_INVALID | \
41251241Sdas			 FE_OVERFLOW | FE_UNDERFLOW)
42251241Sdas#define	OPT_INVALID	(ALL_STD_EXCEPT & ~FE_INVALID)
43251241Sdas#define	OPT_INEXACT	(ALL_STD_EXCEPT & ~FE_INEXACT)
44251241Sdas#define	FLT_ULP()	ldexpl(1.0, 1 - FLT_MANT_DIG)
45251241Sdas#define	DBL_ULP()	ldexpl(1.0, 1 - DBL_MANT_DIG)
46251241Sdas#define	LDBL_ULP()	ldexpl(1.0, 1 - LDBL_MANT_DIG)
47251241Sdas
48251241Sdas/*
49251241Sdas * Flags that control the behavior of various fpequal* functions.
50251241Sdas * XXX This is messy due to merging various notions of "close enough"
51251241Sdas * that are best suited for different functions.
52251241Sdas *
53251241Sdas * CS_REAL
54251241Sdas * CS_IMAG
55251241Sdas * CS_BOTH
56251241Sdas *   (cfpequal_cs, fpequal_tol, cfpequal_tol) Whether to check the sign of
57251241Sdas *   the real part of the result, the imaginary part, or both.
58251241Sdas *
59251241Sdas * FPE_ABS_ZERO
60251241Sdas *   (fpequal_tol, cfpequal_tol) If set, treats the tolerance as an absolute
61251241Sdas *   tolerance when the expected value is 0.  This is useful when there is
62251241Sdas *   round-off error in the input, e.g., cos(Pi/2) ~= 0.
63251241Sdas */
64251241Sdas#define	CS_REAL		0x01
65251241Sdas#define	CS_IMAG		0x02
66251241Sdas#define	CS_BOTH		(CS_REAL | CS_IMAG)
67251241Sdas#define	FPE_ABS_ZERO	0x04
68251241Sdas
69251241Sdas#ifdef	DEBUG
70251241Sdas#define	debug(...)	printf(__VA_ARGS__)
71251241Sdas#else
72251241Sdas#define	debug(...)	(void)0
73251241Sdas#endif
74251241Sdas
75251241Sdas/*
76251241Sdas * XXX The ancient version of gcc in the base system doesn't support CMPLXL,
77251241Sdas * but we can fake it most of the time.
78251241Sdas */
79251241Sdas#ifndef CMPLXL
80251241Sdasstatic inline long double complex
81251241SdasCMPLXL(long double x, long double y)
82251241Sdas{
83251241Sdas	long double complex z;
84251241Sdas
85251241Sdas	__real__ z = x;
86251241Sdas	__imag__ z = y;
87251241Sdas	return (z);
88251241Sdas}
89251241Sdas#endif
90251241Sdas
91315121Sngiestatic int	fpequal(long double, long double) __used;
92315121Sngiestatic int	cfpequal(long double complex, long double complex) __used;
93315121Sngiestatic int	cfpequal_cs(long double complex, long double complex,
94315121Sngie		    int) __used;
95315121Sngiestatic int	cfpequal_tol(long double complex, long double complex,
96315121Sngie		    long double, unsigned int) __used;
97315121Sngie
98251241Sdas/*
99251241Sdas * Compare d1 and d2 using special rules: NaN == NaN and +0 != -0.
100251241Sdas * Fail an assertion if they differ.
101251241Sdas */
102251241Sdasstatic int
103251241Sdasfpequal(long double d1, long double d2)
104251241Sdas{
105251241Sdas
106251241Sdas	if (d1 != d2)
107251241Sdas		return (isnan(d1) && isnan(d2));
108251241Sdas	return (copysignl(1.0, d1) == copysignl(1.0, d2));
109251241Sdas}
110251241Sdas
111251241Sdas/*
112251241Sdas * Determine whether x and y are equal, with two special rules:
113251241Sdas *	+0.0 != -0.0
114251241Sdas *	 NaN == NaN
115251241Sdas * If checksign is 0, we compare the absolute values instead.
116251241Sdas */
117251241Sdasstatic int
118251241Sdasfpequal_cs(long double x, long double y, int checksign)
119251241Sdas{
120251241Sdas	if (isnan(x) && isnan(y))
121251241Sdas		return (1);
122251241Sdas	if (checksign)
123251241Sdas		return (x == y && !signbit(x) == !signbit(y));
124251241Sdas	else
125251241Sdas		return (fabsl(x) == fabsl(y));
126251241Sdas}
127251241Sdas
128251241Sdasstatic int
129315121Sngiefpequal_tol(long double x, long double y, long double tol,
130315121Sngie    unsigned int flags)
131251241Sdas{
132251241Sdas	fenv_t env;
133251241Sdas	int ret;
134251241Sdas
135251241Sdas	if (isnan(x) && isnan(y))
136251241Sdas		return (1);
137251241Sdas	if (!signbit(x) != !signbit(y) && (flags & CS_BOTH))
138251241Sdas		return (0);
139251241Sdas	if (x == y)
140251241Sdas		return (1);
141251241Sdas	if (tol == 0)
142251241Sdas		return (0);
143251241Sdas
144251241Sdas	/* Hard case: need to check the tolerance. */
145251241Sdas	feholdexcept(&env);
146251241Sdas	/*
147251241Sdas	 * For our purposes here, if y=0, we interpret tol as an absolute
148251241Sdas	 * tolerance. This is to account for roundoff in the input, e.g.,
149251241Sdas	 * cos(Pi/2) ~= 0.
150251241Sdas	 */
151251241Sdas	if ((flags & FPE_ABS_ZERO) && y == 0.0)
152251241Sdas		ret = fabsl(x - y) <= fabsl(tol);
153251241Sdas	else
154251241Sdas		ret = fabsl(x - y) <= fabsl(y * tol);
155251241Sdas	fesetenv(&env);
156251241Sdas	return (ret);
157251241Sdas}
158251241Sdas
159251241Sdasstatic int
160251241Sdascfpequal(long double complex d1, long double complex d2)
161251241Sdas{
162251241Sdas
163251241Sdas	return (fpequal(creall(d1), creall(d2)) &&
164315121Sngie	    fpequal(cimagl(d1), cimagl(d2)));
165251241Sdas}
166251241Sdas
167251241Sdasstatic int
168251241Sdascfpequal_cs(long double complex x, long double complex y, int checksign)
169251241Sdas{
170251241Sdas	return (fpequal_cs(creal(x), creal(y), checksign)
171251241Sdas		&& fpequal_cs(cimag(x), cimag(y), checksign));
172251241Sdas}
173251241Sdas
174251241Sdasstatic int
175251241Sdascfpequal_tol(long double complex x, long double complex y, long double tol,
176315121Sngie    unsigned int flags)
177251241Sdas{
178251241Sdas	return (fpequal_tol(creal(x), creal(y), tol, flags)
179251241Sdas		&& fpequal_tol(cimag(x), cimag(y), tol, flags));
180251241Sdas}
181251241Sdas
182251241Sdas#endif /* _TEST_UTILS_H_ */
183