From 8244f62490dd7cc7fd2ea2d15f770f1b025a8927 Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Mon, 15 Apr 2024 23:10:23 +0530 Subject: [PATCH 1/3] refactor: update math/base/special/log2 to freeBSD v12.2.0 --- .../math/base/special/log2/lib/klog.js | 102 -------------- .../math/base/special/log2/lib/main.js | 36 +++-- .../math/base/special/log2/lib/polyval_p.js | 47 ------- .../math/base/special/log2/lib/polyval_q.js | 47 ------- .../math/base/special/log2/manifest.json | 9 +- .../math/base/special/log2/package.json | 1 - .../base/special/log2/scripts/evalpoly.js | 133 ------------------ .../@stdlib/math/base/special/log2/src/main.c | 120 +++------------- .../math/base/special/log2/test/test.js | 7 + .../math/base/special/log2/test/test.klog.js | 50 ------- .../base/special/log2/test/test.native.js | 7 + 11 files changed, 70 insertions(+), 489 deletions(-) delete mode 100644 lib/node_modules/@stdlib/math/base/special/log2/lib/klog.js delete mode 100644 lib/node_modules/@stdlib/math/base/special/log2/lib/polyval_p.js delete mode 100644 lib/node_modules/@stdlib/math/base/special/log2/lib/polyval_q.js delete mode 100644 lib/node_modules/@stdlib/math/base/special/log2/scripts/evalpoly.js delete mode 100644 lib/node_modules/@stdlib/math/base/special/log2/test/test.klog.js diff --git a/lib/node_modules/@stdlib/math/base/special/log2/lib/klog.js b/lib/node_modules/@stdlib/math/base/special/log2/lib/klog.js deleted file mode 100644 index 6445f028e028..000000000000 --- a/lib/node_modules/@stdlib/math/base/special/log2/lib/klog.js +++ /dev/null @@ -1,102 +0,0 @@ -/** -* @license Apache-2.0 -* -* Copyright (c) 2018 The Stdlib Authors. -* -* 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. -* -* -* ## Notice -* -* The following copyright and license were part of the original implementation available as part of [FreeBSD]{@link https://svnweb.freebsd.org/base/release/9.3.0/lib/msun/src/k_log.h}. The implementation follows the original, but has been modified for JavaScript. -* -* ```text -* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. -* -* Developed at SunPro, a Sun Microsystems, Inc. business. -* Permission to use, copy, modify, and distribute this -* software is freely granted, provided that this notice -* is preserved. -* ``` -*/ - -'use strict'; - -// MODULES // - -var getHighWord = require( '@stdlib/number/float64/base/get-high-word' ); -var polyvalP = require( './polyval_p.js' ); -var polyvalQ = require( './polyval_q.js' ); - - -// VARIABLES // - -// 0x000fffff = 1048575 => 0 00000000000 11111111111111111111 -var HIGH_SIGNIFICAND_MASK = 0x000fffff|0; // asm type annotation - -// 1/3 -var ONE_THIRD = 0.33333333333333333; - - -// MAIN // - -/** -* Return `log(x) - (x-1)` for `x` in `~[sqrt(2)/2, sqrt(2)]`. -* -* @private -* @param {number} x - input value -* @returns {number} function value -*/ -function klog( x ) { - var hfsq; - var t1; - var t2; - var hx; - var f; - var s; - var z; - var R; - var w; - var i; - var j; - - hx = getHighWord( x ); - f = x - 1.0; - if ( ( HIGH_SIGNIFICAND_MASK & (2+hx) ) < 3 ) { - // Case: -2**-20 <= f < 2**-20 - if ( f === 0.0 ) { - return 0.0; - } - return f * f * ( ( ONE_THIRD*f )- 0.5 ); - } - s = f / ( 2.0 + f ); - z = s * s; - hx &= HIGH_SIGNIFICAND_MASK; - i = ( hx - 0x6147a )|0; // asm type annotation - w = z * z; - j = ( 0x6b851 - hx )|0; // asm type annotation - t1 = w * polyvalP( w ); - t2 = z * polyvalQ( w ); - i |= j; - R = t2 + t1; - if ( i > 0 ) { - hfsq = 0.5 * f * f; - return ( s * (hfsq+R) ) - hfsq; - } - return s * (R-f); -} - - -// EXPORTS // - -module.exports = klog; diff --git a/lib/node_modules/@stdlib/math/base/special/log2/lib/main.js b/lib/node_modules/@stdlib/math/base/special/log2/lib/main.js index bf01f01118d1..b0216b38d536 100644 --- a/lib/node_modules/@stdlib/math/base/special/log2/lib/main.js +++ b/lib/node_modules/@stdlib/math/base/special/log2/lib/main.js @@ -18,7 +18,7 @@ * * ## Notice * -* The following copyright and license 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_log2.c}. The implementation follows the original, but has been modified for JavaScript. +* The following copyright and license 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_log2.c}. The implementation follows the original, but has been modified for JavaScript. * * ```text * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. @@ -43,7 +43,7 @@ var ABS_MASK = require( '@stdlib/constants/float64/high-word-abs-mask' ); var HIGH_SIGNIFICAND_MASK = require( '@stdlib/constants/float64/high-word-significand-mask' ); var BIAS = require( '@stdlib/constants/float64/exponent-bias' ); var NINF = require( '@stdlib/constants/float64/ninf' ); -var klog = require( './klog.js' ); +var kernelLog1p = require( '@stdlib/math/base/special/kernel-log1p' ); // VARIABLES // @@ -98,11 +98,17 @@ var WORDS = [ 0|0, 0|0 ]; * // returns NaN */ function log2( x ) { - var hi; - var lo; + var valHi; + var valLo; + var hfsq; var hx; var lx; + var hi; + var lo; var f; + var r; + var w; + var y; var i; var k; @@ -127,6 +133,10 @@ function log2( x ) { if ( hx >= HIGH_MAX_NORMAL_EXP ) { return x + x; } + // Case: log(1) = +0 + if ( hx === HIGH_BIASED_EXP_0 && lx === 0 ) { + return 0; + } k += ( (hx>>20) - BIAS )|0; // asm type annotation hx &= HIGH_SIGNIFICAND_MASK; i = ( ( hx+0x95f64 ) & 0x100000 )|0; // asm type annotation @@ -134,11 +144,19 @@ function log2( x ) { // Normalize x or x/2... x = setHighWord( x, hx|(i^HIGH_BIASED_EXP_0) ); k += (i>>20)|0; // asm type annotation - f = klog( x ); - x -= 1; - hi = setLowWord( x, 0 ); - lo = x - hi; - return ( (x+f)*IVLN2LO ) + ( (lo+f)*IVLN2HI ) + ( hi*IVLN2HI ) + k; + y = k; + f = x - 1.0; + hfsq = 0.5 * f * f; + r = kernelLog1p( f ); + hi = f - hfsq; + hi = setLowWord( hi, 0 ); + lo = ( f - hi ) - hfsq + r; + valHi = hi * IVLN2HI; + valLo = ( ( lo + hi ) * IVLN2LO ) + ( lo * IVLN2HI ); + w = y + valHi; + valLo += ( y - w ) + valHi; + valHi = w; + return valLo + valHi; } diff --git a/lib/node_modules/@stdlib/math/base/special/log2/lib/polyval_p.js b/lib/node_modules/@stdlib/math/base/special/log2/lib/polyval_p.js deleted file mode 100644 index 520e03fbd357..000000000000 --- a/lib/node_modules/@stdlib/math/base/special/log2/lib/polyval_p.js +++ /dev/null @@ -1,47 +0,0 @@ -/** -* @license Apache-2.0 -* -* Copyright (c) 2022 The Stdlib Authors. -* -* 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 is a generated file. Do not edit directly. */ -'use strict'; - -// MAIN // - -/** -* Evaluates a polynomial. -* -* ## Notes -* -* - The implementation uses [Horner's rule][horners-method] for efficient computation. -* -* [horners-method]: https://en.wikipedia.org/wiki/Horner%27s_method -* -* @private -* @param {number} x - value at which to evaluate the polynomial -* @returns {number} evaluated polynomial -*/ -function evalpoly( x ) { - if ( x === 0.0 ) { - return 0.3999999999940942; - } - return 0.3999999999940942 + (x * (0.22222198432149784 + (x * 0.15313837699209373))); // eslint-disable-line max-len -} - - -// EXPORTS // - -module.exports = evalpoly; diff --git a/lib/node_modules/@stdlib/math/base/special/log2/lib/polyval_q.js b/lib/node_modules/@stdlib/math/base/special/log2/lib/polyval_q.js deleted file mode 100644 index fb83e6c31d97..000000000000 --- a/lib/node_modules/@stdlib/math/base/special/log2/lib/polyval_q.js +++ /dev/null @@ -1,47 +0,0 @@ -/** -* @license Apache-2.0 -* -* Copyright (c) 2022 The Stdlib Authors. -* -* 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 is a generated file. Do not edit directly. */ -'use strict'; - -// MAIN // - -/** -* Evaluates a polynomial. -* -* ## Notes -* -* - The implementation uses [Horner's rule][horners-method] for efficient computation. -* -* [horners-method]: https://en.wikipedia.org/wiki/Horner%27s_method -* -* @private -* @param {number} x - value at which to evaluate the polynomial -* @returns {number} evaluated polynomial -*/ -function evalpoly( x ) { - if ( x === 0.0 ) { - return 0.6666666666666735; - } - return 0.6666666666666735 + (x * (0.2857142874366239 + (x * (0.1818357216161805 + (x * 0.14798198605116586))))); // eslint-disable-line max-len -} - - -// EXPORTS // - -module.exports = evalpoly; diff --git a/lib/node_modules/@stdlib/math/base/special/log2/manifest.json b/lib/node_modules/@stdlib/math/base/special/log2/manifest.json index 374f446b52c4..9a4bec85ced9 100644 --- a/lib/node_modules/@stdlib/math/base/special/log2/manifest.json +++ b/lib/node_modules/@stdlib/math/base/special/log2/manifest.json @@ -45,7 +45,8 @@ "@stdlib/constants/float64/high-word-abs-mask", "@stdlib/constants/float64/high-word-significand-mask", "@stdlib/constants/float64/exponent-bias", - "@stdlib/constants/float64/ninf" + "@stdlib/constants/float64/ninf", + "@stdlib/math/base/special/kernel-log1p" ] }, { @@ -67,7 +68,8 @@ "@stdlib/constants/float64/high-word-abs-mask", "@stdlib/constants/float64/high-word-significand-mask", "@stdlib/constants/float64/exponent-bias", - "@stdlib/constants/float64/ninf" + "@stdlib/constants/float64/ninf", + "@stdlib/math/base/special/kernel-log1p" ] }, { @@ -89,7 +91,8 @@ "@stdlib/constants/float64/high-word-abs-mask", "@stdlib/constants/float64/high-word-significand-mask", "@stdlib/constants/float64/exponent-bias", - "@stdlib/constants/float64/ninf" + "@stdlib/constants/float64/ninf", + "@stdlib/math/base/special/kernel-log1p" ] } ] diff --git a/lib/node_modules/@stdlib/math/base/special/log2/package.json b/lib/node_modules/@stdlib/math/base/special/log2/package.json index 8861de28d156..af871cee1f5f 100644 --- a/lib/node_modules/@stdlib/math/base/special/log2/package.json +++ b/lib/node_modules/@stdlib/math/base/special/log2/package.json @@ -21,7 +21,6 @@ "example": "./examples", "include": "./include", "lib": "./lib", - "scripts": "./scripts", "src": "./src", "test": "./test" }, diff --git a/lib/node_modules/@stdlib/math/base/special/log2/scripts/evalpoly.js b/lib/node_modules/@stdlib/math/base/special/log2/scripts/evalpoly.js deleted file mode 100644 index e2f8128ed306..000000000000 --- a/lib/node_modules/@stdlib/math/base/special/log2/scripts/evalpoly.js +++ /dev/null @@ -1,133 +0,0 @@ -/** -* @license Apache-2.0 -* -* Copyright (c) 2018 The Stdlib Authors. -* -* 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 script compiles modules for evaluating polynomial functions. If any polynomial coefficients change, this script should be rerun to update the compiled files. -*/ -'use strict'; - -// MODULES // - -var resolve = require( 'path' ).resolve; -var readFileSync = require( '@stdlib/fs/read-file' ).sync; -var writeFileSync = require( '@stdlib/fs/write-file' ).sync; -var currentYear = require( '@stdlib/time/current-year' ); -var substringBefore = require( '@stdlib/string/substring-before' ); -var substringAfter = require( '@stdlib/string/substring-after' ); -var format = require( '@stdlib/string/format' ); -var licenseHeader = require( '@stdlib/_tools/licenses/header' ); -var compile = require( '@stdlib/math/base/tools/evalpoly-compile' ); -var compileC = require( '@stdlib/math/base/tools/evalpoly-compile-c' ); - - -// VARIABLES // - -// Polynomial coefficients ordered in ascending degree... -var P = [ - 3.999999999940941908e-01, // 3FD99999 9997FA04 - 2.222219843214978396e-01, // 3FCC71C5 1D8E78AF - 1.531383769920937332e-01 // 3FC39A09 D078C69F -]; -var Q = [ - 6.666666666666735130e-01, // 3FE55555 55555593 - 2.857142874366239149e-01, // 3FD24924 94229359 - 1.818357216161805012e-01, // 3FC74664 96CB03DE - 1.479819860511658591e-01 // 3FC2F112 DF3E5244 -]; - -// Header to add to output files: -var header = licenseHeader( 'Apache-2.0', 'js', { - 'year': currentYear(), - 'copyright': 'The Stdlib Authors' -}); -header += '\n/* This is a generated file. Do not edit directly. */\n'; - - -// FUNCTIONS // - -/** -* Inserts a compiled function into file content. -* -* @private -* @param {string} text - source content -* @param {string} id - function identifier -* @param {string} str - function string -* @returns {string} updated content -*/ -function insert( text, id, str ) { - var before; - var after; - var begin; - var end; - - begin = '// BEGIN: '+id; - end = '// END: '+id; - - before = substringBefore( text, begin ); - after = substringAfter( text, end ); - - return format( '%s// BEGIN: %s\n\n%s\n%s%s', before, id, str, end, after ); -} - - -// MAIN // - -/** -* Main execution sequence. -* -* @private -*/ -function main() { - var fpath; - var copts; - var opts; - var file; - var str; - - opts = { - 'encoding': 'utf8' - }; - - fpath = resolve( __dirname, '..', 'lib', 'polyval_p.js' ); - str = header + compile( P ); - writeFileSync( fpath, str, opts ); - - fpath = resolve( __dirname, '..', 'lib', 'polyval_q.js' ); - str = header + compile( Q ); - writeFileSync( fpath, str, opts ); - - copts = { - 'dtype': 'double', - 'name': '' - }; - - fpath = resolve( __dirname, '..', 'src', 'main.c' ); - file = readFileSync( fpath, opts ); - - copts.name = 'polyval_p'; - str = compileC( P, copts ); - file = insert( file, copts.name, str ); - - copts.name = 'polyval_q'; - str = compileC( Q, copts ); - file = insert( file, copts.name, str ); - - writeFileSync( fpath, file, opts ); -} - -main(); diff --git a/lib/node_modules/@stdlib/math/base/special/log2/src/main.c b/lib/node_modules/@stdlib/math/base/special/log2/src/main.c index d0d31140683e..d9e8940caa3e 100644 --- a/lib/node_modules/@stdlib/math/base/special/log2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/log2/src/main.c @@ -26,14 +26,12 @@ #include "stdlib/constants/float64/high_word_significand_mask.h" #include "stdlib/constants/float64/exponent_bias.h" #include "stdlib/constants/float64/ninf.h" +#include "stdlib/math/base/special/kernel_log1p.h" #include // 0x000fffff = 1048575 => 0 00000000000 11111111111111111111 static const int32_t HIGH_SIGNIFICAND_MASK = 0x000fffff; -// 1/3 -static const double ONE_THIRD = 0.33333333333333333; - static const double TWO54 = 1.80143985094819840000e+16; // 0x43500000, 0x00000000 static const double IVLN2HI = 1.44269504072144627571e+00; // 0x3ff71547, 0x65200000 @@ -49,95 +47,6 @@ static const int32_t HIGH_MIN_NORMAL_EXP = 0x00100000; // 0x3ff00000 = 1072693248 => 0 01111111111 00000000000000000000 => biased exponent: 1023 = 0+1023 => 2^0 = 1 static const int32_t HIGH_BIASED_EXP_0 = 0x3ff00000; -/* Begin auto-generated functions. The following functions are auto-generated. Do not edit directly. */ - -// BEGIN: polyval_p - -/** -* Evaluates a polynomial. -* -* ## Notes -* -* - The implementation uses [Horner's rule][horners-method] for efficient computation. -* -* [horners-method]: https://en.wikipedia.org/wiki/Horner%27s_method -* -* @param x value at which to evaluate the polynomial -* @return evaluated polynomial -*/ -static double polyval_p( const double x ) { - return 0.3999999999940942 + (x * (0.22222198432149784 + (x * 0.15313837699209373))); -} - -// END: polyval_p - -// BEGIN: polyval_q - -/** -* Evaluates a polynomial. -* -* ## Notes -* -* - The implementation uses [Horner's rule][horners-method] for efficient computation. -* -* [horners-method]: https://en.wikipedia.org/wiki/Horner%27s_method -* -* @param x value at which to evaluate the polynomial -* @return evaluated polynomial -*/ -static double polyval_q( const double x ) { - return 0.6666666666666735 + (x * (0.2857142874366239 + (x * (0.1818357216161805 + (x * 0.14798198605116586))))); -} - -// END: polyval_q - -/* End auto-generated functions. */ - -/** -* Return log(x) - (x-1) for x in ~[sqrt(2)/2, sqrt(2)]. -* -* @param x input value -* @return output value -*/ -double klog( const double x ) { - double hfsq; - uint32_t hx; - int32_t i; - int32_t j; - double t1; - double t2; - double f; - double R; - double s; - double z; - double w; - - stdlib_base_float64_get_high_word( x, &hx ); - f = x - 1.0; - if ( ( HIGH_SIGNIFICAND_MASK & ( 2 + hx ) ) < 3 ) { - // Case: -2**-20 <= f < 2**-20 - if ( f == 0.0 ) { - return 0.0; - } - return f * f * ( ( ONE_THIRD*f )- 0.5 ); - } - s = f / ( 2.0 + f ); - z = s * s; - hx &= HIGH_SIGNIFICAND_MASK; - i = ( hx - 0x6147a ); - w = z * z; - j = ( 0x6b851 - hx ); - t1 = w * polyval_p( w ); - t2 = z * polyval_q( w ); - i |= j; - R = t2 + t1; - if ( i > 0 ) { - hfsq = 0.5 * f * f; - return ( s * (hfsq+R) ) - hfsq; - } - return s * (R-f); -} - /** * Evaluates the binary logarithm (base two). * @@ -145,14 +54,20 @@ double klog( const double x ) { * @return output value */ double stdlib_base_log2( const double x ) { + double val_lo; + double val_hi; uint32_t hx; uint32_t lx; + double hfsq; double hi; int32_t i; int32_t k; double lo; double xc; double f; + double r; + double w; + double y; xc = x; if ( stdlib_base_is_nan( x ) || x < 0.0 ) { @@ -174,6 +89,10 @@ double stdlib_base_log2( const double x ) { if ( hx >= HIGH_MAX_NORMAL_EXP ) { return xc + xc; } + // Case: log(1) = +0 + if ( hx == HIGH_BIASED_EXP_0 && lx == 0 ) { + return 0; + } k += ( ( hx>>20 ) - STDLIB_CONSTANT_FLOAT64_EXPONENT_BIAS ); hx &= HIGH_SIGNIFICAND_MASK; i = ( ( hx + 0x95f64 ) & 0x100000 ); @@ -181,10 +100,17 @@ double stdlib_base_log2( const double x ) { // Normalize x or x/2... stdlib_base_float64_set_high_word( hx|( i^HIGH_BIASED_EXP_0 ), &xc ); k += (i>>20); - f = klog( xc ); - xc -= 1; - hi = xc; + y = k; + f = xc - 1.0; + hfsq = 0.5 * f * f; + r = stdlib_base_kernel_log1p( f ); + hi = f - hfsq; stdlib_base_float64_set_low_word( 0, &hi ); - lo = xc - hi; - return ( ( xc + f ) * IVLN2LO ) + ( ( lo + f ) * IVLN2HI ) + ( hi * IVLN2HI ) + k; + lo = ( f - hi ) - hfsq + r; + val_lo = hi * IVLN2HI; + val_hi = ( ( lo + hi ) * IVLN2LO ) + ( lo * IVLN2HI ); + w = y + val_lo; + val_hi += ( y - w ) + val_lo; + val_lo = w; + return val_hi + val_lo; } diff --git a/lib/node_modules/@stdlib/math/base/special/log2/test/test.js b/lib/node_modules/@stdlib/math/base/special/log2/test/test.js index f7b283dcdc9f..9dc309dda4b8 100644 --- a/lib/node_modules/@stdlib/math/base/special/log2/test/test.js +++ b/lib/node_modules/@stdlib/math/base/special/log2/test/test.js @@ -26,6 +26,7 @@ var PINF = require( '@stdlib/constants/float64/pinf' ); var NINF = require( '@stdlib/constants/float64/ninf' ); var EPS = require( '@stdlib/constants/float64/eps' ); var abs = require( '@stdlib/math/base/special/abs' ); +var isPositiveZero = require( '@stdlib/math/base/assert/is-positive-zero' ); var log2 = require( './../lib' ); @@ -224,3 +225,9 @@ tape( 'the function returns `NaN` if provided a negative number', function test( t.equal( isnan( v ), true, 'returns NaN' ); t.end(); }); + +tape( 'the function returns `0` if provided a `1`', function test( t ) { + var v = log2( 1.0 ); + t.equal( isPositiveZero( v ), true, 'returns +0' ); + t.end(); +}); diff --git a/lib/node_modules/@stdlib/math/base/special/log2/test/test.klog.js b/lib/node_modules/@stdlib/math/base/special/log2/test/test.klog.js deleted file mode 100644 index 69edb949c4d7..000000000000 --- a/lib/node_modules/@stdlib/math/base/special/log2/test/test.klog.js +++ /dev/null @@ -1,50 +0,0 @@ -/** -* @license Apache-2.0 -* -* Copyright (c) 2018 The Stdlib Authors. -* -* 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. -*/ - -'use strict'; - -// MODULES // - -var tape = require( 'tape' ); -var pow = require( '@stdlib/math/base/special/pow' ); -var klog = require( './../lib/klog.js' ); - - -// TESTS // - -tape( 'main export is a function', function test( t ) { - t.ok( true, __filename ); - t.strictEqual( typeof klog, 'function', 'main export is a function' ); - t.end(); -}); - -tape( 'the function returns `0.0` if provided `x = 1`', function test( t ) { - var y = klog( 1.0 ); - t.equal( y, 0.0, 'returns 0.0' ); - t.end(); -}); - -tape( 'the function correctly computes `log(x) - (x-1)` for `x` satisfying `-2**-20 <= x-1 < 2**-20`', function test( t ) { - var y = klog( 1.0 + pow( 2.0, -20.0 ) ); - t.equal( y, -4.547470617660916e-13, 'returns correct value' ); - - y = klog( 1.0 + pow( -2.0, -20.0 ) ); - t.equal( y, -4.547470617660916e-13, 'returns correct value' ); - - t.end(); -}); diff --git a/lib/node_modules/@stdlib/math/base/special/log2/test/test.native.js b/lib/node_modules/@stdlib/math/base/special/log2/test/test.native.js index 22069ef7cf5a..ffd9963546d9 100644 --- a/lib/node_modules/@stdlib/math/base/special/log2/test/test.native.js +++ b/lib/node_modules/@stdlib/math/base/special/log2/test/test.native.js @@ -27,6 +27,7 @@ var PINF = require( '@stdlib/constants/float64/pinf' ); var NINF = require( '@stdlib/constants/float64/ninf' ); var EPS = require( '@stdlib/constants/float64/eps' ); var abs = require( '@stdlib/math/base/special/abs' ); +var isPositiveZero = require( '@stdlib/math/base/assert/is-positive-zero' ); var tryRequire = require( '@stdlib/utils/try-require' ); @@ -233,3 +234,9 @@ tape( 'the function returns `NaN` if provided a negative number', opts, function t.equal( isnan( v ), true, 'returns NaN' ); t.end(); }); + +tape( 'the function returns `0` if provided a `1`', opts, function test( t ) { + var v = log2( 1.0 ); + t.equal( isPositiveZero( v ), true, 'returns +0' ); + t.end(); +}); From 8160c16f1b940867a5be06692e6db0301d22a5ad Mon Sep 17 00:00:00 2001 From: Athan Reines Date: Mon, 15 Apr 2024 22:17:45 -0700 Subject: [PATCH 2/3] refactor: update implementations and add source comments --- .../math/base/special/log2/lib/main.js | 22 ++++++--- .../@stdlib/math/base/special/log2/src/main.c | 45 +++++++++++-------- .../math/base/special/log2/test/test.js | 2 +- .../base/special/log2/test/test.native.js | 2 +- 4 files changed, 46 insertions(+), 25 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/log2/lib/main.js b/lib/node_modules/@stdlib/math/base/special/log2/lib/main.js index b0216b38d536..7860d7af1f6b 100644 --- a/lib/node_modules/@stdlib/math/base/special/log2/lib/main.js +++ b/lib/node_modules/@stdlib/math/base/special/log2/lib/main.js @@ -106,7 +106,7 @@ function log2( x ) { var hi; var lo; var f; - var r; + var R; var w; var y; var i; @@ -135,11 +135,11 @@ function log2( x ) { } // Case: log(1) = +0 if ( hx === HIGH_BIASED_EXP_0 && lx === 0 ) { - return 0; + return 0.0; } k += ( (hx>>20) - BIAS )|0; // asm type annotation hx &= HIGH_SIGNIFICAND_MASK; - i = ( ( hx+0x95f64 ) & 0x100000 )|0; // asm type annotation + i = ( ( hx+0x95f64 ) & HIGH_MIN_NORMAL_EXP )|0; // asm type annotation // Normalize x or x/2... x = setHighWord( x, hx|(i^HIGH_BIASED_EXP_0) ); @@ -147,15 +147,27 @@ function log2( x ) { y = k; f = x - 1.0; hfsq = 0.5 * f * f; - r = kernelLog1p( f ); + R = kernelLog1p( f ); + + /* + * Notes: + * + * - `f-hfsq` must (for args near `1`) be evaluated in extra precision to avoid a large cancellation when `x` is near `sqrt(2)` or `1/sqrt(2)`.This is fairly efficient since `f-hfsq` only depends on `f`, so can be evaluated in parallel with `R`. Not combining `hfsq` with `R` also keeps `R` small (though not as small as a true `lo` term would be), so that extra precision is not needed for terms involving `R`. + * - When implemented in C, compiler bugs involving extra precision used to break Dekker's theorem for spitting `f-hfsq` as `hi+lo`. These problems are now automatically avoided as a side effect of the optimization of combining the Dekker splitting step with the clear-low-bits step. + * - `y` must (for args near `sqrt(2)` and `1/sqrt(2)`) be added in extra precision to avoid a very large cancellation when `x` is very near these values. Unlike the above cancellations, this problem is specific to base `2`. It is strange that adding `+-1` is so much harder than adding `+-ln2` or `+-log10_2`. + * - This implementation uses Dekker's theorem to normalize `y+val_hi`, so, when implemented in C, compiler bugs may reappear in some configurations. + * - The multi-precision calculations for the multiplications are routine. + */ hi = f - hfsq; hi = setLowWord( hi, 0 ); - lo = ( f - hi ) - hfsq + r; + lo = ( f - hi ) - hfsq + R; valHi = hi * IVLN2HI; valLo = ( ( lo + hi ) * IVLN2LO ) + ( lo * IVLN2HI ); + w = y + valHi; valLo += ( y - w ) + valHi; valHi = w; + return valLo + valHi; } diff --git a/lib/node_modules/@stdlib/math/base/special/log2/src/main.c b/lib/node_modules/@stdlib/math/base/special/log2/src/main.c index d9e8940caa3e..890cf99e9bdc 100644 --- a/lib/node_modules/@stdlib/math/base/special/log2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/log2/src/main.c @@ -29,9 +29,6 @@ #include "stdlib/math/base/special/kernel_log1p.h" #include -// 0x000fffff = 1048575 => 0 00000000000 11111111111111111111 -static const int32_t HIGH_SIGNIFICAND_MASK = 0x000fffff; - static const double TWO54 = 1.80143985094819840000e+16; // 0x43500000, 0x00000000 static const double IVLN2HI = 1.44269504072144627571e+00; // 0x3ff71547, 0x65200000 @@ -54,8 +51,8 @@ static const int32_t HIGH_BIASED_EXP_0 = 0x3ff00000; * @return output value */ double stdlib_base_log2( const double x ) { - double val_lo; - double val_hi; + double valLo; + double valHi; uint32_t hx; uint32_t lx; double hfsq; @@ -65,14 +62,14 @@ double stdlib_base_log2( const double x ) { double lo; double xc; double f; - double r; + double R; double w; double y; - xc = x; if ( stdlib_base_is_nan( x ) || x < 0.0 ) { return 0.0 / 0.0; // NaN } + xc = x; stdlib_base_float64_to_words( xc, &hx, &lx ); k = 0; if ( hx < HIGH_MIN_NORMAL_EXP ) { @@ -94,23 +91,35 @@ double stdlib_base_log2( const double x ) { return 0; } k += ( ( hx>>20 ) - STDLIB_CONSTANT_FLOAT64_EXPONENT_BIAS ); - hx &= HIGH_SIGNIFICAND_MASK; - i = ( ( hx + 0x95f64 ) & 0x100000 ); + hx &= STDLIB_CONSTANT_FLOAT64_HIGH_WORD_SIGNIFICAND_MASK; + i = ( hx+0x95f64 ) & HIGH_MIN_NORMAL_EXP; // Normalize x or x/2... stdlib_base_float64_set_high_word( hx|( i^HIGH_BIASED_EXP_0 ), &xc ); k += (i>>20); - y = k; + y = (double)k; f = xc - 1.0; hfsq = 0.5 * f * f; - r = stdlib_base_kernel_log1p( f ); + R = stdlib_base_kernel_log1p( f ); + + /* + * Notes: + * + * - `f-hfsq` must (for args near `1`) be evaluated in extra precision to avoid a large cancellation when `x` is near `sqrt(2)` or `1/sqrt(2)`.This is fairly efficient since `f-hfsq` only depends on `f`, so can be evaluated in parallel with `R`. Not combining `hfsq` with `R` also keeps `R` small (though not as small as a true `lo` term would be), so that extra precision is not needed for terms involving `R`. + * - Compiler bugs involving extra precision used to break Dekker's theorem for spitting `f-hfsq` as `hi+lo`. These problems are now automatically avoided as a side effect of the optimization of combining the Dekker splitting step with the clear-low-bits step. + * - `y` must (for args near `sqrt(2)` and `1/sqrt(2)`) be added in extra precision to avoid a very large cancellation when `x` is very near these values. Unlike the above cancellations, this problem is specific to base `2`. It is strange that adding `+-1` is so much harder than adding `+-ln2` or `+-log10_2`. + * - This implementation uses Dekker's theorem to normalize `y+val_hi`, so compiler bugs may reappear in some configurations. + * - The multi-precision calculations for the multiplications are routine. + */ hi = f - hfsq; stdlib_base_float64_set_low_word( 0, &hi ); - lo = ( f - hi ) - hfsq + r; - val_lo = hi * IVLN2HI; - val_hi = ( ( lo + hi ) * IVLN2LO ) + ( lo * IVLN2HI ); - w = y + val_lo; - val_hi += ( y - w ) + val_lo; - val_lo = w; - return val_hi + val_lo; + lo = ( f - hi ) - hfsq + R; + valLo = hi * IVLN2HI; + valHi = ( ( lo + hi ) * IVLN2LO ) + ( lo * IVLN2HI ); + + w = y + valLo; + valHi += ( y - w ) + valLo; + valLo = w; + + return valHi + valLo; } diff --git a/lib/node_modules/@stdlib/math/base/special/log2/test/test.js b/lib/node_modules/@stdlib/math/base/special/log2/test/test.js index 9dc309dda4b8..6c38bfd320fe 100644 --- a/lib/node_modules/@stdlib/math/base/special/log2/test/test.js +++ b/lib/node_modules/@stdlib/math/base/special/log2/test/test.js @@ -226,7 +226,7 @@ tape( 'the function returns `NaN` if provided a negative number', function test( t.end(); }); -tape( 'the function returns `0` if provided a `1`', function test( t ) { +tape( 'the function returns positive zero if provided `1.0`', function test( t ) { var v = log2( 1.0 ); t.equal( isPositiveZero( v ), true, 'returns +0' ); t.end(); diff --git a/lib/node_modules/@stdlib/math/base/special/log2/test/test.native.js b/lib/node_modules/@stdlib/math/base/special/log2/test/test.native.js index ffd9963546d9..7bcb2226e842 100644 --- a/lib/node_modules/@stdlib/math/base/special/log2/test/test.native.js +++ b/lib/node_modules/@stdlib/math/base/special/log2/test/test.native.js @@ -235,7 +235,7 @@ tape( 'the function returns `NaN` if provided a negative number', opts, function t.end(); }); -tape( 'the function returns `0` if provided a `1`', opts, function test( t ) { +tape( 'the function returns positive zero if provided `1.0`', opts, function test( t ) { var v = log2( 1.0 ); t.equal( isPositiveZero( v ), true, 'returns +0' ); t.end(); From b0b03c594fee9e26004534343d4bc7c094097716 Mon Sep 17 00:00:00 2001 From: Athan Reines Date: Mon, 15 Apr 2024 22:19:35 -0700 Subject: [PATCH 3/3] style: remove blank lines --- lib/node_modules/@stdlib/math/base/special/log2/src/main.c | 2 -- 1 file changed, 2 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/log2/src/main.c b/lib/node_modules/@stdlib/math/base/special/log2/src/main.c index 890cf99e9bdc..45c03e7ca18f 100644 --- a/lib/node_modules/@stdlib/math/base/special/log2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/log2/src/main.c @@ -30,9 +30,7 @@ #include static const double TWO54 = 1.80143985094819840000e+16; // 0x43500000, 0x00000000 - static const double IVLN2HI = 1.44269504072144627571e+00; // 0x3ff71547, 0x65200000 - static const double IVLN2LO = 1.67517131648865118353e-10; // 0x3de705fc, 0x2eefa200 // 0x7ff00000 = 2146435072 => 0 11111111111 00000000000000000000 => biased exponent: 2047 = 1023+1023 => 2^1023