xref: /freebsd/lib/msun/src/s_rint.c (revision 3a8617a83f16ffc9db4f96e1f0f21af94078e6b1)
13a8617a8SJordan K. Hubbard /* @(#)s_rint.c 5.1 93/09/24 */
23a8617a8SJordan K. Hubbard /*
33a8617a8SJordan K. Hubbard  * ====================================================
43a8617a8SJordan K. Hubbard  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
53a8617a8SJordan K. Hubbard  *
63a8617a8SJordan K. Hubbard  * Developed at SunPro, a Sun Microsystems, Inc. business.
73a8617a8SJordan K. Hubbard  * Permission to use, copy, modify, and distribute this
83a8617a8SJordan K. Hubbard  * software is freely granted, provided that this notice
93a8617a8SJordan K. Hubbard  * is preserved.
103a8617a8SJordan K. Hubbard  * ====================================================
113a8617a8SJordan K. Hubbard  */
123a8617a8SJordan K. Hubbard 
133a8617a8SJordan K. Hubbard #ifndef lint
143a8617a8SJordan K. Hubbard static char rcsid[] = "$Id: s_rint.c,v 1.6 1994/08/18 23:07:16 jtc Exp $";
153a8617a8SJordan K. Hubbard #endif
163a8617a8SJordan K. Hubbard 
173a8617a8SJordan K. Hubbard /*
183a8617a8SJordan K. Hubbard  * rint(x)
193a8617a8SJordan K. Hubbard  * Return x rounded to integral value according to the prevailing
203a8617a8SJordan K. Hubbard  * rounding mode.
213a8617a8SJordan K. Hubbard  * Method:
223a8617a8SJordan K. Hubbard  *	Using floating addition.
233a8617a8SJordan K. Hubbard  * Exception:
243a8617a8SJordan K. Hubbard  *	Inexact flag raised if x not equal to rint(x).
253a8617a8SJordan K. Hubbard  */
263a8617a8SJordan K. Hubbard 
273a8617a8SJordan K. Hubbard #include "math.h"
283a8617a8SJordan K. Hubbard #include "math_private.h"
293a8617a8SJordan K. Hubbard 
303a8617a8SJordan K. Hubbard #ifdef __STDC__
313a8617a8SJordan K. Hubbard static const double
323a8617a8SJordan K. Hubbard #else
333a8617a8SJordan K. Hubbard static double
343a8617a8SJordan K. Hubbard #endif
353a8617a8SJordan K. Hubbard TWO52[2]={
363a8617a8SJordan K. Hubbard   4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
373a8617a8SJordan K. Hubbard  -4.50359962737049600000e+15, /* 0xC3300000, 0x00000000 */
383a8617a8SJordan K. Hubbard };
393a8617a8SJordan K. Hubbard 
403a8617a8SJordan K. Hubbard #ifdef __STDC__
413a8617a8SJordan K. Hubbard 	double rint(double x)
423a8617a8SJordan K. Hubbard #else
433a8617a8SJordan K. Hubbard 	double rint(x)
443a8617a8SJordan K. Hubbard 	double x;
453a8617a8SJordan K. Hubbard #endif
463a8617a8SJordan K. Hubbard {
473a8617a8SJordan K. Hubbard 	int32_t i0,j0,sx;
483a8617a8SJordan K. Hubbard 	u_int32_t i,i1;
493a8617a8SJordan K. Hubbard 	double w,t;
503a8617a8SJordan K. Hubbard 	EXTRACT_WORDS(i0,i1,x);
513a8617a8SJordan K. Hubbard 	sx = (i0>>31)&1;
523a8617a8SJordan K. Hubbard 	j0 = ((i0>>20)&0x7ff)-0x3ff;
533a8617a8SJordan K. Hubbard 	if(j0<20) {
543a8617a8SJordan K. Hubbard 	    if(j0<0) {
553a8617a8SJordan K. Hubbard 		if(((i0&0x7fffffff)|i1)==0) return x;
563a8617a8SJordan K. Hubbard 		i1 |= (i0&0x0fffff);
573a8617a8SJordan K. Hubbard 		i0 &= 0xfffe0000;
583a8617a8SJordan K. Hubbard 		i0 |= ((i1|-i1)>>12)&0x80000;
593a8617a8SJordan K. Hubbard 		SET_HIGH_WORD(x,i0);
603a8617a8SJordan K. Hubbard 	        w = TWO52[sx]+x;
613a8617a8SJordan K. Hubbard 	        t =  w-TWO52[sx];
623a8617a8SJordan K. Hubbard 		GET_HIGH_WORD(i0,t);
633a8617a8SJordan K. Hubbard 		SET_HIGH_WORD(t,(i0&0x7fffffff)|(sx<<31));
643a8617a8SJordan K. Hubbard 	        return t;
653a8617a8SJordan K. Hubbard 	    } else {
663a8617a8SJordan K. Hubbard 		i = (0x000fffff)>>j0;
673a8617a8SJordan K. Hubbard 		if(((i0&i)|i1)==0) return x; /* x is integral */
683a8617a8SJordan K. Hubbard 		i>>=1;
693a8617a8SJordan K. Hubbard 		if(((i0&i)|i1)!=0) {
703a8617a8SJordan K. Hubbard 		    if(j0==19) i1 = 0x40000000; else
713a8617a8SJordan K. Hubbard 		    i0 = (i0&(~i))|((0x20000)>>j0);
723a8617a8SJordan K. Hubbard 		}
733a8617a8SJordan K. Hubbard 	    }
743a8617a8SJordan K. Hubbard 	} else if (j0>51) {
753a8617a8SJordan K. Hubbard 	    if(j0==0x400) return x+x;	/* inf or NaN */
763a8617a8SJordan K. Hubbard 	    else return x;		/* x is integral */
773a8617a8SJordan K. Hubbard 	} else {
783a8617a8SJordan K. Hubbard 	    i = ((u_int32_t)(0xffffffff))>>(j0-20);
793a8617a8SJordan K. Hubbard 	    if((i1&i)==0) return x;	/* x is integral */
803a8617a8SJordan K. Hubbard 	    i>>=1;
813a8617a8SJordan K. Hubbard 	    if((i1&i)!=0) i1 = (i1&(~i))|((0x40000000)>>(j0-20));
823a8617a8SJordan K. Hubbard 	}
833a8617a8SJordan K. Hubbard 	INSERT_WORDS(x,i0,i1);
843a8617a8SJordan K. Hubbard 	w = TWO52[sx]+x;
853a8617a8SJordan K. Hubbard 	return w-TWO52[sx];
863a8617a8SJordan K. Hubbard }
87