xref: /openbsd-src/regress/lib/libm/fenv/fenv.c (revision 487aca8602a8b87ea56bd544fa9e22c3e0dcf0dc)
1 /*	$OpenBSD: fenv.c,v 1.8 2021/06/17 12:55:38 kettenis Exp $	*/
2 
3 /*-
4  * Copyright (c) 2004 David Schultz <das@FreeBSD.org>
5  * All rights reserved.
6  *
7  * Redistribution and use in source and binary forms, with or without
8  * modification, are permitted provided that the following conditions
9  * are met:
10  * 1. Redistributions of source code must retain the above copyright
11  *    notice, this list of conditions and the following disclaimer.
12  * 2. Redistributions in binary form must reproduce the above copyright
13  *    notice, this list of conditions and the following disclaimer in the
14  *    documentation and/or other materials provided with the distribution.
15  *
16  * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
17  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
18  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
19  * ARE DISCLAIMED.  IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
20  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
21  * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
22  * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
23  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
24  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
25  * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
26  * SUCH DAMAGE.
27  */
28 
29 /*
30  * Test the correctness and C99-compliance of various fenv.h features.
31  */
32 
33 #include <sys/types.h>
34 #include <sys/wait.h>
35 #include <assert.h>
36 #include <err.h>
37 #include <fenv.h>
38 #include <float.h>
39 #include <math.h>
40 #include <signal.h>
41 #include <stdio.h>
42 #include <string.h>
43 #include <unistd.h>
44 
45 /*
46  * Implementations are permitted to define additional exception flags
47  * not specified in the standard, so it is not necessarily true that
48  * FE_ALL_EXCEPT == ALL_STD_EXCEPT.
49  */
50 #define	ALL_STD_EXCEPT	(FE_DIVBYZERO | FE_INEXACT | FE_INVALID | \
51 			 FE_OVERFLOW | FE_UNDERFLOW)
52 
53 #define	NEXCEPTS	(sizeof(std_excepts) / sizeof(std_excepts[0]))
54 
55 static const int std_excepts[] = {
56 	FE_INVALID,
57 	FE_DIVBYZERO,
58 	FE_OVERFLOW,
59 	FE_UNDERFLOW,
60 	FE_INEXACT,
61 };
62 
63 /* init_exceptsets() initializes this to the power set of std_excepts[] */
64 static int std_except_sets[1 << NEXCEPTS];
65 
66 static void init_exceptsets(void);
67 
68 static void test_dfl_env(void);
69 static void test_fegsetenv(void);
70 static void test_fegsetexceptflag(void);
71 static void test_masking(void);
72 static void test_fegsetround(void);
73 static void test_feholdupdate(void);
74 static void test_feraiseexcept(void);
75 static void test_fetestclearexcept(void);
76 
77 static int getround(void);
78 static void raiseexcept(int excepts);
79 static void trap_handler(int sig);
80 
81 #pragma STDC FENV_ACCESS ON
82 
83 int
main(int argc,char * argv[])84 main(int argc, char *argv[])
85 {
86 
87 	printf("1..8\n");
88 	init_exceptsets();
89 	test_dfl_env();
90 	printf("ok 1 - fenv\n");
91 	test_fetestclearexcept();
92 	printf("ok 2 - fenv\n");
93 	test_fegsetexceptflag();
94 	printf("ok 3 - fenv\n");
95 	test_feraiseexcept();
96 	printf("ok 4 - fenv\n");
97 	test_fegsetround();
98 	printf("ok 5 - fenv\n");
99 	test_fegsetenv();
100 	printf("ok 6 - fenv\n");
101 	test_masking();
102 	printf("ok 7 - fenv\n");
103 	test_feholdupdate();
104 	printf("ok 8 - fenv\n");
105 
106 	return (0);
107 }
108 
109 /*
110  * Initialize std_except_sets[] to the power set of std_excepts[]
111  */
112 void
init_exceptsets(void)113 init_exceptsets(void)
114 {
115 	int i, j, sr;
116 
117 	for (i = 0; i < 1 << NEXCEPTS; i++) {
118 		for (sr = i, j = 0; sr != 0; sr >>= 1, j++)
119 			std_except_sets[i] |= std_excepts[j] & ((~sr & 1) - 1);
120 	}
121 }
122 
123 /*
124  * This tests checks the default FP environment, so it must be first.
125  * The memcmp() test below may be too much to ask for, since there
126  * could be multiple machine-specific default environments.
127  */
128 static void
test_dfl_env(void)129 test_dfl_env(void)
130 {
131 #ifndef NO_STRICT_DFL_ENV
132 	fenv_t env;
133 
134 	fegetenv(&env);
135 #ifdef __amd64
136 	/* Some early amd64 CPUs set fip+fdp for non-x87 instructions */
137 	memset(&env.__x87.__others[0], 0, 14);
138 #endif /* __amd64 */
139 	assert(memcmp(&env, FE_DFL_ENV, sizeof(env)) == 0);
140 #endif
141 	assert(fetestexcept(FE_ALL_EXCEPT) == 0);
142 }
143 
144 /*
145  * Test fetestexcept() and feclearexcept().
146  */
147 static void
test_fetestclearexcept(void)148 test_fetestclearexcept(void)
149 {
150 	int excepts, i;
151 
152 	for (i = 0; i < 1 << NEXCEPTS; i++)
153 		assert(fetestexcept(std_except_sets[i]) == 0);
154 	for (i = 0; i < 1 << NEXCEPTS; i++) {
155 		excepts = std_except_sets[i];
156 
157 		/* FE_ALL_EXCEPT might be special-cased, as on i386. */
158 		raiseexcept(excepts);
159 		assert(fetestexcept(excepts) == excepts);
160 		assert(feclearexcept(FE_ALL_EXCEPT) == 0);
161 		assert(fetestexcept(FE_ALL_EXCEPT) == 0);
162 
163 		raiseexcept(excepts);
164 		assert(fetestexcept(excepts) == excepts);
165 		if ((excepts & (FE_UNDERFLOW | FE_OVERFLOW)) != 0) {
166 			excepts |= FE_INEXACT;
167 			assert((fetestexcept(ALL_STD_EXCEPT) | FE_INEXACT) ==
168 			    excepts);
169 		} else {
170 			assert(fetestexcept(ALL_STD_EXCEPT) == excepts);
171 		}
172 		assert(feclearexcept(excepts) == 0);
173 		assert(fetestexcept(ALL_STD_EXCEPT) == 0);
174 	}
175 }
176 
177 /*
178  * Test fegetexceptflag() and fesetexceptflag().
179  *
180  * Prerequisites: fetestexcept(), feclearexcept()
181  */
182 static void
test_fegsetexceptflag(void)183 test_fegsetexceptflag(void)
184 {
185 	fexcept_t flag;
186 	int excepts, i;
187 
188 	assert(fetestexcept(FE_ALL_EXCEPT) == 0);
189 	for (i = 0; i < 1 << NEXCEPTS; i++) {
190 		excepts = std_except_sets[i];
191 
192 		assert(fegetexceptflag(&flag, excepts) == 0);
193 		raiseexcept(ALL_STD_EXCEPT);
194 		assert(fesetexceptflag(&flag, excepts) == 0);
195 		assert(fetestexcept(ALL_STD_EXCEPT) ==
196 		    (ALL_STD_EXCEPT ^ excepts));
197 
198 		assert(fegetexceptflag(&flag, FE_ALL_EXCEPT) == 0);
199 		assert(feclearexcept(FE_ALL_EXCEPT) == 0);
200 		assert(fesetexceptflag(&flag, excepts) == 0);
201 		assert(fetestexcept(ALL_STD_EXCEPT) == 0);
202 		assert(fesetexceptflag(&flag, ALL_STD_EXCEPT ^ excepts) == 0);
203 		assert(fetestexcept(ALL_STD_EXCEPT) ==
204 		    (ALL_STD_EXCEPT ^ excepts));
205 
206 		assert(feclearexcept(FE_ALL_EXCEPT) == 0);
207 	}
208 }
209 
210 /*
211  * Test feraiseexcept().
212  *
213  * Prerequisites: fetestexcept(), feclearexcept()
214  */
215 static void
test_feraiseexcept(void)216 test_feraiseexcept(void)
217 {
218 	int excepts, i;
219 
220 	for (i = 0; i < 1 << NEXCEPTS; i++) {
221 		excepts = std_except_sets[i];
222 
223 		assert(fetestexcept(FE_ALL_EXCEPT) == 0);
224 		assert(feraiseexcept(excepts) == 0);
225 		if ((excepts & (FE_UNDERFLOW | FE_OVERFLOW)) != 0) {
226 			excepts |= FE_INEXACT;
227 			assert((fetestexcept(ALL_STD_EXCEPT) | FE_INEXACT) ==
228 			    excepts);
229 		} else {
230 			assert(fetestexcept(ALL_STD_EXCEPT) == excepts);
231 		}
232 		assert(feclearexcept(FE_ALL_EXCEPT) == 0);
233 	}
234 	assert(feraiseexcept(FE_INVALID | FE_DIVBYZERO) == 0);
235 	assert(fetestexcept(ALL_STD_EXCEPT) == (FE_INVALID | FE_DIVBYZERO));
236 	assert(feraiseexcept(FE_OVERFLOW | FE_UNDERFLOW | FE_INEXACT) == 0);
237 	assert(fetestexcept(ALL_STD_EXCEPT) == ALL_STD_EXCEPT);
238 	assert(feclearexcept(FE_ALL_EXCEPT) == 0);
239 }
240 
241 /*
242  * Test fegetround() and fesetround().
243  */
244 static void
test_fegsetround(void)245 test_fegsetround(void)
246 {
247 
248 	assert(fegetround() == FE_TONEAREST);
249 	assert(getround() == FE_TONEAREST);
250 	assert(FLT_ROUNDS == 1);
251 
252 	assert(fesetround(FE_DOWNWARD) == 0);
253 	assert(fegetround() == FE_DOWNWARD);
254 	assert(getround() == FE_DOWNWARD);
255 	assert(FLT_ROUNDS == 3);
256 
257 	assert(fesetround(FE_UPWARD) == 0);
258 	assert(getround() == FE_UPWARD);
259 	assert(fegetround() == FE_UPWARD);
260 	assert(FLT_ROUNDS == 2);
261 
262 	assert(fesetround(FE_TOWARDZERO) == 0);
263 	assert(getround() == FE_TOWARDZERO);
264 	assert(fegetround() == FE_TOWARDZERO);
265 	assert(FLT_ROUNDS == 0);
266 
267 	assert(fesetround(FE_TONEAREST) == 0);
268 	assert(getround() == FE_TONEAREST);
269 	assert(FLT_ROUNDS == 1);
270 
271 	assert(feclearexcept(FE_ALL_EXCEPT) == 0);
272 }
273 
274 /*
275  * Test fegetenv() and fesetenv().
276  *
277  * Prerequisites: fetestexcept(), feclearexcept(), fegetround(), fesetround()
278  */
279 static void
test_fegsetenv(void)280 test_fegsetenv(void)
281 {
282 	fenv_t env1, env2;
283 	int excepts, i;
284 
285 	for (i = 0; i < 1 << NEXCEPTS; i++) {
286 		excepts = std_except_sets[i];
287 
288 		assert(fetestexcept(FE_ALL_EXCEPT) == 0);
289 		assert(fegetround() == FE_TONEAREST);
290 		assert(fegetenv(&env1) == 0);
291 
292 		/*
293 		 * fe[gs]etenv() should be able to save and restore
294 		 * exception flags without the spurious inexact
295 		 * exceptions that afflict raiseexcept().
296 		 */
297 		raiseexcept(excepts);
298 		if ((excepts & (FE_UNDERFLOW | FE_OVERFLOW)) != 0 &&
299 		    (excepts & FE_INEXACT) == 0)
300 			assert(feclearexcept(FE_INEXACT) == 0);
301 
302 		fesetround(FE_DOWNWARD);
303 		assert(fegetenv(&env2) == 0);
304 		assert(fesetenv(&env1) == 0);
305 		assert(fetestexcept(FE_ALL_EXCEPT) == 0);
306 		assert(fegetround() == FE_TONEAREST);
307 
308 		assert(fesetenv(&env2) == 0);
309 		assert(fetestexcept(FE_ALL_EXCEPT) == excepts);
310 		assert(fegetround() == FE_DOWNWARD);
311 		assert(fesetenv(&env1) == 0);
312 		assert(fetestexcept(FE_ALL_EXCEPT) == 0);
313 		assert(fegetround() == FE_TONEAREST);
314 	}
315 }
316 
317 /*
318  * Test fegetexcept(), fedisableexcept(), and feenableexcept().
319  *
320  * Prerequisites: fetestexcept(), feraiseexcept()
321  */
322 static void
test_masking(void)323 test_masking(void)
324 {
325 #if !defined(__arm__) && !defined(__aarch64__) && !defined(__riscv)
326 	struct sigaction act;
327 	int except, i, pass, raise, status;
328 
329 	assert((fegetexcept() & ALL_STD_EXCEPT) == 0);
330 	assert((feenableexcept(FE_INVALID|FE_OVERFLOW) & ALL_STD_EXCEPT) == 0);
331 	assert((feenableexcept(FE_UNDERFLOW) & ALL_STD_EXCEPT) ==
332 	    (FE_INVALID | FE_OVERFLOW));
333 	assert((fedisableexcept(FE_OVERFLOW) & ALL_STD_EXCEPT) ==
334 	    (FE_INVALID | FE_OVERFLOW | FE_UNDERFLOW));
335 	assert((fegetexcept() & ALL_STD_EXCEPT) == (FE_INVALID | FE_UNDERFLOW));
336 	assert((fedisableexcept(FE_ALL_EXCEPT) & ALL_STD_EXCEPT) ==
337 	    (FE_INVALID | FE_UNDERFLOW));
338 	assert((fegetexcept() & ALL_STD_EXCEPT) == 0);
339 
340 	sigemptyset(&act.sa_mask);
341 	act.sa_flags = 0;
342 	act.sa_handler = trap_handler;
343 	for (pass = 0; pass < 2; pass++) {
344 		for (i = 0; i < NEXCEPTS; i++) {
345 			except = std_excepts[i];
346 			/* over/underflow may also raise inexact */
347 			if (except == FE_INEXACT)
348 				raise = FE_DIVBYZERO | FE_INVALID;
349 			else
350 				raise = ALL_STD_EXCEPT ^ except;
351 
352 			/*
353 			 * We need to fork a child process because
354 			 * there isn't a portable way to recover from
355 			 * a floating-point exception.
356 			 */
357 			switch(fork()) {
358 			case 0:		/* child */
359 				assert((fegetexcept() & ALL_STD_EXCEPT) == 0);
360 				assert((feenableexcept(except)
361 					   & ALL_STD_EXCEPT) == 0);
362 				assert(fegetexcept() == except);
363 				raiseexcept(raise);
364 				assert(feraiseexcept(raise) == 0);
365 				assert(fetestexcept(ALL_STD_EXCEPT) == raise);
366 
367 				assert(sigaction(SIGFPE, &act, NULL) == 0);
368 				switch (pass) {
369 				case 0:
370 					raiseexcept(except);
371 				case 1:
372 					feraiseexcept(except);
373 				default:
374 					assert(0);
375 				}
376 				assert(0);
377 			default:	/* parent */
378 				assert(wait(&status) > 0);
379 				/*
380 				 * Avoid assert() here so that it's possible
381 				 * to examine a failed child's core dump.
382 				 */
383 				if (!WIFEXITED(status))
384 					errx(1, "child aborted\n");
385 				assert(WEXITSTATUS(status) == 0);
386 				break;
387 			case -1:	/* error */
388 				assert(0);
389 			}
390 		}
391 	}
392 	assert(fetestexcept(FE_ALL_EXCEPT) == 0);
393 #endif
394 }
395 
396 /*
397  * Test feholdexcept() and feupdateenv().
398  *
399  * Prerequisites: fetestexcept(), fegetround(), fesetround(),
400  *	fedisableexcept(), feenableexcept()
401  */
402 static void
test_feholdupdate(void)403 test_feholdupdate(void)
404 {
405 	fenv_t env;
406 
407 	struct sigaction act;
408 	int except, i, pass, status, raise;
409 
410 	sigemptyset(&act.sa_mask);
411 	act.sa_flags = 0;
412 	act.sa_handler = trap_handler;
413 	for (pass = 0; pass < 2; pass++) {
414 		for (i = 0; i < NEXCEPTS; i++) {
415 			except = std_excepts[i];
416 			/* over/underflow may also raise inexact */
417 			if (except == FE_INEXACT)
418 				raise = FE_DIVBYZERO | FE_INVALID;
419 			else
420 				raise = ALL_STD_EXCEPT ^ except;
421 
422 			/*
423 			 * We need to fork a child process because
424 			 * there isn't a portable way to recover from
425 			 * a floating-point exception.
426 			 */
427 			switch(fork()) {
428 			case 0:		/* child */
429 				/*
430 				 * We don't want to cause a fatal exception in
431 				 * the child until the second pass, so we can
432 				 * check other properties of feupdateenv().
433 				 */
434 				if (pass == 1)
435 					assert((feenableexcept(except) &
436 						   ALL_STD_EXCEPT) == 0);
437 				raiseexcept(raise);
438 				assert(fesetround(FE_DOWNWARD) == 0);
439 				assert(feholdexcept(&env) == 0);
440 				assert(fetestexcept(FE_ALL_EXCEPT) == 0);
441 				raiseexcept(except);
442 				assert(fesetround(FE_UPWARD) == 0);
443 
444 				if (pass == 1)
445 					assert(sigaction(SIGFPE, &act, NULL) ==
446 					    0);
447 				assert(feupdateenv(&env) == 0);
448 				assert(fegetround() == FE_DOWNWARD);
449 				assert(fetestexcept(ALL_STD_EXCEPT) ==
450 				    (except | raise));
451 
452 				assert(pass == 0);
453 				_exit(0);
454 			default:	/* parent */
455 				assert(wait(&status) > 0);
456 				/*
457 				 * Avoid assert() here so that it's possible
458 				 * to examine a failed child's core dump.
459 				 */
460 				if (!WIFEXITED(status))
461 					errx(1, "child aborted\n");
462 				assert(WEXITSTATUS(status) == 0);
463 				break;
464 			case -1:	/* error */
465 				assert(0);
466 			}
467 		}
468 #if defined(__arm__) || defined(__aarch64__) || defined(__riscv)
469 		break;
470 #endif
471 	}
472 	assert(fetestexcept(FE_ALL_EXCEPT) == 0);
473 }
474 
475 /*
476  * Raise a floating-point exception without relying on the standard
477  * library routines, which we are trying to test.
478  *
479  * XXX We can't raise an {over,under}flow without also raising an
480  * inexact exception.
481  */
482 static void
raiseexcept(int excepts)483 raiseexcept(int excepts)
484 {
485 	volatile double d;
486 
487 	/*
488 	 * With a compiler that supports the FENV_ACCESS pragma
489 	 * properly, simple expressions like '0.0 / 0.0' should
490 	 * be sufficient to generate traps.  Unfortunately, we
491 	 * need to bring a volatile variable into the equation
492 	 * to prevent incorrect optimizations.
493 	 */
494 	if (excepts & FE_INVALID) {
495 		d = 0.0;
496 		d = 0.0 / d;
497 	}
498 	if (excepts & FE_DIVBYZERO) {
499 		d = 0.0;
500 		d = 1.0 / d;
501 	}
502 	if (excepts & FE_OVERFLOW) {
503 		d = DBL_MAX;
504 		d *= 2.0;
505 	}
506 	if (excepts & FE_UNDERFLOW) {
507 		d = DBL_MIN;
508 		d /= DBL_MAX;
509 	}
510 	if (excepts & FE_INEXACT) {
511 		d = DBL_MIN;
512 		d += 1.0;
513 	}
514 
515 	/*
516 	 * On the x86 (and some other architectures?) the FPU and
517 	 * integer units are decoupled.  We need to execute an FWAIT
518 	 * or a floating-point instruction to get synchronous exceptions.
519 	 */
520 	d = 1.0;
521 	d += 1.0;
522 }
523 
524 /*
525  * Determine the current rounding mode without relying on the fenv
526  * routines.  This function may raise an inexact exception.
527  */
528 static int
getround(void)529 getround(void)
530 {
531 	volatile double d, e;
532 
533 	/*
534 	 * This test works just as well with 0.0 - 0.0, except on ia64
535 	 * where 0.0 - 0.0 gives the wrong sign when rounding downwards.
536 	 * On i386 use two volatile variables d and e to retrieve the
537 	 * value out of the x87 and force rounding.
538 	 * https://bugs.freebsd.org/bugzilla/show_bug.cgi?id=204671
539 	 */
540 	d = 1.0;
541 	d -= 1.0;
542 	e = copysign(1.0, d);
543 	if (e < 0.0)
544 		return (FE_DOWNWARD);
545 
546 	d = 1.0;
547 	e = d + (DBL_EPSILON * 3.0 / 4.0);
548 	if (e == 1.0)
549 		return (FE_TOWARDZERO);
550 	e = d + (DBL_EPSILON * 1.0 / 4.0);
551 	if (e > 1.0)
552 		return (FE_UPWARD);
553 
554 	return (FE_TONEAREST);
555 }
556 
557 static void
trap_handler(int sig)558 trap_handler(int sig)
559 {
560 
561 	assert(sig == SIGFPE);
562 	_exit(0);
563 }
564