xref: /openbsd/lib/libm/src/s_lrint.c (revision cca36db2)
1 /*	$OpenBSD: s_lrint.c,v 1.7 2011/07/06 00:02:42 martynas Exp $	*/
2 /* $NetBSD: lrint.c,v 1.3 2004/10/13 15:18:32 drochner Exp $ */
3 
4 /*-
5  * Copyright (c) 2004
6  *	Matthias Drochner. All rights reserved.
7  *
8  * Redistribution and use in source and binary forms, with or without
9  * modification, are permitted provided that the following conditions
10  * are met:
11  * 1. Redistributions of source code must retain the above copyright
12  *    notice, this list of conditions and the following disclaimer.
13  * 2. Redistributions in binary form must reproduce the above copyright
14  *    notice, this list of conditions and the following disclaimer in the
15  *    documentation and/or other materials provided with the distribution.
16  *
17  * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
18  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
19  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
20  * ARE DISCLAIMED.  IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
21  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
22  * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
23  * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
24  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
25  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
26  * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
27  * SUCH DAMAGE.
28  */
29 
30 /* LINTLIBRARY */
31 
32 #include <sys/types.h>
33 #include <sys/limits.h>
34 #include <float.h>
35 #include <math.h>
36 #include <ieeefp.h>
37 #include <machine/ieee.h>
38 
39 #include "math_private.h"
40 
41 #ifndef LRINTNAME
42 #define LRINTNAME lrint
43 #define RESTYPE long int
44 #define RESTYPE_MIN LONG_MIN
45 #define RESTYPE_MAX LONG_MAX
46 #endif
47 
48 #define RESTYPE_BITS (sizeof(RESTYPE) * 8)
49 
50 static const double
51 TWO52[2]={
52   4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
53  -4.50359962737049600000e+15, /* 0xC3300000, 0x00000000 */
54 };
55 
56 RESTYPE
57 LRINTNAME(double x)
58 {
59 	u_int32_t i0, i1;
60 	int e, s, shift;
61 	RESTYPE res;
62 
63 	GET_HIGH_WORD(i0, x);
64 	e = i0 >> DBL_FRACHBITS;
65 	s = e >> DBL_EXPBITS;
66 	e = (e & 0x7ff) - DBL_EXP_BIAS;
67 
68 	/* 1.0 x 2^31 (or 2^63) is already too large */
69 	if (e >= (int)RESTYPE_BITS - 1)
70 		return (s ? RESTYPE_MIN : RESTYPE_MAX); /* ??? unspecified */
71 
72 	/* >= 2^52 is already an exact integer */
73 	if (e < DBL_FRACBITS) {
74 		volatile double t = x;	/* clip extra precision */
75 		/* round, using current direction */
76 		t += TWO52[s];
77 		t -= TWO52[s];
78 		x = t;
79 	}
80 
81 	EXTRACT_WORDS(i0, i1, x);
82 	e = ((i0 >> DBL_FRACHBITS) & 0x7ff) - DBL_EXP_BIAS;
83 	i0 &= 0xfffff;
84 	i0 |= (1 << DBL_FRACHBITS);
85 
86 	if (e < 0)
87 		return (0);
88 
89 	shift = e - DBL_FRACBITS;
90 	if (shift >=0)
91 		res = (shift < RESTYPE_BITS ? (RESTYPE)i1 << shift : 0);
92 	else
93 		res = (shift > -RESTYPE_BITS ? (RESTYPE)i1 >> -shift : 0);
94 	shift += 32;
95 	if (shift >=0)
96 		res |= (shift < RESTYPE_BITS ? (RESTYPE)i0 << shift : 0);
97 	else
98 		res |= (shift > -RESTYPE_BITS ? (RESTYPE)i0 >> -shift : 0);
99 
100 	return (s ? -res : res);
101 }
102 
103 #if	LDBL_MANT_DIG == 53
104 #ifdef	lint
105 /* PROTOLIB1 */
106 long int lrintl(long double);
107 #else	/* lint */
108 __weak_alias(lrintl, lrint);
109 #endif	/* lint */
110 #endif	/* LDBL_MANT_DIG == 53 */
111