Skip to content

fix: round results near the subnormal limit in math/base/special/ldexp #16297

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
0PrashantYadav0 wants to merge 1 commit into
base: develop
Choose a base branch
from
Open

fix: round results near the subnormal limit in math/base/special/ldexp #16297

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
12 changes: 6 additions & 6 deletions lib/node_modules/@stdlib/math/base/special/ldexp/lib/main.js
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -38,14 +38,14 @@
var copysign = require( '@stdlib/math/base/special/copysign' );
var normalize = require( '@stdlib/number/float64/base/normalize' ).assign;
var floatExp = require( '@stdlib/number/float64/base/exponent' );
var toWords = require( '@stdlib/number/float64/base/to-words' );

Check warning on line 41 in lib/node_modules/@stdlib/math/base/special/ldexp/lib/main.js

View workflow job for this annotation

GitHub Actions / Lint Changed Files 

only the `assign` property of the required module is used; require the property directly
var fromWords = require( '@stdlib/number/float64/base/from-words' );


// VARIABLES //

// 1/(1<<52) = 1/(2**52) = 1/4503599627370496
var TWO52_INV = 2.220446049250313e-16;
// 1/(1<<54) = 1/(2**54) = 1/18014398509481984
var TWO54_INV = 5.551115123125783e-17;

// Exponent all 0s: 1 00000000000 11111111111111111111 => 2148532223
var CLEAR_EXP_MASK = 0x800fffff>>>0; // asm type annotation
Expand Down Expand Up @@ -113,8 +113,8 @@
// Extract the exponent from `frac` and add it to `exp`:
exp += floatExp( frac );

