Home | History | Annotate | Line # | Download | only in libm 
      1 /* $NetBSD: t_round.c,v 1.10 2024/02/02 22:39:09 andvar Exp $ */
      2 
      3 /*-
      4  * Copyright (c) 2011 The NetBSD Foundation, Inc.
      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 NETBSD FOUNDATION, INC. AND CONTRIBUTORS
     17  * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
     18  * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
     19  * PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR CONTRIBUTORS
     20  * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
     21  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
     22  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
     23  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
     24  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
     25  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
     26  * POSSIBILITY OF SUCH DAMAGE.
     27  */
     28 
     29 #include <sys/param.h>
     30 
     31 #include <atf-c.h>
     32 #include <float.h>
     33 #include <math.h>
     34 #include <stdio.h>
     35 #include <stdint.h>
     36 
     37 /*
     38  * This tests for a bug in the initial implementation where
     39  * precision was lost in an internal subtraction, leading to
     40  * rounding into the wrong direction.
     41  */
     42 
     43 /* 0.5 - EPSILON */
     44 #define VAL	0x0.7ffffffffffffcp0
     45 #define VALF	0x0.7fffff8p0
     46 #define VALL	(0.5 - LDBL_EPSILON)
     47 
     48 #ifdef __vax__
     49 #define SMALL_NUM	1.0e-38
     50 #else
     51 #define SMALL_NUM	1.0e-40
     52 #endif
     53 
     54 ATF_TC(round_dir);
     55 ATF_TC_HEAD(round_dir, tc)
     56 {
     57 	atf_tc_set_md_var(tc, "descr","Check for rounding in wrong direction");
     58 }
     59 
     60 ATF_TC_BODY(round_dir, tc)
     61 {
     62 	double a = VAL, b, c;
     63 	float af = VALF, bf, cf;
     64 	long double al = VALL, bl, cl;
     65 
     66 	b = round(a);
     67 	bf = roundf(af);
     68 	bl = roundl(al);
     69 
     70 	ATF_CHECK(fabs(b) < SMALL_NUM);
     71 	ATF_CHECK(fabsf(bf) < SMALL_NUM);
     72 	ATF_CHECK(fabsl(bl) < SMALL_NUM);
     73 
     74 	c = round(-a);
     75 	cf = roundf(-af);
     76 	cl = roundl(-al);
     77 
     78 	ATF_CHECK(fabs(c) < SMALL_NUM);
     79 	ATF_CHECK(fabsf(cf) < SMALL_NUM);
     80 	ATF_CHECK(fabsl(cl) < SMALL_NUM);
     81 }
     82 
     83 ATF_TC(rounding_alpha);
     84 ATF_TC_HEAD(rounding_alpha, tc)
     85 {
     86 	atf_tc_set_md_var(tc, "descr","Checking MPFR's config failure with -mieee on Alpha");
     87 }
     88 
     89 typedef uint64_t gimpy_limb_t;
     90 #define GIMPY_NUMB_BITS 64
     91 
     92 ATF_TC_BODY(rounding_alpha, tc)
     93 {
     94         double d;
     95         gimpy_limb_t u;
     96         int i;
     97 
     98         d = 1.0;
     99         for (i = 0; i < GIMPY_NUMB_BITS - 1; i++)
    100                 d = d + d;
    101 
    102         printf("d = %g\n", d);
    103         u = (gimpy_limb_t) d;
    104 
    105         for (; i > 0; i--) {
    106                 ATF_CHECK_MSG((u % 2 == 0),
    107 		    "%"PRIu64" is not an even number! (iteration %d)", u , i);
    108                 u = u >> 1;
    109         }
    110 }
    111 
    112 ATF_TC(rounding_alpha_simple);
    113 ATF_TC_HEAD(rounding_alpha_simple, tc)
    114 {
    115 	atf_tc_set_md_var(tc, "descr","Checking double to uint64_t edge case");
    116 }
    117 
    118 
    119 double rounding_alpha_simple_even = 9223372036854775808.000000; /* 2^63 */
    120 
    121 ATF_TC_BODY(rounding_alpha_simple, tc)
    122 {
    123 	uint64_t unsigned_even = rounding_alpha_simple_even;
    124 
    125 	ATF_CHECK_MSG(unsigned_even % 2 == 0,
    126 	    "2^63 cast to uint64_t is odd (got %"PRIu64")", unsigned_even);
    127 
    128 }
    129 ATF_TP_ADD_TCS(tp)
    130 {
    131 
    132 	ATF_TP_ADD_TC(tp, round_dir);
    133 	ATF_TP_ADD_TC(tp, rounding_alpha);
    134 	ATF_TP_ADD_TC(tp, rounding_alpha_simple);
    135 
    136 	return atf_no_error();
    137 }
    138