diff --git a/lib/node_modules/@stdlib/math/base/special/atanh/lib/main.js b/lib/node_modules/@stdlib/math/base/special/atanh/lib/main.js index ae3f302923e4..feb81eac2d9d 100644 --- a/lib/node_modules/@stdlib/math/base/special/atanh/lib/main.js +++ b/lib/node_modules/@stdlib/math/base/special/atanh/lib/main.js @@ -18,7 +18,7 @@ * * ## Notice * -* The following copyright, license, and long comment were part of the original implementation available as part of [FreeBSD]{@link https://svnweb.freebsd.org/base/release/9.3.0/lib/msun/src/e_atanh.c?view=markup}. The implementation follows the original, but has been modified for JavaScript. +* The following copyright, license, and long comment were part of the original implementation available as part of [FreeBSD]{@link https://svnweb.freebsd.org/base/release/12.2.0/lib/msun/src/e_atanh.c?view=markup}. The implementation follows the original, but has been modified for JavaScript. * * ```text * Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved. diff --git a/lib/node_modules/@stdlib/math/base/special/atanh/manifest.json b/lib/node_modules/@stdlib/math/base/special/atanh/manifest.json index 4678a491a177..2a6733339c24 100644 --- a/lib/node_modules/@stdlib/math/base/special/atanh/manifest.json +++ b/lib/node_modules/@stdlib/math/base/special/atanh/manifest.json @@ -39,6 +39,9 @@ "@stdlib/math/base/napi/unary", "@stdlib/math/base/special/log1p", "@stdlib/math/base/assert/is-nan", + "@stdlib/constants/float64/high-word-abs-mask", + "@stdlib/number/float64/base/set-high-word", + "@stdlib/number/float64/base/to-words", "@stdlib/constants/float64/pinf", "@stdlib/constants/float64/ninf" ] @@ -56,6 +59,9 @@ "dependencies": [ "@stdlib/math/base/special/log1p", "@stdlib/math/base/assert/is-nan", + "@stdlib/constants/float64/high-word-abs-mask", + "@stdlib/number/float64/base/set-high-word", + "@stdlib/number/float64/base/to-words", "@stdlib/constants/float64/pinf", "@stdlib/constants/float64/ninf" ] @@ -73,6 +79,9 @@ "dependencies": [ "@stdlib/math/base/special/log1p", "@stdlib/math/base/assert/is-nan", + "@stdlib/constants/float64/high-word-abs-mask", + "@stdlib/number/float64/base/set-high-word", + "@stdlib/number/float64/base/to-words", "@stdlib/constants/float64/pinf", "@stdlib/constants/float64/ninf" ] diff --git a/lib/node_modules/@stdlib/math/base/special/atanh/src/main.c b/lib/node_modules/@stdlib/math/base/special/atanh/src/main.c index 954b9f671f17..a932542d7233 100644 --- a/lib/node_modules/@stdlib/math/base/special/atanh/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/atanh/src/main.c @@ -33,6 +33,9 @@ #include "stdlib/math/base/special/atanh.h" #include "stdlib/math/base/assert/is_nan.h" #include "stdlib/math/base/special/log1p.h" +#include "stdlib/constants/float64/high_word_abs_mask.h" +#include "stdlib/number/float64/base/set_high_word.h" +#include "stdlib/number/float64/base/to_words.h" #include "stdlib/constants/float64/pinf.h" #include "stdlib/constants/float64/ninf.h" #include <stdint.h> @@ -79,34 +82,48 @@ static const double NEAR_ZERO = 1.0 / (1 << 28); // 2**-28 * // returns ~1.472 */ double stdlib_base_atanh( const double x ) { - int32_t sgn; + uint32_t ahx; + uint32_t lx; + uint32_t hx; double ax; double t; + if ( stdlib_base_is_nan( x ) || x < -1.0 || x > 1.0 ) { return 0.0 / 0.0; // NaN } + + // Split `x` into high and low words: + stdlib_base_float64_to_words( x, &hx, &lx ); + + ahx = hx & STDLIB_CONSTANT_FLOAT64_HIGH_WORD_ABS_MASK; + + // special cases if ( x == 1.0 ) { return STDLIB_CONSTANT_FLOAT64_PINF; } if ( x == -1.0 ) { return STDLIB_CONSTANT_FLOAT64_NINF; } + if ( x < 0.0 ) { - sgn = 1; + hx = 1; ax = -x; } else { - sgn = 0; + hx = 0; ax = x; } // Case: |x| < 2**-28 if ( ax < NEAR_ZERO ) { - return ( sgn == 1 ) ? -ax : ax; + return ( hx == 1 ) ? -ax : ax; } - if ( ax < 0.5 ) { + + stdlib_base_float64_set_high_word( ahx, &ax ); + + if( ax < 0.5 ) { t = ax + ax; - t = 0.5 * stdlib_base_log1p( t + ( t * ax / ( 1 - ax ) ) ); + t = 0.5 * stdlib_base_log1p( t + ( t * ax / ( 1.0 - ax ) ) ); } else { - t = 0.5 * stdlib_base_log1p( ( ax + ax ) / ( 1 - ax ) ); + t = 0.5 * stdlib_base_log1p( ( ax + ax ) / ( 1.0 - ax ) ); } - return ( sgn == 1 ) ? -t : t; + return ( hx == 1 ) ? -t : t; }