1178140Sdas/*- 2178140Sdas * Copyright (c) 2004-2008 David Schultz <das@FreeBSD.ORG> 3178140Sdas * All rights reserved. 4178140Sdas * 5178140Sdas * Redistribution and use in source and binary forms, with or without 6178140Sdas * modification, are permitted provided that the following conditions 7178140Sdas * are met: 8178140Sdas * 1. Redistributions of source code must retain the above copyright 9178140Sdas * notice, this list of conditions and the following disclaimer. 10178140Sdas * 2. Redistributions in binary form must reproduce the above copyright 11178140Sdas * notice, this list of conditions and the following disclaimer in the 12178140Sdas * documentation and/or other materials provided with the distribution. 13178140Sdas * 14178140Sdas * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND 15178140Sdas * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 16178140Sdas * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 17178140Sdas * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE 18178140Sdas * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL 19178140Sdas * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS 20178140Sdas * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 21178140Sdas * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT 22178140Sdas * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY 23178140Sdas * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF 24178140Sdas * SUCH DAMAGE. 25178140Sdas */ 26178140Sdas 27178140Sdas#include <sys/cdefs.h> 28178140Sdas__FBSDID("$FreeBSD$"); 29178140Sdas 30178140Sdas#include <float.h> 31178140Sdas#include <limits.h> 32178140Sdas#include <math.h> 33178140Sdas#include <stdint.h> 34178140Sdas 35178140Sdas#ifdef __i386__ 36178140Sdas#include <ieeefp.h> 37178140Sdas#endif 38178140Sdas 39178154Sdas#include "../stdio/floatio.h" 40178140Sdas#include "fpmath.h" 41178140Sdas#include "gdtoaimp.h" 42178140Sdas 43178140Sdas#if (LDBL_MANT_DIG > DBL_MANT_DIG) 44178140Sdas 45178140Sdas/* Strings values used by dtoa() */ 46178140Sdas#define INFSTR "Infinity" 47178140Sdas#define NANSTR "NaN" 48178140Sdas 49178140Sdas#ifdef LDBL_IMPLICIT_NBIT 50178140Sdas#define MANH_SIZE LDBL_MANH_SIZE 51178140Sdas#else 52178140Sdas#define MANH_SIZE (LDBL_MANH_SIZE - 1) 53178140Sdas#endif 54178140Sdas 55178140Sdas#if MANH_SIZE > 32 56178140Sdastypedef uint64_t manh_t; 57178140Sdas#else 58178140Sdastypedef uint32_t manh_t; 59178140Sdas#endif 60178140Sdas 61178140Sdas#if LDBL_MANL_SIZE > 32 62178140Sdastypedef uint64_t manl_t; 63178140Sdas#else 64178140Sdastypedef uint32_t manl_t; 65178140Sdas#endif 66178140Sdas 67178140Sdas#define LDBL_ADJ (LDBL_MAX_EXP - 2) 68178140Sdas#define SIGFIGS ((LDBL_MANT_DIG + 3) / 4 + 1) 69178140Sdas 70178140Sdasstatic const float one[] = { 1.0f, -1.0f }; 71178140Sdas 72178140Sdas/* 73178140Sdas * This is the long double version of __hdtoa(). 74178140Sdas */ 75178140Sdaschar * 76178140Sdas__hldtoa(long double e, const char *xdigs, int ndigits, int *decpt, int *sign, 77178140Sdas char **rve) 78178140Sdas{ 79178140Sdas union IEEEl2bits u; 80178140Sdas char *s, *s0; 81178140Sdas manh_t manh; 82178140Sdas manl_t manl; 83178140Sdas int bufsize; 84178140Sdas#ifdef __i386__ 85178140Sdas fp_prec_t oldprec; 86178140Sdas#endif 87178140Sdas 88178140Sdas u.e = e; 89178140Sdas *sign = u.bits.sign; 90178140Sdas 91178140Sdas switch (fpclassify(e)) { 92178140Sdas case FP_NORMAL: 93178140Sdas *decpt = u.bits.exp - LDBL_ADJ; 94178140Sdas break; 95178140Sdas case FP_ZERO: 96178140Sdas *decpt = 1; 97178140Sdas return (nrv_alloc("0", rve, 1)); 98178140Sdas case FP_SUBNORMAL: 99178140Sdas#ifdef __i386__ 100178140Sdas oldprec = fpsetprec(FP_PE); 101178140Sdas#endif 102178140Sdas u.e *= 0x1p514L; 103178140Sdas *decpt = u.bits.exp - (514 + LDBL_ADJ); 104178140Sdas#ifdef __i386__ 105178140Sdas fpsetprec(oldprec); 106178140Sdas#endif 107178140Sdas break; 108178140Sdas case FP_INFINITE: 109178140Sdas *decpt = INT_MAX; 110178140Sdas return (nrv_alloc(INFSTR, rve, sizeof(INFSTR) - 1)); 111178140Sdas default: /* FP_NAN or unrecognized */ 112178140Sdas *decpt = INT_MAX; 113178140Sdas return (nrv_alloc(NANSTR, rve, sizeof(NANSTR) - 1)); 114178140Sdas } 115178140Sdas 116178140Sdas /* FP_NORMAL or FP_SUBNORMAL */ 117178140Sdas 118178140Sdas if (ndigits == 0) /* dtoa() compatibility */ 119178140Sdas ndigits = 1; 120178140Sdas 121178140Sdas /* 122178140Sdas * If ndigits < 0, we are expected to auto-size, so we allocate 123178140Sdas * enough space for all the digits. 124178140Sdas */ 125178140Sdas bufsize = (ndigits > 0) ? ndigits : SIGFIGS; 126178140Sdas s0 = rv_alloc(bufsize); 127178140Sdas 128178140Sdas /* Round to the desired number of digits. */ 129178140Sdas if (SIGFIGS > ndigits && ndigits > 0) { 130178140Sdas float redux = one[u.bits.sign]; 131178140Sdas int offset = 4 * ndigits + LDBL_MAX_EXP - 4 - LDBL_MANT_DIG; 132178140Sdas#ifdef __i386__ 133178140Sdas oldprec = fpsetprec(FP_PE); 134178140Sdas#endif 135178140Sdas u.bits.exp = offset; 136178140Sdas u.e += redux; 137178140Sdas u.e -= redux; 138178140Sdas *decpt += u.bits.exp - offset; 139178140Sdas#ifdef __i386__ 140178140Sdas fpsetprec(oldprec); 141178140Sdas#endif 142178140Sdas } 143178140Sdas 144178140Sdas mask_nbit_l(u); 145178140Sdas manh = u.bits.manh; 146178140Sdas manl = u.bits.manl; 147178140Sdas *s0 = '1'; 148178140Sdas for (s = s0 + 1; s < s0 + bufsize; s++) { 149178140Sdas *s = xdigs[(manh >> (MANH_SIZE - 4)) & 0xf]; 150178140Sdas manh = (manh << 4) | (manl >> (LDBL_MANL_SIZE - 4)); 151178140Sdas manl <<= 4; 152178140Sdas } 153178140Sdas 154178140Sdas /* If ndigits < 0, we are expected to auto-size the precision. */ 155178140Sdas if (ndigits < 0) { 156178140Sdas for (ndigits = SIGFIGS; s0[ndigits - 1] == '0'; ndigits--) 157178140Sdas ; 158178140Sdas } 159178140Sdas 160178140Sdas s = s0 + ndigits; 161178140Sdas *s = '\0'; 162178140Sdas if (rve != NULL) 163178140Sdas *rve = s; 164178140Sdas return (s0); 165178140Sdas} 166178140Sdas 167178140Sdas#else /* (LDBL_MANT_DIG == DBL_MANT_DIG) */ 168178140Sdas 169178140Sdaschar * 170178140Sdas__hldtoa(long double e, const char *xdigs, int ndigits, int *decpt, int *sign, 171178140Sdas char **rve) 172178140Sdas{ 173178140Sdas 174178140Sdas return (__hdtoa((double)e, xdigs, ndigits, decpt, sign, rve)); 175178140Sdas} 176178140Sdas 177178140Sdas#endif /* (LDBL_MANT_DIG == DBL_MANT_DIG) */ 178