t_exp.c revision 274600
1/* $NetBSD: t_exp.c,v 1.7 2014/03/17 11:08:11 martin Exp $ */
2
3/*-
4 * Copyright (c) 2011 The NetBSD Foundation, Inc.
5 * All rights reserved.
6 *
7 * This code is derived from software contributed to The NetBSD Foundation
8 * by Jukka Ruohonen.
9 *
10 * Redistribution and use in source and binary forms, with or without
11 * modification, are permitted provided that the following conditions
12 * are met:
13 * 1. Redistributions of source code must retain the above copyright
14 *    notice, this list of conditions and the following disclaimer.
15 * 2. Redistributions in binary form must reproduce the above copyright
16 *    notice, this list of conditions and the following disclaimer in the
17 *    documentation and/or other materials provided with the distribution.
18 *
19 * THIS SOFTWARE IS PROVIDED BY THE NETBSD FOUNDATION, INC. AND CONTRIBUTORS
20 * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
21 * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
22 * PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR CONTRIBUTORS
23 * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
24 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
25 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
26 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
27 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
28 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
29 * POSSIBILITY OF SUCH DAMAGE.
30 */
31
32#include <atf-c.h>
33#include <math.h>
34#include "t_libm.h"
35
36/* y = exp(x) */
37static const struct {
38	double x;
39	double y;
40	double e;
41} exp_values[] = {
42	{  -10, 0.4539992976248485e-4, 1e-4, },
43	{   -5, 0.6737946999085467e-2, 1e-2, },
44	{   -1, 0.3678794411714423,    1e-1, },
45	{ -0.1, 0.9048374180359595,    1e-1, },
46	{    0, 1.0000000000000000,    1,    },
47	{  0.1, 1.1051709180756477,    1,    },
48	{    1, 2.7182818284590452,    1,    },
49	{    5, 148.41315910257660,    1e2, },
50	{   10, 22026.465794806718,    1e4, },
51};
52
53/*
54 * exp2/exp2f(3)
55 */
56ATF_LIBM_TEST(exp2_is_nan, "Test exp2(x) == NaN")
57{
58#ifdef T_LIBM_NAN
59	T_LIBM_CHECK_NAN(0, exp2, T_LIBM_NAN);
60	T_LIBM_CHECK_NAN(0, exp2f, T_LIBM_NAN);
61#else
62	atf_tc_skip("no NaN on this machine");
63#endif
64}
65
66ATF_LIBM_TEST(exp2_is_plus_zero, "Test exp2(x) == +0.0")
67{
68#ifdef T_LIBM_MINUS_INF
69	T_LIBM_CHECK_PLUS_ZERO(0, exp2, T_LIBM_MINUS_INF);
70	T_LIBM_CHECK_PLUS_ZERO(0, exp2f, T_LIBM_MINUS_INF);
71#else
72	atf_tc_skip("no +/-Inf on this machine");
73#endif
74}
75
76ATF_LIBM_TEST(exp2_powers, "Test exp2(x) is correct for some integer x")
77{
78	static const struct {
79		double	x;
80		double	d_y;
81		double	f_y;
82	} v[] = {
83	    { +0.0,	1.0,	1.0 },
84	    { -0.0,	1.0,	1.0 },
85	    {    1,	0x1p1,	0x1p1 },
86	    {    2,	0x1p2,	0x1p2 },
87	    {  100,	0x1p100,	0x1p100 },
88	    {  125,	0x1p125,	0x1p125 },
89	    {  126,	0x1p126,	0x1p126 },
90#if __DBL_MAX_EXP__ > 129
91	    {  127,	0x1p127,	0x1p127 },
92#endif
93#ifdef T_LIBM_PLUS_INF
94	    {  128,	0x1p128,	T_LIBM_PLUS_INF },
95	    {  129,	0x1p129,	T_LIBM_PLUS_INF },
96	    { 1000,	0x1p1000,	T_LIBM_PLUS_INF },
97	    { 1020,	0x1p1020,	T_LIBM_PLUS_INF },
98	    { 1023,	0x1p1023,	T_LIBM_PLUS_INF },
99	    { 1024,	T_LIBM_PLUS_INF,	T_LIBM_PLUS_INF },
100	    { 1030,	T_LIBM_PLUS_INF,	T_LIBM_PLUS_INF },
101	    { 1050,	T_LIBM_PLUS_INF,	T_LIBM_PLUS_INF },
102	    { 2000,	T_LIBM_PLUS_INF,	T_LIBM_PLUS_INF },
103	    { 16383,	T_LIBM_PLUS_INF,	T_LIBM_PLUS_INF },
104	    { 16384,	T_LIBM_PLUS_INF,	T_LIBM_PLUS_INF },
105	    { 16385,	T_LIBM_PLUS_INF,	T_LIBM_PLUS_INF },
106#endif
107	    {   -1,	0x1p-1,	0x1p-1 },
108	    {   -2,	0x1p-2,	0x1p-2 },
109	    { -100,	0x1p-100,	0x1p-100 },
110	    { -127,	0x1p-127,	0x1p-127 },
111	    { -128,	0x1p-128,	0x1p-128 },
112#if __LDBL_MIN_EXP__ < -129
113	    { -300,	0x1p-300,	0.0},
114	    { -400,	0x1p-400,	0.0},
115	    {-1000,	0x1p-1000,	0.0},
116	    {-1022,	0x1p-1022,	0.0},
117	    /* These should be denormal numbers */
118	    {-1023,	0x1p-1023,	0.0},
119	    {-1024,	0x1p-1024,	0.0},
120	    {-1040,	0x1p-1040,	0.0},
121	    {-1060,	0x1p-1060,	0.0},
122	    /* This is the smallest result gcc will allow */
123	    {-1074,	0x1p-1074,	0.0},
124#endif
125	    {-1075,	0x0,	0.0},
126	    {-1080,	0x0,	0.0},
127	    {-2000,	0x0,	0.0},
128	    {-16382,	0x0,	0.0},
129	    {-16383,	0x0,	0.0},
130	    {-16384,	0x0,	0.0},
131	};
132	unsigned int i;
133
134#if defined(__FreeBSD__) && defined(__i386__)
135	atf_tc_expect_fail("a number of the assertions fail on i386");
136#endif
137
138	for (i = 0; i < __arraycount(v); i++) {
139		T_LIBM_CHECK(i, exp2, v[i].x, v[i].d_y, 0.0);
140		T_LIBM_CHECK(i, exp2f, v[i].x, v[i].f_y, 0.0);
141	}
142}
143
144ATF_LIBM_TEST(exp2_values, "Test exp2(x) is correct for some x")
145{
146	static const struct {
147		double	x;
148		double	y;
149		double	d_eps;
150		double	f_eps;
151	} v[] = {
152#if __DBL_MAX_EXP__ > 128
153	    /* The largest double constant */
154	    { 0x1.fffffffffffffp9,	0x1.ffffffffffd3ap1023,
155		0x1p969,	0.0 },
156	    /* The largest float constant */
157	    { 0x1.fffffep6,	0x1.ffff4ep+127,	6e30,	0.0 },
158#endif
159#ifdef T_LIBM_PLUS_INF
160	    { T_LIBM_PLUS_INF,	T_LIBM_PLUS_INF,	0.0,	0.0 },
161#endif
162
163	    /* The few values from the old tests */
164	    /* Results from i386/amd64, d_eps needed on i386 */
165	    {  1.1,	0x1.125fbee250664p+1,	0x1p-52,	0x1.8p-22 },
166	    {  2.2,	0x1.2611186bae675p+2,	0x1p-51,	0x1.8p-21 },
167	    {  3.3,	0x1.3b2c47bff8328p+3,	0x1p-50,	0x1.8p-20 },
168	    {  4.4,	0x1.51cb453b9536ep+4,	0x1p-49,	0x1.8p-19 },
169	    {  5.5,	0x1.6a09e667f3bcdp+5,	0x1p-48,	0x1.8p-18 },
170	    {  6.6,	0x1.8406003b2ae5bp+6,	0x1p-47,	0x1.8p-17 },
171	    /*
172	     * These two currently fail for 'float'.
173	     * 8.8 is definitely out by more than it should be.
174	     */
175	    {  7.7,	0x1.9fdf8bcce533ep+7,	0x1p-46,	0x1.8p-16 },
176	    {  8.8,	0x1.bdb8cdadbe124p+8,	0x1p-45,	0x1.8p-15 },
177	};
178	unsigned int i;
179
180#ifdef __FreeBSD__
181	atf_tc_expect_fail("Some of the cases produce failures on FreeBSD "
182	    "due to the error epsilon being so small");
183#endif
184
185	for (i = 0; i < __arraycount(v); i++) {
186		T_LIBM_CHECK(i, exp2, v[i].x, v[i].y, v[i].d_eps);
187		if (i > 1)
188			T_LIBM_CHECK(i, exp2f, v[i].x, v[i].y, v[i].f_eps);
189	}
190}
191
192
193/*
194 * exp(3)
195 */
196ATF_TC(exp_nan);
197ATF_TC_HEAD(exp_nan, tc)
198{
199	atf_tc_set_md_var(tc, "descr", "Test exp(NaN) == NaN");
200}
201
202ATF_TC_BODY(exp_nan, tc)
203{
204	const double x = 0.0L / 0.0L;
205
206	if (isnan(exp(x)) == 0)
207		atf_tc_fail_nonfatal("exp(NaN) != NaN");
208}
209
210ATF_TC(exp_inf_neg);
211ATF_TC_HEAD(exp_inf_neg, tc)
212{
213	atf_tc_set_md_var(tc, "descr", "Test exp(-Inf) == +0.0");
214}
215
216ATF_TC_BODY(exp_inf_neg, tc)
217{
218	const double x = -1.0L / 0.0L;
219	double y = exp(x);
220
221	if (fabs(y) > 0.0 || signbit(y) != 0)
222		atf_tc_fail_nonfatal("exp(-Inf) != +0.0");
223}
224
225ATF_TC(exp_inf_pos);
226ATF_TC_HEAD(exp_inf_pos, tc)
227{
228	atf_tc_set_md_var(tc, "descr", "Test exp(+Inf) == +Inf");
229}
230
231ATF_TC_BODY(exp_inf_pos, tc)
232{
233	const double x = 1.0L / 0.0L;
234	double y = exp(x);
235
236	if (isinf(y) == 0 || signbit(y) != 0)
237		atf_tc_fail_nonfatal("exp(+Inf) != +Inf");
238}
239
240ATF_TC(exp_product);
241ATF_TC_HEAD(exp_product, tc)
242{
243	atf_tc_set_md_var(tc, "descr", "Test some selected exp(x)");
244}
245
246ATF_TC_BODY(exp_product, tc)
247{
248	double eps;
249	double x;
250	double y;
251	size_t i;
252
253	for (i = 0; i < __arraycount(exp_values); i++) {
254		x = exp_values[i].x;
255		y = exp_values[i].y;
256		eps = 1e-15 * exp_values[i].e;
257
258		if (fabs(exp(x) - y) > eps)
259			atf_tc_fail_nonfatal("exp(%0.01f) != %18.18e", x, y);
260	}
261}
262
263ATF_TC(exp_zero_neg);
264ATF_TC_HEAD(exp_zero_neg, tc)
265{
266	atf_tc_set_md_var(tc, "descr", "Test exp(-0.0) == 1.0");
267}
268
269ATF_TC_BODY(exp_zero_neg, tc)
270{
271	const double x = -0.0L;
272
273	if (fabs(exp(x) - 1.0) > 0.0)
274		atf_tc_fail_nonfatal("exp(-0.0) != 1.0");
275}
276
277ATF_TC(exp_zero_pos);
278ATF_TC_HEAD(exp_zero_pos, tc)
279{
280	atf_tc_set_md_var(tc, "descr", "Test exp(+0.0) == 1.0");
281}
282
283ATF_TC_BODY(exp_zero_pos, tc)
284{
285	const double x = 0.0L;
286
287	if (fabs(exp(x) - 1.0) > 0.0)
288		atf_tc_fail_nonfatal("exp(+0.0) != 1.0");
289}
290
291/*
292 * expf(3)
293 */
294ATF_TC(expf_nan);
295ATF_TC_HEAD(expf_nan, tc)
296{
297	atf_tc_set_md_var(tc, "descr", "Test expf(NaN) == NaN");
298}
299
300ATF_TC_BODY(expf_nan, tc)
301{
302	const float x = 0.0L / 0.0L;
303
304	if (isnan(expf(x)) == 0)
305		atf_tc_fail_nonfatal("expf(NaN) != NaN");
306}
307
308ATF_TC(expf_inf_neg);
309ATF_TC_HEAD(expf_inf_neg, tc)
310{
311	atf_tc_set_md_var(tc, "descr", "Test expf(-Inf) == +0.0");
312}
313
314ATF_TC_BODY(expf_inf_neg, tc)
315{
316	const float x = -1.0L / 0.0L;
317	float y = expf(x);
318
319	if (fabsf(y) > 0.0 || signbit(y) != 0)
320		atf_tc_fail_nonfatal("expf(-Inf) != +0.0");
321}
322
323ATF_TC(expf_inf_pos);
324ATF_TC_HEAD(expf_inf_pos, tc)
325{
326	atf_tc_set_md_var(tc, "descr", "Test expf(+Inf) == +Inf");
327}
328
329ATF_TC_BODY(expf_inf_pos, tc)
330{
331	const float x = 1.0L / 0.0L;
332	float y = expf(x);
333
334	if (isinf(y) == 0 || signbit(y) != 0)
335		atf_tc_fail_nonfatal("expf(+Inf) != +Inf");
336}
337
338ATF_TC(expf_product);
339ATF_TC_HEAD(expf_product, tc)
340{
341	atf_tc_set_md_var(tc, "descr", "Test some selected expf(x)");
342}
343
344ATF_TC_BODY(expf_product, tc)
345{
346	float eps;
347	float x;
348	float y;
349	size_t i;
350
351	for (i = 0; i < __arraycount(exp_values); i++) {
352		x = exp_values[i].x;
353		y = exp_values[i].y;
354		eps = 1e-6 * exp_values[i].e;
355
356		if (fabsf(expf(x) - y) > eps)
357			atf_tc_fail_nonfatal("expf(%0.01f) != %18.18e", x, y);
358	}
359}
360
361ATF_TC(expf_zero_neg);
362ATF_TC_HEAD(expf_zero_neg, tc)
363{
364	atf_tc_set_md_var(tc, "descr", "Test expf(-0.0) == 1.0");
365}
366
367ATF_TC_BODY(expf_zero_neg, tc)
368{
369	const float x = -0.0L;
370
371	if (fabsf(expf(x) - 1.0f) > 0.0)
372		atf_tc_fail_nonfatal("expf(-0.0) != 1.0");
373}
374
375ATF_TC(expf_zero_pos);
376ATF_TC_HEAD(expf_zero_pos, tc)
377{
378	atf_tc_set_md_var(tc, "descr", "Test expf(+0.0) == 1.0");
379}
380
381ATF_TC_BODY(expf_zero_pos, tc)
382{
383	const float x = 0.0L;
384
385	if (fabsf(expf(x) - 1.0f) > 0.0)
386		atf_tc_fail_nonfatal("expf(+0.0) != 1.0");
387}
388
389/*
390 * expm1(3)
391 */
392ATF_TC(expm1_nan);
393ATF_TC_HEAD(expm1_nan, tc)
394{
395	atf_tc_set_md_var(tc, "descr", "Test expm1(NaN) == NaN");
396}
397
398ATF_TC_BODY(expm1_nan, tc)
399{
400	const double x = 0.0L / 0.0L;
401
402	if (isnan(expm1(x)) == 0)
403		atf_tc_fail_nonfatal("expm1(NaN) != NaN");
404}
405
406ATF_TC(expm1_inf_neg);
407ATF_TC_HEAD(expm1_inf_neg, tc)
408{
409	atf_tc_set_md_var(tc, "descr", "Test expm1(-Inf) == -1");
410}
411
412ATF_TC_BODY(expm1_inf_neg, tc)
413{
414	const double x = -1.0L / 0.0L;
415
416	if (expm1(x) != -1.0)
417		atf_tc_fail_nonfatal("expm1(-Inf) != -1.0");
418}
419
420ATF_TC(expm1_inf_pos);
421ATF_TC_HEAD(expm1_inf_pos, tc)
422{
423	atf_tc_set_md_var(tc, "descr", "Test expm1(+Inf) == +Inf");
424}
425
426ATF_TC_BODY(expm1_inf_pos, tc)
427{
428	const double x = 1.0L / 0.0L;
429	double y = expm1(x);
430
431	if (isinf(y) == 0 || signbit(y) != 0)
432		atf_tc_fail_nonfatal("expm1(+Inf) != +Inf");
433}
434
435ATF_TC(expm1_zero_neg);
436ATF_TC_HEAD(expm1_zero_neg, tc)
437{
438	atf_tc_set_md_var(tc, "descr", "Test expm1(-0.0) == -0.0");
439}
440
441ATF_TC_BODY(expm1_zero_neg, tc)
442{
443	const double x = -0.0L;
444	double y = expm1(x);
445
446	if (fabs(y) > 0.0 || signbit(y) == 0)
447		atf_tc_fail_nonfatal("expm1(-0.0) != -0.0");
448}
449
450ATF_TC(expm1_zero_pos);
451ATF_TC_HEAD(expm1_zero_pos, tc)
452{
453	atf_tc_set_md_var(tc, "descr", "Test expm1(+0.0) == 1.0");
454}
455
456ATF_TC_BODY(expm1_zero_pos, tc)
457{
458	const double x = 0.0L;
459	double y = expm1(x);
460
461	if (fabs(y) > 0.0 || signbit(y) != 0)
462		atf_tc_fail_nonfatal("expm1(+0.0) != +0.0");
463}
464
465/*
466 * expm1f(3)
467 */
468ATF_TC(expm1f_nan);
469ATF_TC_HEAD(expm1f_nan, tc)
470{
471	atf_tc_set_md_var(tc, "descr", "Test expm1f(NaN) == NaN");
472}
473
474ATF_TC_BODY(expm1f_nan, tc)
475{
476	const float x = 0.0L / 0.0L;
477
478	if (isnan(expm1f(x)) == 0)
479		atf_tc_fail_nonfatal("expm1f(NaN) != NaN");
480}
481
482ATF_TC(expm1f_inf_neg);
483ATF_TC_HEAD(expm1f_inf_neg, tc)
484{
485	atf_tc_set_md_var(tc, "descr", "Test expm1f(-Inf) == -1");
486}
487
488ATF_TC_BODY(expm1f_inf_neg, tc)
489{
490	const float x = -1.0L / 0.0L;
491
492	if (expm1f(x) != -1.0)
493		atf_tc_fail_nonfatal("expm1f(-Inf) != -1.0");
494}
495
496ATF_TC(expm1f_inf_pos);
497ATF_TC_HEAD(expm1f_inf_pos, tc)
498{
499	atf_tc_set_md_var(tc, "descr", "Test expm1f(+Inf) == +Inf");
500}
501
502ATF_TC_BODY(expm1f_inf_pos, tc)
503{
504	const float x = 1.0L / 0.0L;
505	float y = expm1f(x);
506
507	if (isinf(y) == 0 || signbit(y) != 0)
508		atf_tc_fail_nonfatal("expm1f(+Inf) != +Inf");
509}
510
511ATF_TC(expm1f_zero_neg);
512ATF_TC_HEAD(expm1f_zero_neg, tc)
513{
514	atf_tc_set_md_var(tc, "descr", "Test expm1f(-0.0) == -0.0");
515}
516
517ATF_TC_BODY(expm1f_zero_neg, tc)
518{
519	const float x = -0.0L;
520	float y = expm1f(x);
521
522	if (fabsf(y) > 0.0 || signbit(y) == 0)
523		atf_tc_fail_nonfatal("expm1f(-0.0) != -0.0");
524}
525
526ATF_TC(expm1f_zero_pos);
527ATF_TC_HEAD(expm1f_zero_pos, tc)
528{
529	atf_tc_set_md_var(tc, "descr", "Test expm1f(+0.0) == 1.0");
530}
531
532ATF_TC_BODY(expm1f_zero_pos, tc)
533{
534	const float x = 0.0L;
535	float y = expm1f(x);
536
537	if (fabsf(y) > 0.0 || signbit(y) != 0)
538		atf_tc_fail_nonfatal("expm1f(+0.0) != +0.0");
539}
540
541ATF_TP_ADD_TCS(tp)
542{
543
544	ATF_TP_ADD_TC(tp, exp2_is_nan);
545	ATF_TP_ADD_TC(tp, exp2_is_plus_zero);
546	ATF_TP_ADD_TC(tp, exp2_values);
547	ATF_TP_ADD_TC(tp, exp2_powers);
548
549	ATF_TP_ADD_TC(tp, exp_nan);
550	ATF_TP_ADD_TC(tp, exp_inf_neg);
551	ATF_TP_ADD_TC(tp, exp_inf_pos);
552	ATF_TP_ADD_TC(tp, exp_product);
553	ATF_TP_ADD_TC(tp, exp_zero_neg);
554	ATF_TP_ADD_TC(tp, exp_zero_pos);
555
556	ATF_TP_ADD_TC(tp, expf_nan);
557	ATF_TP_ADD_TC(tp, expf_inf_neg);
558	ATF_TP_ADD_TC(tp, expf_inf_pos);
559	ATF_TP_ADD_TC(tp, expf_product);
560	ATF_TP_ADD_TC(tp, expf_zero_neg);
561	ATF_TP_ADD_TC(tp, expf_zero_pos);
562
563	ATF_TP_ADD_TC(tp, expm1_nan);
564	ATF_TP_ADD_TC(tp, expm1_inf_neg);
565	ATF_TP_ADD_TC(tp, expm1_inf_pos);
566	ATF_TP_ADD_TC(tp, expm1_zero_neg);
567	ATF_TP_ADD_TC(tp, expm1_zero_pos);
568
569	ATF_TP_ADD_TC(tp, expm1f_nan);
570	ATF_TP_ADD_TC(tp, expm1f_inf_neg);
571	ATF_TP_ADD_TC(tp, expm1f_inf_pos);
572	ATF_TP_ADD_TC(tp, expm1f_zero_neg);
573	ATF_TP_ADD_TC(tp, expm1f_zero_pos);
574
575	return atf_no_error();
576}
577