1*0a6a1f1dSLionel Sambuc /* $NetBSD: lrint.c,v 1.5 2015/07/09 06:17:13 nat Exp $ */
22fe8fb19SBen Gras
32fe8fb19SBen Gras /*-
42fe8fb19SBen Gras * Copyright (c) 2004
52fe8fb19SBen Gras * Matthias Drochner. All rights reserved.
62fe8fb19SBen Gras *
72fe8fb19SBen Gras * Redistribution and use in source and binary forms, with or without
82fe8fb19SBen Gras * modification, are permitted provided that the following conditions
92fe8fb19SBen Gras * are met:
102fe8fb19SBen Gras * 1. Redistributions of source code must retain the above copyright
112fe8fb19SBen Gras * notice, this list of conditions and the following disclaimer.
122fe8fb19SBen Gras * 2. Redistributions in binary form must reproduce the above copyright
132fe8fb19SBen Gras * notice, this list of conditions and the following disclaimer in the
142fe8fb19SBen Gras * documentation and/or other materials provided with the distribution.
152fe8fb19SBen Gras *
162fe8fb19SBen Gras * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
172fe8fb19SBen Gras * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
182fe8fb19SBen Gras * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
192fe8fb19SBen Gras * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
202fe8fb19SBen Gras * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
212fe8fb19SBen Gras * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
222fe8fb19SBen Gras * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
232fe8fb19SBen Gras * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
242fe8fb19SBen Gras * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
252fe8fb19SBen Gras * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
262fe8fb19SBen Gras * SUCH DAMAGE.
272fe8fb19SBen Gras */
282fe8fb19SBen Gras
292fe8fb19SBen Gras #include <math.h>
302fe8fb19SBen Gras #include <sys/ieee754.h>
312fe8fb19SBen Gras #include <machine/limits.h>
322fe8fb19SBen Gras #include "math_private.h"
332fe8fb19SBen Gras
342fe8fb19SBen Gras #ifndef LRINTNAME
352fe8fb19SBen Gras #define LRINTNAME lrint
362fe8fb19SBen Gras #define RESTYPE long int
372fe8fb19SBen Gras #define RESTYPE_MIN LONG_MIN
382fe8fb19SBen Gras #define RESTYPE_MAX LONG_MAX
392fe8fb19SBen Gras #endif
402fe8fb19SBen Gras
412fe8fb19SBen Gras #define RESTYPE_BITS (sizeof(RESTYPE) * 8)
422fe8fb19SBen Gras
432fe8fb19SBen Gras static const double
442fe8fb19SBen Gras TWO52[2]={
452fe8fb19SBen Gras 4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
462fe8fb19SBen Gras -4.50359962737049600000e+15, /* 0xC3300000, 0x00000000 */
472fe8fb19SBen Gras };
482fe8fb19SBen Gras
492fe8fb19SBen Gras RESTYPE
LRINTNAME(double x)502fe8fb19SBen Gras LRINTNAME(double x)
512fe8fb19SBen Gras {
522fe8fb19SBen Gras u_int32_t i0, i1;
532fe8fb19SBen Gras int e, s, shift;
542fe8fb19SBen Gras RESTYPE res;
552fe8fb19SBen Gras
562fe8fb19SBen Gras GET_HIGH_WORD(i0, x);
572fe8fb19SBen Gras e = i0 >> 20;
582fe8fb19SBen Gras s = (uint32_t)e >> DBL_EXPBITS;
592fe8fb19SBen Gras e = (e & 0x7ff) - DBL_EXP_BIAS;
602fe8fb19SBen Gras
612fe8fb19SBen Gras /* 1.0 x 2^-1 is the smallest number which can be rounded to 1 */
622fe8fb19SBen Gras if (e < -1)
632fe8fb19SBen Gras return (0);
642fe8fb19SBen Gras /* 1.0 x 2^31 (or 2^63) is already too large */
652fe8fb19SBen Gras if (e >= (int)RESTYPE_BITS - 1)
662fe8fb19SBen Gras return (s ? RESTYPE_MIN : RESTYPE_MAX); /* ??? unspecified */
672fe8fb19SBen Gras
682fe8fb19SBen Gras /* >= 2^52 is already an exact integer */
692fe8fb19SBen Gras if (e < DBL_FRACBITS) {
702fe8fb19SBen Gras /* round, using current direction */
712fe8fb19SBen Gras x += TWO52[s];
722fe8fb19SBen Gras x -= TWO52[s];
73*0a6a1f1dSLionel Sambuc } else
74*0a6a1f1dSLionel Sambuc return x;
752fe8fb19SBen Gras
762fe8fb19SBen Gras EXTRACT_WORDS(i0, i1, x);
772fe8fb19SBen Gras e = ((i0 >> 20) & 0x7ff) - DBL_EXP_BIAS;
782fe8fb19SBen Gras i0 &= 0xfffff;
792fe8fb19SBen Gras i0 |= (1 << 20);
802fe8fb19SBen Gras
812fe8fb19SBen Gras shift = e - DBL_FRACBITS;
822fe8fb19SBen Gras if (shift >=0)
832fe8fb19SBen Gras res = (shift < 32 ? (RESTYPE)i1 << shift : 0);
842fe8fb19SBen Gras else
852fe8fb19SBen Gras res = (shift > -32 ? i1 >> -shift : 0);
862fe8fb19SBen Gras shift += 32;
872fe8fb19SBen Gras if (shift >=0)
882fe8fb19SBen Gras res |= (shift < 32 ? (RESTYPE)i0 << shift : 0);
892fe8fb19SBen Gras else
902fe8fb19SBen Gras res |= (shift > -32 ? i0 >> -shift : 0);
912fe8fb19SBen Gras
922fe8fb19SBen Gras return (s ? -res : res);
932fe8fb19SBen Gras }
94