Home | History | Annotate | Line # | Download | only in spmath
      1  1.5     skrll /*	$NetBSD: frnd.c,v 1.5 2012/02/04 17:03:10 skrll Exp $	*/
      2  1.1  fredette 
      3  1.1  fredette /*	$OpenBSD: frnd.c,v 1.5 2001/03/29 03:58:18 mickey Exp $	*/
      4  1.1  fredette 
      5  1.1  fredette /*
      6  1.1  fredette  * Copyright 1996 1995 by Open Software Foundation, Inc.
      7  1.1  fredette  *              All Rights Reserved
      8  1.1  fredette  *
      9  1.1  fredette  * Permission to use, copy, modify, and distribute this software and
     10  1.1  fredette  * its documentation for any purpose and without fee is hereby granted,
     11  1.1  fredette  * provided that the above copyright notice appears in all copies and
     12  1.1  fredette  * that both the copyright notice and this permission notice appear in
     13  1.1  fredette  * supporting documentation.
     14  1.1  fredette  *
     15  1.1  fredette  * OSF DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE
     16  1.1  fredette  * INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
     17  1.1  fredette  * FOR A PARTICULAR PURPOSE.
     18  1.1  fredette  *
     19  1.1  fredette  * IN NO EVENT SHALL OSF BE LIABLE FOR ANY SPECIAL, INDIRECT, OR
     20  1.1  fredette  * CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM
     21  1.1  fredette  * LOSS OF USE, DATA OR PROFITS, WHETHER IN ACTION OF CONTRACT,
     22  1.1  fredette  * NEGLIGENCE, OR OTHER TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION
     23  1.1  fredette  * WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
     24  1.1  fredette  *
     25  1.1  fredette  */
     26  1.1  fredette /*
     27  1.1  fredette  * pmk1.1
     28  1.1  fredette  */
     29  1.1  fredette /*
     30  1.1  fredette  * (c) Copyright 1986 HEWLETT-PACKARD COMPANY
     31  1.1  fredette  *
     32  1.1  fredette  * To anyone who acknowledges that this file is provided "AS IS"
     33  1.1  fredette  * without any express or implied warranty:
     34  1.1  fredette  *     permission to use, copy, modify, and distribute this file
     35  1.1  fredette  * for any purpose is hereby granted without fee, provided that
     36  1.1  fredette  * the above copyright notice and this notice appears in all
     37  1.1  fredette  * copies, and that the name of Hewlett-Packard Company not be
     38  1.1  fredette  * used in advertising or publicity pertaining to distribution
     39  1.1  fredette  * of the software without specific, written prior permission.
     40  1.1  fredette  * Hewlett-Packard Company makes no representations about the
     41  1.1  fredette  * suitability of this software for any purpose.
     42  1.1  fredette  */
     43  1.2     lukem 
     44  1.2     lukem #include <sys/cdefs.h>
     45  1.5     skrll __KERNEL_RCSID(0, "$NetBSD: frnd.c,v 1.5 2012/02/04 17:03:10 skrll Exp $");
     46  1.1  fredette 
     47  1.1  fredette #include "../spmath/float.h"
     48  1.1  fredette #include "../spmath/sgl_float.h"
     49  1.1  fredette #include "../spmath/dbl_float.h"
     50  1.1  fredette #include "../spmath/quad_float.h"
     51  1.1  fredette #include "../spmath/cnv_float.h"
     52  1.1  fredette 
     53  1.1  fredette /*
     54  1.1  fredette  *  Single Floating-point Round to Integer
     55  1.1  fredette  */
     56  1.1  fredette 
     57  1.1  fredette /*ARGSUSED*/
     58  1.1  fredette int
     59  1.5     skrll sgl_frnd(sgl_floating_point *srcptr, sgl_floating_point *dstptr,
     60  1.5     skrll     unsigned int *status)
     61  1.1  fredette {
     62  1.1  fredette 	register unsigned int src, result;
     63  1.1  fredette 	register int src_exponent;
     64  1.4   thorpej 	register int inexact = false;
     65  1.1  fredette 
     66  1.1  fredette 	src = *srcptr;
     67  1.1  fredette 	/*
     68  1.1  fredette 	 * check source operand for NaN or infinity
     69  1.1  fredette 	 */
     70  1.1  fredette 	if ((src_exponent = Sgl_exponent(src)) == SGL_INFINITY_EXPONENT) {
     71  1.1  fredette 		/*
     72  1.1  fredette 		 * is signaling NaN?
     73  1.1  fredette 		 */
     74  1.1  fredette 		if (Sgl_isone_signaling(src)) {
     75  1.1  fredette 			/* trap if INVALIDTRAP enabled */
     76  1.1  fredette 			if (Is_invalidtrap_enabled()) return(INVALIDEXCEPTION);
     77  1.1  fredette 			/* make NaN quiet */
     78  1.1  fredette 			Set_invalidflag();
     79  1.1  fredette 			Sgl_set_quiet(src);
     80  1.1  fredette 		}
     81  1.1  fredette 		/*
     82  1.1  fredette 		 * return quiet NaN or infinity
     83  1.1  fredette 		 */
     84  1.1  fredette 		*dstptr = src;
     85  1.1  fredette 		return(NOEXCEPTION);
     86  1.1  fredette 	}
     87  1.1  fredette 	/*
     88  1.1  fredette 	 * Need to round?
     89  1.1  fredette 	 */
     90  1.1  fredette 	if ((src_exponent -= SGL_BIAS) >= SGL_P - 1) {
     91  1.1  fredette 		*dstptr = src;
     92  1.1  fredette 		return(NOEXCEPTION);
     93  1.1  fredette 	}
     94  1.1  fredette 	/*
     95  1.1  fredette 	 * Generate result
     96  1.1  fredette 	 */
     97  1.1  fredette 	if (src_exponent >= 0) {
     98  1.1  fredette 		Sgl_clear_exponent_set_hidden(src);
     99  1.1  fredette 		result = src;
    100  1.1  fredette 		Sgl_rightshift(result,(SGL_P-1) - (src_exponent));
    101  1.1  fredette 		/* check for inexact */
    102  1.1  fredette 		if (Sgl_isinexact_to_fix(src,src_exponent)) {
    103  1.4   thorpej 			inexact = true;
    104  1.1  fredette 			/*  round result  */
    105  1.1  fredette 			switch (Rounding_mode()) {
    106  1.1  fredette 			case ROUNDPLUS:
    107  1.1  fredette 			     if (Sgl_iszero_sign(src)) Sgl_increment(result);
    108  1.1  fredette 			     break;
    109  1.1  fredette 			case ROUNDMINUS:
    110  1.1  fredette 			     if (Sgl_isone_sign(src)) Sgl_increment(result);
    111  1.1  fredette 			     break;
    112  1.1  fredette 			case ROUNDNEAREST:
    113  1.1  fredette 			     if (Sgl_isone_roundbit(src,src_exponent))
    114  1.1  fredette 				if (Sgl_isone_stickybit(src,src_exponent)
    115  1.1  fredette 				|| (Sgl_isone_lowmantissa(result)))
    116  1.1  fredette 					Sgl_increment(result);
    117  1.1  fredette 			}
    118  1.1  fredette 		}
    119  1.1  fredette 		Sgl_leftshift(result,(SGL_P-1) - (src_exponent));
    120  1.1  fredette 		if (Sgl_isone_hiddenoverflow(result))
    121  1.1  fredette 			Sgl_set_exponent(result,src_exponent + (SGL_BIAS+1));
    122  1.1  fredette 		else Sgl_set_exponent(result,src_exponent + SGL_BIAS);
    123  1.1  fredette 	}
    124  1.1  fredette 	else {
    125  1.1  fredette 		result = src;		/* set sign */
    126  1.1  fredette 		Sgl_setzero_exponentmantissa(result);
    127  1.1  fredette 		/* check for inexact */
    128  1.1  fredette 		if (Sgl_isnotzero_exponentmantissa(src)) {
    129  1.4   thorpej 			inexact = true;
    130  1.1  fredette 			/*  round result  */
    131  1.1  fredette 			switch (Rounding_mode()) {
    132  1.1  fredette 			case ROUNDPLUS:
    133  1.1  fredette 			     if (Sgl_iszero_sign(src))
    134  1.1  fredette 				Sgl_set_exponent(result,SGL_BIAS);
    135  1.1  fredette 			     break;
    136  1.1  fredette 			case ROUNDMINUS:
    137  1.1  fredette 			     if (Sgl_isone_sign(src))
    138  1.1  fredette 				Sgl_set_exponent(result,SGL_BIAS);
    139  1.1  fredette 			     break;
    140  1.1  fredette 			case ROUNDNEAREST:
    141  1.1  fredette 			     if (src_exponent == -1)
    142  1.1  fredette 				if (Sgl_isnotzero_mantissa(src))
    143  1.1  fredette 				   Sgl_set_exponent(result,SGL_BIAS);
    144  1.1  fredette 			}
    145  1.1  fredette 		}
    146  1.1  fredette 	}
    147  1.1  fredette 	*dstptr = result;
    148  1.1  fredette 	if (inexact) {
    149  1.1  fredette 		if (Is_inexacttrap_enabled()) return(INEXACTEXCEPTION);
    150  1.1  fredette 		else Set_inexactflag();
    151  1.1  fredette 	}
    152  1.1  fredette 	return(NOEXCEPTION);
    153  1.1  fredette }
    154  1.1  fredette 
    155  1.1  fredette /*
    156  1.1  fredette  *  Double Floating-point Round to Integer
    157  1.1  fredette  */
    158  1.1  fredette 
    159  1.1  fredette /*ARGSUSED*/
    160  1.1  fredette int
    161  1.5     skrll dbl_frnd(dbl_floating_point *srcptr, dbl_floating_point *dstptr,
    162  1.5     skrll     unsigned int *status)
    163  1.1  fredette {
    164  1.1  fredette 	register unsigned int srcp1, srcp2, resultp1, resultp2;
    165  1.1  fredette 	register int src_exponent;
    166  1.4   thorpej 	register int inexact = false;
    167  1.1  fredette 
    168  1.1  fredette 	Dbl_copyfromptr(srcptr,srcp1,srcp2);
    169  1.1  fredette 	/*
    170  1.1  fredette 	 * check source operand for NaN or infinity
    171  1.1  fredette 	 */
    172  1.1  fredette 	if ((src_exponent = Dbl_exponent(srcp1)) == DBL_INFINITY_EXPONENT) {
    173  1.1  fredette 		/*
    174  1.1  fredette 		 * is signaling NaN?
    175  1.1  fredette 		 */
    176  1.1  fredette 		if (Dbl_isone_signaling(srcp1)) {
    177  1.1  fredette 			/* trap if INVALIDTRAP enabled */
    178  1.1  fredette 			if (Is_invalidtrap_enabled()) return(INVALIDEXCEPTION);
    179  1.1  fredette 			/* make NaN quiet */
    180  1.1  fredette 			Set_invalidflag();
    181  1.1  fredette 			Dbl_set_quiet(srcp1);
    182  1.1  fredette 		}
    183  1.1  fredette 		/*
    184  1.1  fredette 		 * return quiet NaN or infinity
    185  1.1  fredette 		 */
    186  1.1  fredette 		Dbl_copytoptr(srcp1,srcp2,dstptr);
    187  1.1  fredette 		return(NOEXCEPTION);
    188  1.1  fredette 	}
    189  1.1  fredette 	/*
    190  1.1  fredette 	 * Need to round?
    191  1.1  fredette 	 */
    192  1.1  fredette 	if ((src_exponent -= DBL_BIAS) >= DBL_P - 1) {
    193  1.1  fredette 		Dbl_copytoptr(srcp1,srcp2,dstptr);
    194  1.1  fredette 		return(NOEXCEPTION);
    195  1.1  fredette 	}
    196  1.1  fredette 	/*
    197  1.1  fredette 	 * Generate result
    198  1.1  fredette 	 */
    199  1.1  fredette 	if (src_exponent >= 0) {
    200  1.1  fredette 		Dbl_clear_exponent_set_hidden(srcp1);
    201  1.1  fredette 		resultp1 = srcp1;
    202  1.1  fredette 		resultp2 = srcp2;
    203  1.1  fredette 		Dbl_rightshift(resultp1,resultp2,(DBL_P-1) - (src_exponent));
    204  1.1  fredette 		/* check for inexact */
    205  1.1  fredette 		if (Dbl_isinexact_to_fix(srcp1,srcp2,src_exponent)) {
    206  1.4   thorpej 			inexact = true;
    207  1.1  fredette 			/*  round result  */
    208  1.1  fredette 			switch (Rounding_mode()) {
    209  1.1  fredette 			case ROUNDPLUS:
    210  1.1  fredette 			     if (Dbl_iszero_sign(srcp1))
    211  1.1  fredette 				Dbl_increment(resultp1,resultp2);
    212  1.1  fredette 			     break;
    213  1.1  fredette 			case ROUNDMINUS:
    214  1.1  fredette 			     if (Dbl_isone_sign(srcp1))
    215  1.1  fredette 				Dbl_increment(resultp1,resultp2);
    216  1.1  fredette 			     break;
    217  1.1  fredette 			case ROUNDNEAREST:
    218  1.1  fredette 			     if (Dbl_isone_roundbit(srcp1,srcp2,src_exponent))
    219  1.1  fredette 			      if (Dbl_isone_stickybit(srcp1,srcp2,src_exponent)
    220  1.1  fredette 				  || (Dbl_isone_lowmantissap2(resultp2)))
    221  1.1  fredette 					Dbl_increment(resultp1,resultp2);
    222  1.1  fredette 			}
    223  1.1  fredette 		}
    224  1.1  fredette 		Dbl_leftshift(resultp1,resultp2,(DBL_P-1) - (src_exponent));
    225  1.1  fredette 		if (Dbl_isone_hiddenoverflow(resultp1))
    226  1.1  fredette 			Dbl_set_exponent(resultp1,src_exponent + (DBL_BIAS+1));
    227  1.1  fredette 		else Dbl_set_exponent(resultp1,src_exponent + DBL_BIAS);
    228  1.1  fredette 	}
    229  1.1  fredette 	else {
    230  1.1  fredette 		resultp1 = srcp1;  /* set sign */
    231  1.1  fredette 		Dbl_setzero_exponentmantissa(resultp1,resultp2);
    232  1.1  fredette 		/* check for inexact */
    233  1.1  fredette 		if (Dbl_isnotzero_exponentmantissa(srcp1,srcp2)) {
    234  1.4   thorpej 			inexact = true;
    235  1.1  fredette 			/*  round result  */
    236  1.1  fredette 			switch (Rounding_mode()) {
    237  1.1  fredette 			case ROUNDPLUS:
    238  1.1  fredette 			     if (Dbl_iszero_sign(srcp1))
    239  1.1  fredette 				Dbl_set_exponent(resultp1,DBL_BIAS);
    240  1.1  fredette 			     break;
    241  1.1  fredette 			case ROUNDMINUS:
    242  1.1  fredette 			     if (Dbl_isone_sign(srcp1))
    243  1.1  fredette 				Dbl_set_exponent(resultp1,DBL_BIAS);
    244  1.1  fredette 			     break;
    245  1.1  fredette 			case ROUNDNEAREST:
    246  1.1  fredette 			     if (src_exponent == -1)
    247  1.1  fredette 				if (Dbl_isnotzero_mantissa(srcp1,srcp2))
    248  1.1  fredette 				   Dbl_set_exponent(resultp1,DBL_BIAS);
    249  1.1  fredette 			}
    250  1.1  fredette 		}
    251  1.1  fredette 	}
    252  1.1  fredette 	Dbl_copytoptr(resultp1,resultp2,dstptr);
    253  1.1  fredette 	if (inexact) {
    254  1.1  fredette 		if (Is_inexacttrap_enabled()) return(INEXACTEXCEPTION);
    255  1.1  fredette 		else Set_inexactflag();
    256  1.1  fredette 	}
    257  1.1  fredette 	return(NOEXCEPTION);
    258  1.1  fredette }
    259  1.1  fredette 
    260  1.1  fredette /*ARGSUSED*/
    261  1.1  fredette int
    262  1.5     skrll quad_frnd(quad_floating_point *srcptr, quad_floating_point *dstptr,
    263  1.5     skrll     unsigned int *status)
    264  1.1  fredette {
    265  1.1  fredette 	return(UNIMPLEMENTEDEXCEPTION);
    266  1.1  fredette }
    267  1.1  fredette 
    268