// Check for underflow/overflow...
if ( exp < MIN_SUBNORMAL_EXPONENT ) {
// Check for underflow/overflow (note: values on the interval `[2^-1075, 2^-1074)` may round to the smallest subnormal)...
if ( exp < MIN_SUBNORMAL_EXPONENT-1 ) {
return copysign( 0.0, frac );
}
if ( exp > MAX_EXPONENT ) {
Expand All @@ -125,8 +125,8 @@
}
// Check for a subnormal and scale accordingly to retain precision...
if ( exp <= MAX_SUBNORMAL_EXPONENT ) {
exp += 52;
m = TWO52_INV;
exp += 54;
m = TWO54_INV;
} else {
m = 1.0;
}
Expand Down
12 changes: 6 additions & 6 deletions lib/node_modules/@stdlib/math/base/special/ldexp/src/main.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -32,8 +32,8 @@
#include "stdlib/constants/float64/exponent_bias.h"
#include <stdint.h>

// 1/(1<<52) = 1/(2**52) = 1/4503599627370496
static const double TWO52_INV = 2.220446049250313e-16;
// 1/(1<<54) = 1/(2**54) = 1/18014398509481984
static const double TWO54_INV = 5.551115123125783e-17;

// Exponent all 0s: 1 00000000000 11111111111111111111 => 2148532223
static const uint32_t CLEAR_EXP_MASK = 2148532223;
Expand Down Expand Up @@ -64,8 +64,8 @@ double stdlib_base_ldexp( const double frac, const int32_t exp ) {
// Extract the exponent from `frac` and add it to `exp`:
e += exp + stdlib_base_float64_exponent( y );

// Check for underflow/overflow...
if ( e < STDLIB_CONSTANT_FLOAT64_MIN_BASE2_EXPONENT_SUBNORMAL ) {
// Check for underflow/overflow (note: values on the interval `[2^-1075, 2^-1074)` may round to the smallest subnormal)...
if ( e < STDLIB_CONSTANT_FLOAT64_MIN_BASE2_EXPONENT_SUBNORMAL-1 ) {
return stdlib_base_copysign( 0.0, y );
}
if ( e > STDLIB_CONSTANT_FLOAT64_MAX_BASE2_EXPONENT ) {
Expand All @@ -76,8 +76,8 @@ double stdlib_base_ldexp( const double frac, const int32_t exp ) {
}
// Check for a subnormal and scale accordingly to retain precision...
if ( e <= STDLIB_CONSTANT_FLOAT64_MAX_BASE2_EXPONENT_SUBNORMAL ) {
e += 52;
m = TWO52_INV;
e += 54;
m = TWO54_INV;
} else {
m = 1.0;
}
Expand Down
24 changes: 24 additions & 0 deletions lib/node_modules/@stdlib/math/base/special/ldexp/test/test.js
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -113,6 +113,30 @@ tape( 'the function multiplies a number by an integer power of two (subnormals)'
t.end();
});

tape( 'the function rounds results on the interval `[2^-1075, 2^-1074)` to the nearest representable value', function test( t ) {
var v;

v = ldexp( 0.75, -1074 );
t.strictEqual( v, 5.0e-324, 'returns expected value' );

v = ldexp( -0.75, -1074 );
t.strictEqual( v, -5.0e-324, 'returns expected value' );

v = ldexp( 3.0, -1076 );
t.strictEqual( v, 5.0e-324, 'returns expected value' );

v = ldexp( 0.5000000000000001, -1074 );
t.strictEqual( v, 5.0e-324, 'returns expected value' );

// Ties round to even:
v = ldexp( 0.5, -1074 );
t.strictEqual( isPositiveZero( v ), true, 'returns expected value' );

v = ldexp( -0.5, -1074 );
t.strictEqual( isNegativeZero( v ), true, 'returns expected value' );
t.end();
});

tape( 'if provided a fraction equal to `+0`, the function returns `+0`', function test( t ) {
var v = ldexp( 0.0, 10 );
t.strictEqual( isPositiveZero( v ), true, 'returns expected value' );
Expand Down
24 changes: 24 additions & 0 deletions lib/node_modules/@stdlib/math/base/special/ldexp/test/test.native.js
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -122,6 +122,30 @@ tape( 'the function multiplies a number by an integer power of two (subnormals)'
t.end();
});

tape( 'the function rounds results on the interval `[2^-1075, 2^-1074)` to the nearest representable value', opts, function test( t ) {
var v;

v = ldexp( 0.75, -1074 );
t.strictEqual( v, 5.0e-324, 'returns expected value' );

v = ldexp( -0.75, -1074 );
t.strictEqual( v, -5.0e-324, 'returns expected value' );

v = ldexp( 3.0, -1076 );
t.strictEqual( v, 5.0e-324, 'returns expected value' );

v = ldexp( 0.5000000000000001, -1074 );
t.strictEqual( v, 5.0e-324, 'returns expected value' );

// Ties round to even:
v = ldexp( 0.5, -1074 );
t.strictEqual( isPositiveZero( v ), true, 'returns expected value' );

v = ldexp( -0.5, -1074 );
t.strictEqual( isNegativeZero( v ), true, 'returns expected value' );
t.end();
});

tape( 'if provided a fraction equal to `+0`, the function returns `+0`', opts, function test( t ) {
var v = ldexp( 0.0, 10 );
t.strictEqual( isPositiveZero( v ), true, 'returns expected value' );
Expand Down
2 changes: 1 addition & 1 deletion lib/node_modules/@stdlib/math/base/special/pow/test/test.js
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -389,7 +389,7 @@ tape( 'the function evaluates the exponential function (near underflow)', functi
x = 1.9999999999999998;
y = -1075;
v = pow( x, y );
t.strictEqual( v, 0.0, 'pow('+x+','+y+') returns 0' );
t.strictEqual( v, 5.0e-324, 'pow('+x+','+y+') returns 5e-324' );

x = 1.9999999999999999;
y = -1075;
Expand Down
2 changes: 1 addition & 1 deletion lib/node_modules/@stdlib/math/base/special/pow/test/test.native.js
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -401,7 +401,7 @@ tape( 'the function evaluates the exponential function (near underflow)', opts,
x = 1.9999999999999998;
y = -1075;
v = pow( x, y );
t.strictEqual( v, 0.0, 'pow('+x+','+y+') returns 0' );
t.strictEqual( v, 5.0e-324, 'pow('+x+','+y+') returns 5e-324' );

x = 1.9999999999999999;
y = -1075;
Expand Down
Loading