s_truncl.c revision 145637
1178479Sjb/* 2178479Sjb * ==================================================== 3178479Sjb * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. 4178479Sjb * 5178479Sjb * Developed at SunPro, a Sun Microsystems, Inc. business. 6178479Sjb * Permission to use, copy, modify, and distribute this 7178479Sjb * software is freely granted, provided that this notice 8178479Sjb * is preserved. 9178479Sjb * ==================================================== 10178479Sjb * 11178479Sjb * From: @(#)s_floor.c 5.1 93/09/24 12178479Sjb */ 13178479Sjb 14178479Sjb#ifndef lint 15178479Sjbstatic char rcsid[] = "$FreeBSD: head/lib/msun/src/s_truncl.c 145637 2005-04-28 19:45:55Z stefanf $"; 16178479Sjb#endif 17178479Sjb 18178479Sjb/* 19178479Sjb * truncl(x) 20178479Sjb * Return x rounded toward 0 to integral value 21178479Sjb * Method: 22178479Sjb * Bit twiddling. 23178479Sjb * Exception: 24267941Srpaulo * Inexact flag raised if x not equal to truncl(x). 25267941Srpaulo */ 26178479Sjb 27178479Sjb#include <float.h> 28178479Sjb#include <math.h> 29178479Sjb#include <stdint.h> 30178479Sjb 31178479Sjb#include "fpmath.h" 32178479Sjb 33178479Sjb#ifdef LDBL_IMPLICIT_NBIT 34178479Sjb#define MANH_SIZE (LDBL_MANH_SIZE + 1) 35178479Sjb#else 36178479Sjb#define MANH_SIZE LDBL_MANH_SIZE 37178479Sjb#endif 38178479Sjb 39178479Sjbstatic const long double huge = 1.0e300; 40178479Sjb 41178479Sjblong double 42178479Sjbtruncl(long double x) 43178479Sjb{ 44178479Sjb union IEEEl2bits u = { .e = x }; 45178479Sjb int e = u.bits.exp - LDBL_MAX_EXP + 1; 46178479Sjb 47178479Sjb if (e < MANH_SIZE - 1) { 48178479Sjb if (e < 0) { /* raise inexact if x != 0 */ 49178479Sjb if (huge + x > 0.0) 50178479Sjb u.e = 0.0; 51178479Sjb } else { 52178479Sjb uint64_t m = ((1llu << MANH_SIZE) - 1) >> (e + 1); 53178479Sjb if (((u.bits.manh & m) | u.bits.manl) == 0) 54178479Sjb return (x); /* x is integral */ 55178479Sjb if (huge + x > 0.0) { /* raise inexact flag */ 56178479Sjb u.bits.manh &= ~m; 57178479Sjb u.bits.manl = 0; 58178479Sjb } 59178479Sjb } 60178479Sjb } else if (e < LDBL_MANT_DIG - 1) { 61178479Sjb uint64_t m = (uint64_t)-1 >> (64 - LDBL_MANT_DIG + e + 1); 62178479Sjb if ((u.bits.manl & m) == 0) 63178479Sjb return (x); /* x is integral */ 64178479Sjb if (huge + x > 0.0) /* raise inexact flag */ 65178479Sjb u.bits.manl &= ~m; 66178479Sjb } 67178479Sjb return (u.e); 68178479Sjb} 69178479Sjb