1 /*-
2  * Copyright (c) 2008-2010 David Schultz <das@FreeBSD.org>
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, this list of conditions and the following disclaimer.
10  * 2. Redistributions in binary form must reproduce the above copyright
11  *    notice, this list of conditions and the following disclaimer in the
12  *    documentation and/or other materials provided with the distribution.
13  *
14  * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
15  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
16  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
17  * ARE DISCLAIMED.  IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
18  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
19  * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
20  * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
21  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
22  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
23  * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
24  * SUCH DAMAGE.
25  *
26  * $FreeBSD: src/tools/regression/lib/msun/test-logarithm.c,v 1.3 2011/10/15 05:28:13 das Exp $
27  */
28 
29 /*
30  * Tests for corner cases in log*().
31  */
32 
33 #include <assert.h>
34 #include <fenv.h>
35 #include <float.h>
36 #include <math.h>
37 #include <stdio.h>
38 
39 #ifdef __i386__
40 #include <ieeefp.h>
41 #endif
42 
43 #define	ALL_STD_EXCEPT	(FE_DIVBYZERO | FE_INEXACT | FE_INVALID | \
44 			 FE_OVERFLOW | FE_UNDERFLOW)
45 
46 #pragma STDC FENV_ACCESS ON
47 
48 /*
49  * Test that a function returns the correct value and sets the
50  * exception flags correctly. The exceptmask specifies which
51  * exceptions we should check. We need to be lenient for several
52  * reasoons, but mainly because on some architectures it's impossible
53  * to raise FE_OVERFLOW without raising FE_INEXACT.
54  *
55  * These are macros instead of functions so that assert provides more
56  * meaningful error messages.
57  *
58  * XXX The volatile here is to avoid gcc's bogus constant folding and work
59  *     around the lack of support for the FENV_ACCESS pragma.
60  */
61 #define	test(func, x, result, exceptmask, excepts)	do {		\
62 	volatile long double _d = x;					\
63 	assert(feclearexcept(FE_ALL_EXCEPT) == 0);			\
64 	assert(fpequal((func)(_d), (result)));				 \
65 	assert(((func), fetestexcept(exceptmask) == (excepts)));	\
66 } while (0)
67 
68 /* Test all the functions that compute log(x). */
69 #define	testall0(x, result, exceptmask, excepts)	do {		\
70 	test(log, x, result, exceptmask, excepts);			\
71 	test(logf, x, result, exceptmask, excepts);			\
72 	test(log2, x, result, exceptmask, excepts);			\
73 	test(log2f, x, result, exceptmask, excepts);			\
74 	test(log10, x, result, exceptmask, excepts);			\
75 	test(log10f, x, result, exceptmask, excepts);			\
76 } while (0)
77 
78 /* Test all the functions that compute log(1+x). */
79 #define	testall1(x, result, exceptmask, excepts)	do {		\
80 	test(log1p, x, result, exceptmask, excepts);			\
81 	test(log1pf, x, result, exceptmask, excepts);			\
82 } while (0)
83 
84 /*
85  * Determine whether x and y are equal, with two special rules:
86  *	+0.0 != -0.0
87  *	 NaN == NaN
88  */
89 int
90 fpequal(long double x, long double y)
91 {
92 	return ((x == y && !signbit(x) == !signbit(y)) || isnan(x) && isnan(y));
93 }
94 
95 void
96 run_generic_tests(void)
97 {
98 
99 	/* log(1) == 0, no exceptions raised */
100 	testall0(1.0, 0.0, ALL_STD_EXCEPT, 0);
101 	testall1(0.0, 0.0, ALL_STD_EXCEPT, 0);
102 	testall1(-0.0, -0.0, ALL_STD_EXCEPT, 0);
103 
104 	/* log(NaN) == NaN, no exceptions raised */
105 	testall0(NAN, NAN, ALL_STD_EXCEPT, 0);
106 	testall1(NAN, NAN, ALL_STD_EXCEPT, 0);
107 
108 	/* log(Inf) == Inf, no exceptions raised */
109 	testall0(INFINITY, INFINITY, ALL_STD_EXCEPT, 0);
110 	testall1(INFINITY, INFINITY, ALL_STD_EXCEPT, 0);
111 
112 	/* log(x) == NaN for x < 0, invalid exception raised */
113 	testall0(-INFINITY, NAN, ALL_STD_EXCEPT, FE_INVALID);
114 	testall1(-INFINITY, NAN, ALL_STD_EXCEPT, FE_INVALID);
115 	testall0(-1.0, NAN, ALL_STD_EXCEPT, FE_INVALID);
116 	testall1(-1.5, NAN, ALL_STD_EXCEPT, FE_INVALID);
117 
118 	/* log(0) == -Inf, divide-by-zero exception */
119 	testall0(0.0, -INFINITY, ALL_STD_EXCEPT & ~FE_INEXACT, FE_DIVBYZERO);
120 	testall0(-0.0, -INFINITY, ALL_STD_EXCEPT & ~FE_INEXACT, FE_DIVBYZERO);
121 	testall1(-1.0, -INFINITY, ALL_STD_EXCEPT & ~FE_INEXACT, FE_DIVBYZERO);
122 }
123 
124 void
125 run_log2_tests(void)
126 {
127 	int i;
128 
129 	/*
130 	 * We should insist that log2() return exactly the correct
131 	 * result and not raise an inexact exception for powers of 2.
132 	 */
133 	feclearexcept(FE_ALL_EXCEPT);
134 	for (i = FLT_MIN_EXP - FLT_MANT_DIG; i < FLT_MAX_EXP; i++) {
135 		assert(log2f(ldexpf(1.0, i)) == i);
136 		assert(fetestexcept(ALL_STD_EXCEPT) == 0);
137 	}
138 	for (i = DBL_MIN_EXP - DBL_MANT_DIG; i < DBL_MAX_EXP; i++) {
139 		assert(log2(ldexp(1.0, i)) == i);
140 		assert(fetestexcept(ALL_STD_EXCEPT) == 0);
141 	}
142 }
143 
144 void
145 run_roundingmode_tests(void)
146 {
147 
148 	/*
149 	 * Corner cases in other rounding modes.
150 	 */
151 	fesetround(FE_DOWNWARD);
152 	/* These are still positive per IEEE 754R */
153 	testall0(1.0, 0.0, ALL_STD_EXCEPT, 0);
154 	testall1(0.0, 0.0, ALL_STD_EXCEPT, 0);
155 	fesetround(FE_TOWARDZERO);
156 	testall0(1.0, 0.0, ALL_STD_EXCEPT, 0);
157 	testall1(0.0, 0.0, ALL_STD_EXCEPT, 0);
158 
159 	fesetround(FE_UPWARD);
160 	testall0(1.0, 0.0, ALL_STD_EXCEPT, 0);
161 	testall1(0.0, 0.0, ALL_STD_EXCEPT, 0);
162 	/* log1p(-0.0) == -0.0 even when rounding upwards */
163 	testall1(-0.0, -0.0, ALL_STD_EXCEPT, 0);
164 
165 	fesetround(FE_TONEAREST);
166 }
167 
168 int
169 main(int argc, char *argv[])
170 {
171 
172 	printf("1..3\n");
173 
174 	run_generic_tests();
175 	printf("ok 1 - logarithm\n");
176 
177 	run_log2_tests();
178 	printf("ok 2 - logarithm\n");
179 
180 	run_roundingmode_tests();
181 	printf("ok 3 - logarithm\n");
182 
183 	return (0);
184 }
185