1 /* @(#)s_floor.c 5.1 93/09/24 */
2 /*
3  * ====================================================
4  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
5  *
6  * Developed at SunPro, a Sun Microsystems, Inc. business.
7  * Permission to use, copy, modify, and distribute this
8  * software is freely granted, provided that this notice
9  * is preserved.
10  * ====================================================
11  */
12 
13 /*
14  * floorl(x)
15  * Return x rounded toward -inf to integral value
16  * Method:
17  *	Bit twiddling.
18  * Exception:
19  *	Inexact flag raised if x not equal to floor(x).
20  */
21 
22 #include <math.h>
23 
24 #include "math_private.h"
25 
26 static const long double huge = 1.0e4930L;
27 
28 long double
29 floorl(long double x)
30 {
31 	int32_t i1,jj0;
32 	u_int32_t i,j,se,i0,sx;
33 	GET_LDOUBLE_WORDS(se,i0,i1,x);
34 	sx = (se>>15)&1;
35 	jj0 = (se&0x7fff)-0x3fff;
36 	if(jj0<31) {
37 	    if(jj0<0) {	/* raise inexact if x != 0 */
38 		if(huge+x>0.0) {
39 		    if(sx==0)
40 			return 0.0L;
41 		    else if(((se&0x7fff)|i0|i1)!=0)
42 			return -1.0L;
43 		}
44 	    } else {
45 		i = (0x7fffffff)>>jj0;
46 		if(((i0&i)|i1)==0) return x; /* x is integral */
47 		if(huge+x>0.0) {	/* raise inexact flag */
48 		    if(sx) {
49 			if (jj0>0 && (i0+(0x80000000>>jj0))>i0)
50 			  i0 += (0x80000000)>>jj0;
51 			else
52 			  {
53 			    i = 0x7fffffff;
54 			    ++se;
55 			  }
56 		    }
57 		    i0 &= (~i); i1=0;
58 		}
59 	    }
60 	} else if (jj0>62) {
61 	    if(jj0==0x4000) return x+x;	/* inf or NaN */
62 	    else return x;		/* x is integral */
63 	} else {
64 	    i = ((u_int32_t)(0xffffffff))>>(jj0-31);
65 	    if((i1&i)==0) return x;	/* x is integral */
66 	    if(huge+x>0.0) {		/* raise inexact flag */
67 		if(sx) {
68 		    if(jj0==31) i0+=1;
69 		    else {
70 			j = i1+(1<<(63-jj0));
71 			if(j<(u_int32_t)i1) i0 +=1 ;	/* got a carry */
72 			i1=j;
73 		    }
74 		}
75 		i1 &= (~i);
76 	    }
77 	}
78 	SET_LDOUBLE_WORDS(x,se,i0,i1);
79 	return x;
80 }
81