From 2a5468a2ac65cff2e132d9547a1fbbdbafa100ca Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A1szl=C3=B3=20Lang=C3=B3?= Date: Tue, 3 May 2016 12:40:46 +0200 Subject: [PATCH] Improve double to string conversion MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Added errol implementation. Fixed '3d-raytrace.js', because the new algorithm has better precision and also reduce the code size. JerryScript-DCO-1.0-Signed-off-by: László Langó llango.u-szeged@partner.samsung.com --- .../ecma/base/ecma-helpers-conversion.c | 205 +-------------- jerry-core/ecma/base/ecma-helpers-errol.c | 238 ++++++++++++++++++ jerry-core/ecma/base/ecma-helpers.h | 5 + jerry-libm/include/math.h | 2 + jerry-libm/jerry-libm-internal.h | 2 + jerry-libm/nextafter.c | 135 ++++++++++ 6 files changed, 386 insertions(+), 201 deletions(-) create mode 100644 jerry-core/ecma/base/ecma-helpers-errol.c create mode 100644 jerry-libm/nextafter.c diff --git a/jerry-core/ecma/base/ecma-helpers-conversion.c b/jerry-core/ecma/base/ecma-helpers-conversion.c index 0639f8114..da7a2856a 100644 --- a/jerry-core/ecma/base/ecma-helpers-conversion.c +++ b/jerry-core/ecma/base/ecma-helpers-conversion.c @@ -14,6 +14,8 @@ * limitations under the License. */ +#include + #include "ecma-globals.h" #include "ecma-helpers.h" #include "jrt-libc-includes.h" @@ -916,74 +918,6 @@ ecma_number_to_int32 (ecma_number_t num) /**< ecma-number */ return ret; } /* ecma_number_to_int32 */ -#if CONFIG_ECMA_NUMBER_TYPE == CONFIG_ECMA_NUMBER_FLOAT64 -/** - * Perform conversion of 128-bit binary representation of number - * to decimal representation with decimal exponent. - */ -static void -ecma_number_helper_binary_to_decimal (ECMA_NUMBER_CONVERSION_128BIT_INTEGER_ARG (fraction_uint128), /**< mantissa */ - int32_t binary_exponent, /**< binary exponent */ - int32_t *out_decimal_exp_p) /**< [out] decimal exponent */ -{ - int32_t decimal_exp = 0; - - if (binary_exponent > 0) - { - while (binary_exponent > 0) - { - if (!ECMA_NUMBER_CONVERSION_128BIT_INTEGER_IS_HIGH_BIT_MASK_ZERO (fraction_uint128, 124)) - { - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_INC (fraction_uint128); - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_RIGHT_SHIFT (fraction_uint128); - binary_exponent++; - } - else - { - ECMA_NUMBER_CONVERSION_128BIT_INTEGER (fraction_uint128_tmp); - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_COPY (fraction_uint128_tmp, fraction_uint128); - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_DIV_10 (fraction_uint128_tmp); - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_MUL_10 (fraction_uint128_tmp); - - if (!ECMA_NUMBER_CONVERSION_128BIT_INTEGER_ARE_EQUAL (fraction_uint128, fraction_uint128_tmp) - && ECMA_NUMBER_CONVERSION_128BIT_INTEGER_IS_HIGH_BIT_MASK_ZERO (fraction_uint128, 123)) - { - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_LEFT_SHIFT (fraction_uint128); - binary_exponent--; - } - else - { - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_DIV_10 (fraction_uint128); - decimal_exp++; - } - } - } - } - else if (binary_exponent < 0) - { - while (binary_exponent < 0) - { - if (ECMA_NUMBER_CONVERSION_128BIT_INTEGER_IS_LOW_BIT_MASK_ZERO (fraction_uint128, 0) - || !ECMA_NUMBER_CONVERSION_128BIT_INTEGER_IS_HIGH_BIT_MASK_ZERO (fraction_uint128, 124)) - { - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_RIGHT_SHIFT (fraction_uint128); - - binary_exponent++; - } - else - { - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_MUL_10 (fraction_uint128); - - decimal_exp--; - } - } - } - - *out_decimal_exp_p = decimal_exp; -} /* ecma_number_helper_binary_to_decimal */ - -#endif /* CONFIG_ECMA_NUMBER_TYPE == CONFIG_ECMA_NUMBER_FLOAT64 */ - /** * Perform conversion of ecma-number to decimal representation with decimal exponent * @@ -1002,142 +936,11 @@ ecma_number_to_decimal (ecma_number_t num, /**< ecma-number */ JERRY_ASSERT (!ecma_number_is_nan (num)); JERRY_ASSERT (!ecma_number_is_zero (num)); JERRY_ASSERT (!ecma_number_is_infinity (num)); + JERRY_ASSERT (!ecma_number_is_negative (num)); #if CONFIG_ECMA_NUMBER_TYPE == CONFIG_ECMA_NUMBER_FLOAT64 - ecma_number_t num_m1 = ecma_number_get_prev (num); - ecma_number_t num_p1 = ecma_number_get_next (num); - ECMA_NUMBER_CONVERSION_128BIT_INTEGER (fraction_uint128); - ECMA_NUMBER_CONVERSION_128BIT_INTEGER (fraction_uint128_m1); - ECMA_NUMBER_CONVERSION_128BIT_INTEGER (fraction_uint128_p1); - - uint64_t fraction_uint64, fraction_uint64_m1, fraction_uint64_p1; - int32_t binary_exponent, binary_exponent_m1, binary_exponent_p1; - int32_t decimal_exp, decimal_exp_m1, decimal_exp_p1; - int32_t dot_shift, dot_shift_m1, dot_shift_p1; - - dot_shift_m1 = ecma_number_get_fraction_and_exponent (num_m1, &fraction_uint64_m1, &binary_exponent_m1); - dot_shift = ecma_number_get_fraction_and_exponent (num, &fraction_uint64, &binary_exponent); - dot_shift_p1 = ecma_number_get_fraction_and_exponent (num_p1, &fraction_uint64_p1, &binary_exponent_p1); - - binary_exponent_m1 -= dot_shift_m1; - binary_exponent -= dot_shift; - binary_exponent_p1 -= dot_shift_p1; - - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_INIT (fraction_uint128, - 0ull, - 0ull, - (fraction_uint64) >> 32u, - ((fraction_uint64) << 32u) >> 32u); - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_INIT (fraction_uint128_m1, - 0ull, - 0ull, - (fraction_uint64_m1) >> 32u, - ((fraction_uint64_m1) << 32u) >> 32u); - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_INIT (fraction_uint128_p1, - 0ull, - 0ull, - (fraction_uint64_p1) >> 32u, - ((fraction_uint64_p1) << 32u) >> 32u); - - ecma_number_helper_binary_to_decimal (fraction_uint128, binary_exponent, &decimal_exp); - ecma_number_helper_binary_to_decimal (fraction_uint128_m1, binary_exponent_m1, &decimal_exp_m1); - ecma_number_helper_binary_to_decimal (fraction_uint128_p1, binary_exponent_p1, &decimal_exp_p1); - - if (ECMA_NUMBER_CONVERSION_128BIT_INTEGER_IS_ZERO (fraction_uint128_m1)) - { - decimal_exp_m1 = decimal_exp; - } - - while (decimal_exp != decimal_exp_m1 - || decimal_exp != decimal_exp_p1) - { - while (decimal_exp > decimal_exp_m1 - || decimal_exp > decimal_exp_p1) - { - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_MUL_10 (fraction_uint128); - decimal_exp--; - } - while (decimal_exp_m1 > decimal_exp - || decimal_exp_m1 > decimal_exp_p1) - { - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_MUL_10 (fraction_uint128_m1); - decimal_exp_m1--; - } - while (decimal_exp_p1 > decimal_exp - || decimal_exp_p1 > decimal_exp_m1) - { - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_MUL_10 (fraction_uint128_p1); - decimal_exp_p1--; - } - } - - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_ADD (fraction_uint128_m1, fraction_uint128); - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_RIGHT_SHIFT (fraction_uint128_m1); - - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_ADD (fraction_uint128_p1, fraction_uint128); - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_RIGHT_SHIFT (fraction_uint128_p1); - - /* While fraction doesn't fit to integer, divide it by 10 - and simultaneously increment decimal exponent */ - uint64_t digits_min, digits_max; - - while (!ECMA_NUMBER_CONVERSION_128BIT_INTEGER_IS_HIGH_BIT_MASK_ZERO (fraction_uint128_m1, 63)) - { - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_DIV_10 (fraction_uint128_m1); - decimal_exp_m1++; - } - while (!ECMA_NUMBER_CONVERSION_128BIT_INTEGER_IS_HIGH_BIT_MASK_ZERO (fraction_uint128_p1, 63)) - { - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_DIV_10 (fraction_uint128_p1); - decimal_exp_p1++; - } - - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_ROUND_MIDDLE_AND_LOW_TO_UINT64 (fraction_uint128_m1, digits_min); - ECMA_NUMBER_CONVERSION_128BIT_INTEGER_ROUND_MIDDLE_AND_LOW_TO_UINT64 (fraction_uint128_p1, digits_max); - - digits_min++; - - if (decimal_exp_m1 < decimal_exp_p1) - { - JERRY_ASSERT (decimal_exp_m1 == decimal_exp_p1 - 1); - - digits_min /= 10; - decimal_exp_m1++; - } - else if (decimal_exp_m1 > decimal_exp_p1) - { - JERRY_ASSERT (decimal_exp_m1 == decimal_exp_p1 + 1); - - digits_max /= 10; - decimal_exp_p1++; - } - - JERRY_ASSERT (digits_max >= digits_min); - - while (digits_min / 10 != digits_max / 10) - { - digits_min /= 10; - digits_max /= 10; - decimal_exp_m1++; - decimal_exp_p1++; - } - - uint64_t digits = (digits_min + digits_max + 1) / 2; - int32_t digits_num = 0; - uint64_t t = digits; - - while (t != 0) - { - t /= 10; - digits_num++; - } - - JERRY_ASSERT (digits_num > 0); - - *out_digits_p = digits; - *out_digits_num_p = digits_num; - *out_decimal_exp_p = decimal_exp_p1 + digits_num; + *out_digits_p = ecma_errol0_dtoa (num, out_digits_num_p, out_decimal_exp_p); #elif CONFIG_ECMA_NUMBER_TYPE == CONFIG_ECMA_NUMBER_FLOAT32 /* Less precise conversion */ diff --git a/jerry-core/ecma/base/ecma-helpers-errol.c b/jerry-core/ecma/base/ecma-helpers-errol.c new file mode 100644 index 000000000..149f8182b --- /dev/null +++ b/jerry-core/ecma/base/ecma-helpers-errol.c @@ -0,0 +1,238 @@ +/* Copyright 2016 Samsung Electronics Co., Ltd. + * Copyright 2016 University of Szeged + * + * Licensed under the Apache License, Version 2.0 (the "License"); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0 + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + * + * This file is based on work under the following copyright and permission + * notice: + * + * Copyright (c) 2016 Marc Andrysco + * + * Permission is hereby granted, free of charge, to any person obtaining a copy + * of this software and associated documentation files (the "Software"), to deal + * in the Software without restriction, including without limitation the rights + * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell + * copies of the Software, and to permit persons to whom the Software is + * furnished to do so, subject to the following conditions: + * + * The above copyright notice and this permission notice shall be included in all + * copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR + * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, + * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE + * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER + * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, + * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include + +#include "config.h" +#include "ecma-helpers.h" + +/** \addtogroup ecma ECMA + * @{ + * + * \addtogroup ecmahelpers Helpers for operations with ECMA data types + * @{ + */ + +#if CONFIG_ECMA_NUMBER_TYPE == CONFIG_ECMA_NUMBER_FLOAT64 + +/** + * Printing Floating-Point Numbers + * + * available at http://cseweb.ucsd.edu/~mandrysc/pub/dtoa.pdf + */ + +/** + * Floating point format definitions + */ +#define ECMA_NEXT_FLOAT(value) (nextafter ((value), INFINITY)) +#define ECMA_PREV_FLOAT(value) (nextafter ((value), -INFINITY)) + +#define ERROL0_EPSILON 0.0000001 + +/** + * High-precision data structure. + */ +typedef struct +{ + ecma_number_t value; /**< value */ + ecma_number_t offset; /**< offset */ +} ecma_high_prec_t; + +/** + * Normalize the number by factoring in the error. + */ +static inline void __attr_always_inline___ +ecma_normalize_high_prec_data (ecma_high_prec_t *hp_data_p) /**< [in, out] float pair */ +{ + ecma_number_t val = hp_data_p->value; + + hp_data_p->value += hp_data_p->offset; + hp_data_p->offset += val - hp_data_p->value; +} /* ecma_normalize_high_prec_data */ + +/** + * Multiply the high-precision number by ten. + */ +static inline void __attr_always_inline___ +ecma_multiply_high_prec_by_10 (ecma_high_prec_t *hp_data_p) /**< [in, out] high-precision number */ +{ + ecma_number_t value = hp_data_p->value; + + hp_data_p->value *= 10.0; + hp_data_p->offset *= 10.0; + + ecma_number_t offset = hp_data_p->value; + + offset -= value * 8.0; + offset -= value * 2.0; + + hp_data_p->offset -= offset; + + ecma_normalize_high_prec_data (hp_data_p); +} /* ecma_multiply_high_prec_by_10 */ + +/** + * Divide the high-precision number by ten. + */ +static void +ecma_divide_high_prec_by_10 (ecma_high_prec_t *hp_data_p) /**< [in, out] high-precision number */ +{ + ecma_number_t value = hp_data_p->value; + + hp_data_p->value /= 10.0; + hp_data_p->offset /= 10.0; + + value -= hp_data_p->value * 8.0; + value -= hp_data_p->value * 2.0; + + hp_data_p->offset += value / 10.0; + + ecma_normalize_high_prec_data (hp_data_p); +} /* ecma_divide_high_prec_by_10 */ + +/** + * Errol0 double to ASCII conversion, guaranteed correct but possibly not optimal. + * + * @return digits + */ +inline uint64_t __attr_always_inline___ +ecma_errol0_dtoa (ecma_number_t val, /**< ecma number */ + int32_t *num_of_digits_p, /**< [out] number of digits */ + int32_t *exp_p) /**< [out] exponent */ +{ + uint64_t digits = 0u; + int32_t num_of_digits = 0; + double power_of_10 = 1.0; + int32_t exp = 1; + + /* normalize the midpoint */ + ecma_high_prec_t mid; + + mid.value = val; + mid.offset = 0.0; + + while (((mid.value > 10.0) || ((mid.value == 10.0) && (mid.offset >= 0.0))) && (exp < 308)) + { + exp++; + ecma_divide_high_prec_by_10 (&mid); + power_of_10 /= 10.0; + } + + while (((mid.value < 1.0) || ((mid.value == 1.0) && (mid.offset < 0.0))) && (exp > -307)) + { + exp--; + ecma_multiply_high_prec_by_10 (&mid); + power_of_10 *= 10.0; + } + + ecma_high_prec_t high_bound, low_bound; + + high_bound.value = mid.value; + high_bound.offset = mid.offset + (ECMA_NEXT_FLOAT (val) - val) * power_of_10 / (2.0 + ERROL0_EPSILON); + low_bound.value = mid.value; + low_bound.offset = mid.offset + (ECMA_PREV_FLOAT (val) - val) * power_of_10 / (2.0 + ERROL0_EPSILON); + + ecma_normalize_high_prec_data (&high_bound); + ecma_normalize_high_prec_data (&low_bound); + + /* normalized boundaries */ + + while (high_bound.value > 10.0 || (high_bound.value == 10.0 && (high_bound.offset >= 0.0))) + { + exp++; + ecma_divide_high_prec_by_10 (&high_bound); + ecma_divide_high_prec_by_10 (&low_bound); + } + + while (high_bound.value < 1.0 || (high_bound.value == 1.0 && (high_bound.offset < 0.0))) + { + exp--; + ecma_multiply_high_prec_by_10 (&high_bound); + ecma_multiply_high_prec_by_10 (&low_bound); + } + + /* digit generation */ + + while (high_bound.value != 0.0 || high_bound.offset != 0.0) + { + uint8_t high_digit = (uint8_t) high_bound.value; + + if ((high_bound.value == high_digit) && (high_bound.offset < 0)) + { + high_digit = (uint8_t) (high_digit - 1u); + } + + uint8_t low_digit = (uint8_t) low_bound.value; + + if ((low_bound.value == low_digit) && (low_bound.offset < 0)) + { + low_digit = (uint8_t) (low_digit - 1u); + } + + if (low_digit != high_digit) + { + break; + } + + digits *= 10; + digits += (uint64_t) high_digit; + num_of_digits++; + + high_bound.value -= high_digit; + ecma_multiply_high_prec_by_10 (&high_bound); + + low_bound.value -= low_digit; + ecma_multiply_high_prec_by_10 (&low_bound); + } + + double mdig = (high_bound.value + low_bound.value) / 2.0 + 0.5; + + *num_of_digits_p = num_of_digits + 1; + *exp_p = exp; + digits *= 10; + + return digits + (uint64_t) mdig; +} /* ecma_errol0_dtoa */ + +#endif /* CONFIG_ECMA_NUMBER_TYPE == CONFIG_ECMA_NUMBER_FLOAT64 */ + +/** + * @} + * @} + */ diff --git a/jerry-core/ecma/base/ecma-helpers.h b/jerry-core/ecma/base/ecma-helpers.h index 36257949a..40807e313 100644 --- a/jerry-core/ecma/base/ecma-helpers.h +++ b/jerry-core/ecma/base/ecma-helpers.h @@ -308,6 +308,11 @@ extern uint32_t ecma_number_to_uint32 (ecma_number_t); extern int32_t ecma_number_to_int32 (ecma_number_t); extern lit_utf8_size_t ecma_number_to_utf8_string (ecma_number_t, lit_utf8_byte_t *, lit_utf8_size_t); +/* ecma-helpers-errol.c */ +#if CONFIG_ECMA_NUMBER_TYPE == CONFIG_ECMA_NUMBER_FLOAT64 +extern uint64_t ecma_errol0_dtoa (ecma_number_t, int32_t *, int32_t *); +#endif /* CONFIG_ECMA_NUMBER_TYPE == CONFIG_ECMA_NUMBER_FLOAT64 */ + /** * @} * @} diff --git a/jerry-libm/include/math.h b/jerry-libm/include/math.h index fc714432c..53dc17d30 100644 --- a/jerry-libm/include/math.h +++ b/jerry-libm/include/math.h @@ -76,6 +76,8 @@ double floor (double); double fabs (double); double fmod (double, double); +double nextafter (double, double); + #ifdef __cplusplus } #endif /* __cplusplus */ diff --git a/jerry-libm/jerry-libm-internal.h b/jerry-libm/jerry-libm-internal.h index a30f421f4..3551f1441 100644 --- a/jerry-libm/jerry-libm-internal.h +++ b/jerry-libm/jerry-libm-internal.h @@ -72,6 +72,8 @@ extern double fmod (double, double); extern int isnan (double); extern int finite (double); +double nextafter (double, double); + /* * Functions callable from C, intended to support IEEE arithmetic. */ diff --git a/jerry-libm/nextafter.c b/jerry-libm/nextafter.c new file mode 100644 index 000000000..a6b8fabe1 --- /dev/null +++ b/jerry-libm/nextafter.c @@ -0,0 +1,135 @@ +/* Copyright 2016 Samsung Electronics Co., Ltd. + * Copyright 2016 University of Szeged + * + * Licensed under the Apache License, Version 2.0 (the "License"); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0 + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + * + * This file is based on work under the following copyright and permission + * notice: + * + * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. + * + * Developed at SunSoft, a Sun Microsystems, Inc. business. + * Permission to use, copy, modify, and distribute this + * software is freely granted, provided that this notice + * is preserved. + * + * @(#)s_nextafter.c 1.3 95/01/18 + */ + +#include "jerry-libm-internal.h" + +double +nextafter (double x, + double y) +{ + int hx, hy, ix, iy; + unsigned lx, ly; + + hx = __HI (x); /* high word of x */ + lx = __LO (x); /* low word of x */ + hy = __HI (y); /* high word of y */ + ly = __LO (y); /* low word of y */ + ix = hx & 0x7fffffff; /* |x| */ + iy = hy & 0x7fffffff; /* |y| */ + + if (((ix >= 0x7ff00000) && ((ix - 0x7ff00000) | lx) != 0) /* x is nan */ + || ((iy >= 0x7ff00000) && ((iy - 0x7ff00000) | ly) != 0)) /* y is nan */ + { + return x + y; + } + + if (x == y) + { + return x; /* x=y, return x */ + } + + if ((ix | lx) == 0) + { /* x == 0 */ + __HI (x) = hy & 0x80000000; /* return +-minsubnormal */ + __LO (x) = 1; + y = x * x; + if (y == x) + { + return y; + } + else + { + return x; /* raise underflow flag */ + } + } + + if (hx >= 0) + { /* x > 0 */ + if (hx > hy || ((hx == hy) && (lx > ly))) + { /* x > y, x -= ulp */ + if (lx == 0) + { + hx -= 1; + } + + lx -= 1; + } + else + { /* x < y, x += ulp */ + lx += 1; + + if (lx == 0) + { + hx += 1; + } + } + } + else + { /* x < 0 */ + if (hy >= 0 || hx > hy || ((hx == hy) && (lx > ly))) + { /* x < y, x -= ulp */ + if (lx == 0) + { + hx -= 1; + } + + lx -= 1; + } + else + { /* x > y, x += ulp */ + lx += 1; + + if (lx == 0) + { + hx += 1; + } + } + } + + hy = hx & 0x7ff00000; + + if (hy >= 0x7ff00000) + { + return x + x; /* overflow */ + } + + if (hy < 0x00100000) + { /* underflow */ + y = x * x; + if (y != x) + { /* raise underflow flag */ + __HI (y) = hx; + __LO (y) = lx; + return y; + } + } + + __HI (x) = hx; + __LO (x) = lx; + return x; +} /* nextafter */