From 7e8a907a3a6f72da5bf8560a99b38c9d1e63b39c Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 17:38:10 +0530 Subject: [PATCH 01/28] feat: draft commit --- .../math/base/special/rempio2/README.md | 91 +++ .../rempio2/benchmark/benchmark.native.js | 63 +++ .../rempio2/benchmark/c/native/Makefile | 146 +++++ .../rempio2/benchmark/c/native/benchmark.c | 138 +++++ .../math/base/special/rempio2/binding.gyp | 170 ++++++ .../base/special/rempio2/examples/c/Makefile | 146 +++++ .../base/special/rempio2/examples/c/example.c | 33 ++ .../math/base/special/rempio2/include.gypi | 53 ++ .../stdlib/math/base/special/rempio2.h | 41 ++ .../math/base/special/rempio2/lib/native.js | 56 ++ .../math/base/special/rempio2/manifest.json | 86 +++ .../math/base/special/rempio2/src/Makefile | 70 +++ .../math/base/special/rempio2/src/addon.c | 119 ++++ .../math/base/special/rempio2/src/main.c | 530 ++++++++++++++++++ .../base/special/rempio2/test/test.native.js | 281 ++++++++++ 15 files changed, 2023 insertions(+) create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/benchmark.native.js create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/Makefile create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/benchmark.c create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/binding.gyp create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/Makefile create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/include.gypi create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/include/stdlib/math/base/special/rempio2.h create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/lib/native.js create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/src/Makefile create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/src/addon.c create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/README.md b/lib/node_modules/@stdlib/math/base/special/rempio2/README.md index 32551ba1d644..8c32bb0e53c5 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/README.md +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/README.md @@ -111,6 +111,97 @@ for ( i = 0; i < x.length; i++ ) { + + +* * * + +
+ +## C APIs + + + +
+ +
+ + + + + +
+ +### Usage + +```c +#include "stdlib/math/base/special/rempio2.h" +``` + +#### stdlib_base_rempio2( x, &rem1, &rem2 ) + +Computes `x - nπ/2 = r`. The function returns `n` and stores the remainder `r` as two numbers in `y`, such that `y[0]+y[1] = r`. + +```c +double rem1; +double rem2; + +stdlib_base_rempio2( 4.0, &rem1, &rem2 ); +``` + +The function accepts the following arguments: + +- **x**: `[in] double` input value. +- **rem1**: `[out] double*` destination for rem1. +- **rem2**: `[out] double*` destination for rem2. + +```c +void stdlib_base_rempio2( const double x, double *rem1, double *rem2 ); +``` + +
+ + + + + +
+ +
+ + + + + +
+ +### Examples + +```c +#include "stdlib/math/base/special/rempio2.h" +#include + +int main( void ) { + const double x[] = { 0.0, 1.0, 4.0, 128.0 }; + + double rem1; + double rem2; + double n; + int i; + for ( i = 0; i < 4; i++ ) { + stdlib_base_rempio2( x[ i ], &rem1, &rem2 ); + printf( "%lf - %lfπ/2 = %lf + %lf", x[ i ], n, rem1, rem2 ); + } +} +``` + +
+ + + +
+ + +
diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/benchmark.native.js b/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/benchmark.native.js new file mode 100644 index 000000000000..ba0fa56577dd --- /dev/null +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/benchmark.native.js @@ -0,0 +1,63 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2024 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 resolve = require( 'path' ).resolve; +var bench = require( '@stdlib/bench' ); +var randu = require( '@stdlib/random/base/randu' ); +var isnan = require( '@stdlib/math/base/assert/is-nan' ); +var tryRequire = require( '@stdlib/utils/try-require' ); +var pkg = require( './../package.json' ).name; + + +// VARIABLES // + +var rempio2 = tryRequire( resolve( __dirname, './../lib/native.js' ) ); +var opts = { + 'skip': ( rempio2 instanceof Error ) +}; + + +// MAIN // + +bench( pkg+'::native', opts, function benchmark( b ) { + var x; + var y; + var n; + var i; + + y = [ 0.0, 0.0 ]; + + b.tic(); + for ( i = 0; i < b.iterations; i++ ) { + x = ( randu() * 200.0 ) - 100.0; + n = rempio2( x, y ); + if ( isnan( n ) ) { + b.fail( 'should not return NaN' ); + } + } + b.toc(); + if ( isnan( n ) || isnan( y[0] ) || isnan( y[1] ) ) { + b.fail( 'should not return NaN' ); + } + b.pass( 'benchmark finished' ); + b.end(); +}); diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/Makefile b/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/Makefile new file mode 100644 index 000000000000..f69e9da2b4d3 --- /dev/null +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/Makefile @@ -0,0 +1,146 @@ +#/ +# @license Apache-2.0 +# +# Copyright (c) 2024 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. +#/ + +# VARIABLES # + +ifndef VERBOSE + QUIET := @ +else + QUIET := +endif + +# Determine the OS ([1][1], [2][2]). +# +# [1]: https://en.wikipedia.org/wiki/Uname#Examples +# [2]: http://stackoverflow.com/a/27776822/2225624 +OS ?= $(shell uname) +ifneq (, $(findstring MINGW,$(OS))) + OS := WINNT +else +ifneq (, $(findstring MSYS,$(OS))) + OS := WINNT +else +ifneq (, $(findstring CYGWIN,$(OS))) + OS := WINNT +else +ifneq (, $(findstring Windows_NT,$(OS))) + OS := WINNT +endif +endif +endif +endif + +# Define the program used for compiling C source files: +ifdef C_COMPILER + CC := $(C_COMPILER) +else + CC := gcc +endif + +# Define the command-line options when compiling C files: +CFLAGS ?= \ + -std=c99 \ + -O3 \ + -Wall \ + -pedantic + +# Determine whether to generate position independent code ([1][1], [2][2]). +# +# [1]: https://gcc.gnu.org/onlinedocs/gcc/Code-Gen-Options.html#Code-Gen-Options +# [2]: http://stackoverflow.com/questions/5311515/gcc-fpic-option +ifeq ($(OS), WINNT) + fPIC ?= +else + fPIC ?= -fPIC +endif + +# List of includes (e.g., `-I /foo/bar -I /beep/boop/include`): +INCLUDE ?= + +# List of source files: +SOURCE_FILES ?= + +# List of libraries (e.g., `-lopenblas -lpthread`): +LIBRARIES ?= + +# List of library paths (e.g., `-L /foo/bar -L /beep/boop`): +LIBPATH ?= + +# List of C targets: +c_targets := benchmark.out + + +# RULES # + +#/ +# Compiles source files. +# +# @param {string} [C_COMPILER] - C compiler (e.g., `gcc`) +# @param {string} [CFLAGS] - C compiler options +# @param {(string|void)} [fPIC] - compiler flag determining whether to generate position independent code (e.g., `-fPIC`) +# @param {string} [INCLUDE] - list of includes (e.g., `-I /foo/bar -I /beep/boop/include`) +# @param {string} [SOURCE_FILES] - list of source files +# @param {string} [LIBPATH] - list of library paths (e.g., `-L /foo/bar -L /beep/boop`) +# @param {string} [LIBRARIES] - list of libraries (e.g., `-lopenblas -lpthread`) +# +# @example +# make +# +# @example +# make all +#/ +all: $(c_targets) + +.PHONY: all + +#/ +# Compiles C source files. +# +# @private +# @param {string} CC - C compiler (e.g., `gcc`) +# @param {string} CFLAGS - C compiler options +# @param {(string|void)} fPIC - compiler flag determining whether to generate position independent code (e.g., `-fPIC`) +# @param {string} INCLUDE - list of includes (e.g., `-I /foo/bar`) +# @param {string} SOURCE_FILES - list of source files +# @param {string} LIBPATH - list of library paths (e.g., `-L /foo/bar`) +# @param {string} LIBRARIES - list of libraries (e.g., `-lopenblas`) +#/ +$(c_targets): %.out: %.c + $(QUIET) $(CC) $(CFLAGS) $(fPIC) $(INCLUDE) -o $@ $(SOURCE_FILES) $< $(LIBPATH) -lm $(LIBRARIES) + +#/ +# Runs compiled benchmarks. +# +# @example +# make run +#/ +run: $(c_targets) + $(QUIET) ./$< + +.PHONY: run + +#/ +# Removes generated files. +# +# @example +# make clean +#/ +clean: + $(QUIET) -rm -f *.o *.out + +.PHONY: clean diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/benchmark.c b/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/benchmark.c new file mode 100644 index 000000000000..78439da428ff --- /dev/null +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/benchmark.c @@ -0,0 +1,138 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2024 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. +*/ + +/** +* Benchmark `rempio2`. +*/ +#include "stdlib/math/base/special/rempio2.h" +#include +#include +#include +#include +#include + +#define NAME "rempio2" +#define ITERATIONS 1000000 +#define REPEATS 3 + +/** +* Prints the TAP version. +*/ +void print_version() { + printf( "TAP version 13\n" ); +} + +/** +* Prints the TAP summary. +* +* @param total total number of tests +* @param passing total number of passing tests +*/ +void print_summary( int total, int passing ) { + printf( "#\n" ); + printf( "1..%d\n", total ); // TAP plan + printf( "# total %d\n", total ); + printf( "# pass %d\n", passing ); + printf( "#\n" ); + printf( "# ok\n" ); +} + +/** +* Prints benchmarks results. +* +* @param elapsed elapsed time in seconds +*/ +void print_results( double elapsed ) { + double rate = (double)ITERATIONS / elapsed; + printf( " ---\n" ); + printf( " iterations: %d\n", ITERATIONS ); + printf( " elapsed: %0.9f\n", elapsed ); + printf( " rate: %0.9f\n", rate ); + printf( " ...\n" ); +} + +/** +* Returns a clock time. +* +* @return clock time +*/ +double tic() { + struct timeval now; + gettimeofday( &now, NULL ); + return (double)now.tv_sec + (double)now.tv_usec/1.0e6; +} + +/** +* Generates a random number on the interval [0,1]. +* +* @return random number +*/ +double rand_double() { + int r = rand(); + return (double)r / ( (double)RAND_MAX + 1.0 ); +} + +/** +* Runs a benchmark. +* +* @return elapsed time in seconds +*/ +double benchmark() { + double elapsed; + double rem1; + double rem2; + double x; + double y; + double t; + int i; + + t = tic(); + for ( i = 0; i < ITERATIONS; i++ ) { + x = ( rand_double() * 200.0 ) - 100.0; + y = stdlib_base_rempio2( x, &rem1, &rem2 ); + if ( y != y ) { + printf( "should not return NaN\n" ); + break; + } + } + elapsed = tic() - t; + if ( y != y ) { + printf( "should not return NaN\n" ); + } + return elapsed; +} + +/** +* Main execution sequence. +*/ +int main( void ) { + double elapsed; + int i; + + // Use the current time to seed the random number generator: + srand( time( NULL ) ); + + print_version(); + for ( i = 0; i < REPEATS; i++ ) { + printf( "# c::native::%s\n", NAME ); + elapsed = benchmark(); + print_results( elapsed ); + printf( "ok %d benchmark finished\n", i+1 ); + } + print_summary( REPEATS, REPEATS ); +} diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/binding.gyp b/lib/node_modules/@stdlib/math/base/special/rempio2/binding.gyp new file mode 100644 index 000000000000..ec3992233442 --- /dev/null +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/binding.gyp @@ -0,0 +1,170 @@ +# @license Apache-2.0 +# +# Copyright (c) 2024 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. + +# A `.gyp` file for building a Node.js native add-on. +# +# [1]: https://gyp.gsrc.io/docs/InputFormatReference.md +# [2]: https://gyp.gsrc.io/docs/UserDocumentation.md +{ + # List of files to include in this file: + 'includes': [ + './include.gypi', + ], + + # Define variables to be used throughout the configuration for all targets: + 'variables': { + # Target name should match the add-on export name: + 'addon_target_name%': 'addon', + + # Set variables based on the host OS: + 'conditions': [ + [ + 'OS=="win"', + { + # Define the object file suffix: + 'obj': 'obj', + }, + { + # Define the object file suffix: + 'obj': 'o', + } + ], # end condition (OS=="win") + ], # end conditions + }, # end variables + + # Define compile targets: + 'targets': [ + + # Target to generate an add-on: + { + # The target name should match the add-on export name: + 'target_name': '<(addon_target_name)', + + # Define dependencies: + 'dependencies': [], + + # Define directories which contain relevant include headers: + 'include_dirs': [ + # Local include directory: + '<@(include_dirs)', + ], + + # List of source files: + 'sources': [ + '<@(src_files)', + ], + + # Settings which should be applied when a target's object files are used as linker input: + 'link_settings': { + # Define libraries: + 'libraries': [ + '<@(libraries)', + ], + + # Define library directories: + 'library_dirs': [ + '<@(library_dirs)', + ], + }, + + # C/C++ compiler flags: + 'cflags': [ + # Enable commonly used warning options: + '-Wall', + + # Aggressive optimization: + '-O3', + ], + + # C specific compiler flags: + 'cflags_c': [ + # Specify the C standard to which a program is expected to conform: + '-std=c99', + ], + + # C++ specific compiler flags: + 'cflags_cpp': [ + # Specify the C++ standard to which a program is expected to conform: + '-std=c++11', + ], + + # Linker flags: + 'ldflags': [], + + # Apply conditions based on the host OS: + 'conditions': [ + [ + 'OS=="mac"', + { + # Linker flags: + 'ldflags': [ + '-undefined dynamic_lookup', + '-Wl,-no-pie', + '-Wl,-search_paths_first', + ], + }, + ], # end condition (OS=="mac") + [ + 'OS!="win"', + { + # C/C++ flags: + 'cflags': [ + # Generate platform-independent code: + '-fPIC', + ], + }, + ], # end condition (OS!="win") + ], # end conditions + }, # end target <(addon_target_name) + + # Target to copy a generated add-on to a standard location: + { + 'target_name': 'copy_addon', + + # Declare that the output of this target is not linked: + 'type': 'none', + + # Define dependencies: + 'dependencies': [ + # Require that the add-on be generated before building this target: + '<(addon_target_name)', + ], + + # Define a list of actions: + 'actions': [ + { + 'action_name': 'copy_addon', + 'message': 'Copying addon...', + + # Explicitly list the inputs in the command-line invocation below: + 'inputs': [], + + # Declare the expected outputs: + 'outputs': [ + '<(addon_output_dir)/<(addon_target_name).node', + ], + + # Define the command-line invocation: + 'action': [ + 'cp', + '<(PRODUCT_DIR)/<(addon_target_name).node', + '<(addon_output_dir)/<(addon_target_name).node', + ], + }, + ], # end actions + }, # end target copy_addon + ], # end targets +} diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/Makefile b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/Makefile new file mode 100644 index 000000000000..6aed70daf167 --- /dev/null +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/Makefile @@ -0,0 +1,146 @@ +#/ +# @license Apache-2.0 +# +# Copyright (c) 2024 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. +#/ + +# VARIABLES # + +ifndef VERBOSE + QUIET := @ +else + QUIET := +endif + +# Determine the OS ([1][1], [2][2]). +# +# [1]: https://en.wikipedia.org/wiki/Uname#Examples +# [2]: http://stackoverflow.com/a/27776822/2225624 +OS ?= $(shell uname) +ifneq (, $(findstring MINGW,$(OS))) + OS := WINNT +else +ifneq (, $(findstring MSYS,$(OS))) + OS := WINNT +else +ifneq (, $(findstring CYGWIN,$(OS))) + OS := WINNT +else +ifneq (, $(findstring Windows_NT,$(OS))) + OS := WINNT +endif +endif +endif +endif + +# Define the program used for compiling C source files: +ifdef C_COMPILER + CC := $(C_COMPILER) +else + CC := gcc +endif + +# Define the command-line options when compiling C files: +CFLAGS ?= \ + -std=c99 \ + -O3 \ + -Wall \ + -pedantic + +# Determine whether to generate position independent code ([1][1], [2][2]). +# +# [1]: https://gcc.gnu.org/onlinedocs/gcc/Code-Gen-Options.html#Code-Gen-Options +# [2]: http://stackoverflow.com/questions/5311515/gcc-fpic-option +ifeq ($(OS), WINNT) + fPIC ?= +else + fPIC ?= -fPIC +endif + +# List of includes (e.g., `-I /foo/bar -I /beep/boop/include`): +INCLUDE ?= + +# List of source files: +SOURCE_FILES ?= + +# List of libraries (e.g., `-lopenblas -lpthread`): +LIBRARIES ?= + +# List of library paths (e.g., `-L /foo/bar -L /beep/boop`): +LIBPATH ?= + +# List of C targets: +c_targets := example.out + + +# RULES # + +#/ +# Compiles source files. +# +# @param {string} [C_COMPILER] - C compiler (e.g., `gcc`) +# @param {string} [CFLAGS] - C compiler options +# @param {(string|void)} [fPIC] - compiler flag determining whether to generate position independent code (e.g., `-fPIC`) +# @param {string} [INCLUDE] - list of includes (e.g., `-I /foo/bar -I /beep/boop/include`) +# @param {string} [SOURCE_FILES] - list of source files +# @param {string} [LIBPATH] - list of library paths (e.g., `-L /foo/bar -L /beep/boop`) +# @param {string} [LIBRARIES] - list of libraries (e.g., `-lopenblas -lpthread`) +# +# @example +# make +# +# @example +# make all +#/ +all: $(c_targets) + +.PHONY: all + +#/ +# Compiles C source files. +# +# @private +# @param {string} CC - C compiler (e.g., `gcc`) +# @param {string} CFLAGS - C compiler options +# @param {(string|void)} fPIC - compiler flag determining whether to generate position independent code (e.g., `-fPIC`) +# @param {string} INCLUDE - list of includes (e.g., `-I /foo/bar`) +# @param {string} SOURCE_FILES - list of source files +# @param {string} LIBPATH - list of library paths (e.g., `-L /foo/bar`) +# @param {string} LIBRARIES - list of libraries (e.g., `-lopenblas`) +#/ +$(c_targets): %.out: %.c + $(QUIET) $(CC) $(CFLAGS) $(fPIC) $(INCLUDE) -o $@ $(SOURCE_FILES) $< $(LIBPATH) -lm $(LIBRARIES) + +#/ +# Runs compiled examples. +# +# @example +# make run +#/ +run: $(c_targets) + $(QUIET) ./$< + +.PHONY: run + +#/ +# Removes generated files. +# +# @example +# make clean +#/ +clean: + $(QUIET) -rm -f *.o *.out + +.PHONY: clean diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c new file mode 100644 index 000000000000..8286358afd0a --- /dev/null +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c @@ -0,0 +1,33 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2024 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. +*/ + +#include "stdlib/math/base/special/rempio2.h" +#include + +int main( void ) { + const double x[] = { 0.0, 1.0, 4.0, 128.0 }; + + double rem1; + double rem2; + double n; + int i; + for ( i = 0; i < 4; i++ ) { + n = stdlib_base_rempio2( x[ i ], &rem1, &rem2 ); + printf( "%lf - %lfπ/2 = %lf + %lf", x[ i ], n, rem1, rem2 ); + } +} diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/include.gypi b/lib/node_modules/@stdlib/math/base/special/rempio2/include.gypi new file mode 100644 index 000000000000..575cb043c0bf --- /dev/null +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/include.gypi @@ -0,0 +1,53 @@ +# @license Apache-2.0 +# +# Copyright (c) 2024 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. + +# A GYP include file for building a Node.js native add-on. +# +# Main documentation: +# +# [1]: https://gyp.gsrc.io/docs/InputFormatReference.md +# [2]: https://gyp.gsrc.io/docs/UserDocumentation.md +{ + # Define variables to be used throughout the configuration for all targets: + 'variables': { + # Source directory: + 'src_dir': './src', + + # Include directories: + 'include_dirs': [ + ' +#include +#include + +/** +* Receives JavaScript callback invocation data. +* +* @private +* @param env environment under which the function is invoked +* @param info callback data +* @return Node-API value +*/ +static napi_value addon( napi_env env, napi_callback_info info ) { + napi_status status; + + // Get callback arguments: + size_t argc = 2; + napi_value argv[ 2 ]; + status = napi_get_cb_info( env, info, &argc, argv, NULL, NULL ); + assert( status == napi_ok ); + + // Check whether we were provided the correct number of arguments: + if ( argc < 2 ) { + status = napi_throw_error( env, NULL, "invalid invocation. Insufficient arguments." ); + assert( status == napi_ok ); + return NULL; + } + if ( argc > 2 ) { + status = napi_throw_error( env, NULL, "invalid invocation. Too many arguments." ); + assert( status == napi_ok ); + return NULL; + } + + bool res; + status = napi_is_typedarray( env, argv[ 0 ], &res ); + assert( status == napi_ok ); + if ( res == false ) { + status = napi_throw_type_error( env, NULL, "invalid argument. First argument must be a Float64Array." ); + assert( status == napi_ok ); + return NULL; + } + + napi_valuetype vtype1; + status = napi_typeof( env, argv[ 1 ], &vtype1 ); + assert( status == napi_ok ); + if ( vtype1 != napi_number ) { + status = napi_throw_type_error( env, NULL, "invalid argument. Second argument must be a number." ); + assert( status == napi_ok ); + return NULL; + } + + napi_typedarray_type vtype0; + size_t len; + void *Out; + status = napi_get_typedarray_info( env, argv[ 0 ], &vtype0, &len, &Out, NULL, NULL ); + assert( status == napi_ok ); + if ( vtype0 != napi_float64_array ) { + status = napi_throw_type_error( env, NULL, "invalid argument. First argument must be a Float64Array." ); + assert( status == napi_ok ); + return NULL; + } + if ( len != 2 ) { + status = napi_throw_range_error( env, NULL, "invalid argument. First argument must have 2 elements." ); + assert( status == napi_ok ); + return NULL; + } + + double value; + status = napi_get_value_double( env, argv[ 1 ], &value ); + assert( status == napi_ok ); + + double rem1; + double rem2; + napi_value v; + status = napi_create_double( env, stdlib_base_rempio2( value, &rem1, &rem2 ), &v ); + assert( status == napi_ok ); + + double *op = (double *)Out; + op[ 0 ] = rem1; + op[ 1 ] = rem2; + + return v; +} + +/** +* Initializes a Node-API module. +* +* @private +* @param env environment under which the function is invoked +* @param exports exports object +* @return main export +*/ +static napi_value init( napi_env env, napi_value exports ) { + napi_value fcn; + napi_status status = napi_create_function( env, "exports", NAPI_AUTO_LENGTH, addon, NULL, &fcn ); + assert( status == napi_ok ); + return fcn; +} + +NAPI_MODULE( NODE_GYP_MODULE_NAME, init ) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c new file mode 100644 index 000000000000..9ca6f82ad658 --- /dev/null +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c @@ -0,0 +1,530 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2024 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, 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_rem_pio2.c}. The implementation follows the original, but has been modified according to project conventions. +* +* ```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. +* ``` +*/ + +#include "stdlib/math/base/special/rempio2.h" +#include "stdlib/math/base/special/round.h" +#include "stdlib/number/float64/base/get_high_word.h" +#include "stdlib/number/float64/base/get_low_word.h" +#include "stdlib/number/float64/base/from_words.h" +#include "stdlib/math/base/special/floor.h" +#include "stdlib/math/base/special/ldexp.h" +#include "stdlib/constants/float64/high_word_abs_mask.h" +#include "stdlib/constants/float64/high_word_exponent_mask.h" +#include "stdlib/constants/float64/high_word_sign_mask.h" +#include + +static const double ZERO = 0.00000000000000000000e+00; // 0x00000000, 0x00000000 +static const double TWO24 = 1.67772160000000000000e+07; // 0x41700000, 0x00000000 +static const double TWON24 = 5.96046447753906250000e-08; // 0x3E700000, 0x00000000 +static const double INVPIO2 = 6.36619772367581382433e-01; // 0x3FE45F30, 0x6DC9C883 +static const double PIO2_1 = 1.57079632673412561417e+00; // 0x3FF921FB, 0x54400000 +static const double PIO2_1T = 6.07710050650619224932e-11; // 0x3DD0B461, 0x1A626331 +static const double PIO2_2 = 6.07710050630396597660e-11; // 0x3DD0B461, 0x1A600000 +static const double PIO2_2T = 2.02226624879595063154e-21; // 0x3BA3198A, 0x2E037073 +static const double PIO2_3 = 2.02226624871116645580e-21; // 0x3BA3198A, 0x2E000000 +static const double PIO2_3T = 8.47842766036889956997e-32; // 0x397B839A, 0x252049C1 + +// High word significand for π and π/2: 0x921fb = 598523 => 00000000000010010010000111111011 +static const int32_t PI_HIGH_WORD_SIGNIFICAND = 0x921fb; + +// High word for π/4: 0x3fe921fb = 1072243195 => 00111111111010010010000111111011 +static const int32_t PIO4_HIGH_WORD = 0x3fe921fb; + +// High word for 3π/4: 0x4002d97c = 1073928572 => 01000000000000101101100101111100 +static const int32_t THREE_PIO4_HIGH_WORD = 0x4002d97c; + +// High word for 5π/4: 0x400f6a7a = 1074752122 => 01000000000011110110101001111010 +static const int32_t FIVE_PIO4_HIGH_WORD = 0x400f6a7a; + +// High word for 6π/4: 0x4012d97c = 1074977148 => 01000000000100101101100101111100 +static const int32_t THREE_PIO2_HIGH_WORD = 0x4012d97c; + +// High word for 7π/4: 0x4015fdbc = 1075183036 => 01000000000101011111110110111100 +static const int32_t SEVEN_PIO4_HIGH_WORD = 0x4015fdbc; + +// High word for 8π/4: 0x401921fb = 1075388923 => 01000000000110010010000111111011 +static const int32_t TWO_PI_HIGH_WORD = 0x401921fb; + +// High word for 9π/4: 0x401c463b = 1075594811 => 01000000000111000100011000111011 +static const int32_t NINE_PIO4_HIGH_WORD = 0x401c463b; + +// 2^20*π/2 = 1647099.3291652855 => 0100000100111001001000011111101101010100010001000010110100011000 => high word => 0x413921fb = 1094263291 => 01000001001110010010000111111011 +static const int32_t MEDIUM = 0x413921fb; + + +/* +* Table of constants for `2/π` (`396` hex digits, `476` decimal). +* +* Integer array which contains the (`24*i`)-th to (`24*i+23`)-th bit of `2/π` after binary point. The corresponding floating value is +* +* ```tex +* \operatorname{ipio2}[i] \cdot 2^{-24(i+1)} +* ``` +* +* This table must have at least `(e0-3)/24 + jk` terms. For quad precision (`e0 <= 16360`, `jk = 6`), this is `686`. +*/ +static const int32_t IPIO2[] = { + 0xA2F983, 0x6E4E44, 0x1529FC, 0x2757D1, 0xF534DD, 0xC0DB62, + 0x95993C, 0x439041, 0xFE5163, 0xABDEBB, 0xC561B7, 0x246E3A, + 0x424DD2, 0xE00649, 0x2EEA09, 0xD1921C, 0xFE1DEB, 0x1CB129, + 0xA73EE8, 0x8235F5, 0x2EBB44, 0x84E99C, 0x7026B4, 0x5F7E41, + 0x3991D6, 0x398353, 0x39F49C, 0x845F8B, 0xBDF928, 0x3B1FF8, + 0x97FFDE, 0x05980F, 0xEF2F11, 0x8B5A0A, 0x6D1F6D, 0x367ECF, + 0x27CB09, 0xB74F46, 0x3F669E, 0x5FEA2D, 0x7527BA, 0xC7EBE5, + 0xF17B3D, 0x0739F7, 0x8A5292, 0xEA6BFB, 0x5FB11F, 0x8D5D08, + 0x560330, 0x46FC7B, 0x6BABF0, 0xCFBC20, 0x9AF436, 0x1DA9E3, + 0x91615E, 0xE61B08, 0x659985, 0x5F14A0, 0x68408D, 0xFFD880, + 0x4D7327, 0x310606, 0x1556CA, 0x73A8C9, 0x60E27B, 0xC08C6B +}; + +// Double precision array, obtained by cutting `π/2` into `24` bits chunks... +static const double PIO2[] = { + 1.57079625129699707031e+00, // 0x3FF921FB, 0x40000000 + 7.54978941586159635335e-08, // 0x3E74442D, 0x00000000 + 5.39030252995776476554e-15, // 0x3CF84698, 0x80000000 + 3.28200341580791294123e-22, // 0x3B78CC51, 0x60000000 + 1.27065575308067607349e-29, // 0x39F01B83, 0x80000000 + 1.22933308981111328932e-36, // 0x387A2520, 0x40000000 + 2.73370053816464559624e-44, // 0x36E38222, 0x80000000 + 2.16741683877804819444e-51 // 0x3569F31D, 0x00000000 +}; + +// Arrays for storing temporary values: +double TX[ 3 ] = { 0.0 }; +double TY[ 2 ] = { 0.0 }; +double F[ 20 ]; +double Q[ 20 ]; +double FQ[ 20 ]; +int32_t IQ[ 20 ]; + +/** +* Returns the last three binary digits of `N` with `y = x - Nπ/2` so that `|y| < π/2`. +* +* ## Method +* +* - The method is to compute the integer (`mod 8`) and fraction parts of `2x/π` without doing the full multiplication. In general, we skip the part of the product that is known to be a huge integer (more accurately, equals `0 mod 8` ). Thus, the number of operations is independent of the exponent of the input. +* +* @private +* @param x input value +* @param rem1 remainder element 1 +* @param rem2 remainder element 2 +* @param e0 the exponent of `x[0]` (must be <= 16360) +* @param nx dimension of `x[]` +*/ +double rempio2Kernel( const double* x, const double* y, int32_t e0, int32_t nx ) { + int32_t carry; + int32_t jz; + int32_t jx; + int32_t jv; + int32_t jp; + int32_t jk; + int32_t q0; + int32_t ih; + int32_t n; + int32_t i; + int32_t j; + int32_t k; + int32_t m; + double fw; + double z; + + // Initialize `jk` for double-precision floating-point numbers: + jk = 4; + jp = jk; + + // Determine `jx`, `jv`, `q0` (note that `q0 < 3`): + jx = nx - 1; + jv = ( e0 - 3 ) / 24; + if ( jv < 0 ) { + jv = 0; + } + q0 = e0 - ( 24 * ( jv + 1 ) ); + + // Set up `F[0]` to `F[jx+jk]` where `F[jx+jk] = IPIO2[jv+jk]`: + j = jv - jx; + m = jx + jk; + for ( i = 0; i <= m; i++ ) { + if ( j < 0 ) { + F[ i ] = ZERO; + } else { + F[ i ] = ( double )IPIO2[ j ]; + } + j++; + } + + // Compute `Q[0],Q[1],...,Q[jk]`: + for ( i = 0; i <= jk; i++ ) { + fw = 0.0; + for ( j = 0; j <= jx; j++ ) { + fw += x[ j ] * F[ jx + ( i - j ) ]; + } + Q[ i ] = fw; + } + + jz = jk; +recompute: + // Distill `q[]` into `IQ[]` in reverse order... + z = Q[ jz ]; + j = jz; + for ( i = 0; j > 0; i++ ) { + fw = ( double )( ( int32_t )( TWON24 * z ) ); + IQ[ i ] = (int32_t)( z - ( TWO24 * fw ) ); + z = Q[ j - 1 ] + fw; + j--; + } + // Compute `n`... + z = stdlib_base_ldexp( z, q0 ); + z -= 8.0 * stdlib_base_floor( z * 0.125 ); // Trim off integer >= 8 + n = ( int32_t )z; + z -= ( double )n; + ih = 0; + if ( q0 > 0 ) { + // Need `IQ[jz-1]` to determine `n`... + i = ( IQ[ jz - 1 ] >> ( 24 - q0 ) ); + n += i; + IQ[ jz - 1 ] -= ( i << ( 24 - q0 ) ); + ih = ( IQ[ jz - 1 ] >> ( 23 - q0 ) ); + } else if ( q0 == 0 ) { + ih = ( IQ[ jz - 1 ] >> 23 ); + } else if ( z >= 0.5 ) { + ih = 2; + } + + // Case: q > 0.5 + if ( ih > 0 ) { + n += 1; + carry = 0; + // Compute `1-q`: + for ( i = 0; i < jz; i++ ) { + j = IQ[ i ]; + if ( carry == 0 ) { + if ( j != 0 ) { + carry = 1; + IQ[ i ] = 0x1000000 - j; + } + } else { + IQ[ i ] = 0xffffff - j; + } + } + if ( q0 > 0 ) { + // Rare case: chance is 1 in 12... + switch ( q0 ) { + case 1: + IQ[ jz - 1 ] &= 0x7fffff; + break; + case 2: + IQ[ jz - 1 ] &= 0x3fffff; + break; + } + } + if ( ih == 2 ) { + z = 1.0 - z; + if ( carry != 0 ) { + z -= stdlib_base_ldexp( 1.0, q0 ); + } + } + } + + // Check if re-computation is needed... + if ( z == ZERO ) { + j = 0; + for ( i = jz - 1; i >= jk; i-- ) { + j |= IQ[ i ]; + } + if ( j == 0 ) { + // Need re-computation... + for ( k = 1; IQ[ jk - k ] == 0; k++ ) { + // `k` is the number of terms needed... + } + + for ( i = jz + 1; i <= jz + k; i++ ) { + // Add `q[jz+1]` to `q[jz+k]`... + F[ jx + i ] = ( double )IPIO2[ jv + i ]; + fw = 0.0; + for ( j = 0; j <= jx; j++ ) { + fw += x[ j ] * F[ jx + ( i - j ) ]; + } + Q[ i ] = fw; + } + jz += k; + goto recompute; + } + } + + // Chop off zero terms... + if ( z == 1.0 ) { + jz -= 1; + q0 -= 24; + while ( IQ[ jz ] == 0 ) { + jz--; + q0 -= 24; + } + } else { + // Break `z` into 24-bit if necessary... + z = stdlib_base_ldexp( z, -q0 ); + if ( z >= TWO24 ) { + fw = ( double )( ( int32_t )( TWON24 * z ) ); + IQ[ jz ] = ( int32_t )( z - ( TWO24 * fw ) ); + jz += 1; + q0 += 24; + IQ[ jz ] = ( int32_t )fw; + } else { + IQ[ jz ] = ( int32_t )z; + } + } + + // Convert integer "bit" chunk to floating-point value... + fw = stdlib_base_ldexp( 1.0, q0 ); + for ( i = jz; i >= 0; i-- ) { + Q[ i ] = fw * ( double )IQ[ i ]; + fw *= TWON24; + } + + // Compute `PIO2[0,...,jp]*q[jz,...,0]`... + for ( i = jz; i >= 0; i-- ) { + fw = 0.0; + for ( k = 0; k <= jp && k <= jz - i; k++ ) { + fw += PIO2[ k ] * Q[ i + k ]; + } + FQ[ jz - i ] = fw; + } + // Compress `FQ[]` into `y[]`... + fw = 0.0; + for ( i = jz; i >= 0; i-- ) { + fw += FQ[ i ]; + } + if ( ih == 0 ) { + y[ 0 ] = fw; + } else { + y[ 0 ] = -fw; + } + fw = FQ[ 0 ] - fw; + for ( i = 1; i <= jz; i++ ) { + fw += FQ[i]; + } + if ( ih == 0 ) { + y[ 1 ] = fw; + } else { + y[ 1 ] = -fw; + } + return n & 7; +} + + +/** +* Computes `x - nπ/2 = r`. +* +* ## Notes +* +* - Returns `n` and stores the remainder `r` as two numbers `rem1` and `rem2`, such that `rem1+rem2 = r`. +* +* @param x input value +* @param rem1 remainder element 1 +* @param rem2 remainder element 2 +* +* @example +* double x = 128.0; +* +* double rem1; +* double rem2; +* stdlib_base_rempio2( x, &rem1, &rem2 ); +*/ +double stdlib_base_rempio2( const double x, double* rem1, double* rem2 ) { + uint32_t uhigh; + uint32_t ulow; + uint32_t uhx; + int32_t high; + int32_t low; + int32_t ix; + int32_t hx; + int32_t e0; + int32_t nx; + int32_t i; + int32_t j; + int32_t n; + double z; + double w; + double t; + double r; + + // high word of x + stdlib_base_float64_get_high_word( x, &uhx ); + hx = ( int32_t )uhx; + ix = hx & STDLIB_CONSTANT_FLOAT64_HIGH_WORD_ABS_MASK; + + // Case: |x| ~<= π/4 (no need for reduction) + if ( ix <= PIO4_HIGH_WORD ) { + *rem1 = x; + *rem2 = 0; + return 0; + } + // Case: |x| ~<= 5π/4 + if ( ix <= FIVE_PIO4_HIGH_WORD ) { + // Case: |x| ~= π/2 or π + if ( ( ix & STDLIB_CONSTANT_FLOAT64_HIGH_WORD_SIGN_MASK ) == PI_HIGH_WORD_SIGNIFICAND ) { + // Cancellation => use medium case + goto medium; + } + // Case: |x| ~<= 3π/4 + if ( ix <= THREE_PIO4_HIGH_WORD ) { + if ( hx > 0 ) { + // one round good to 85 bits + z = x - PIO2_1; + *rem1 = z - PIO2_1T; + *rem2 = ( z - *rem1 ) - PIO2_1T; + return 1; + } else { + z = x + PIO2_1; + *rem1 = z + PIO2_1T; + *rem2 = ( z - *rem1 ) + PIO2_1T; + return -1; + } + } else { + if ( hx > 0 ) { + z = x - ( 2.0 * PIO2_1 ); + *rem1 = z - ( 2 * PIO2_1T ); + *rem2 = ( z - *rem1 ) - ( 2 * PIO2_1T ); + return 2; + } else { + z = x + ( 2.0 * PIO2_1 ); + *rem1 = z + ( 2.0 * PIO2_1T ); + *rem2 = ( z - *rem1 ) + ( 2.0 * PIO2_1T ); + return -2; + } + } + } + // Case: |x| ~<= 9π/4 + if ( ix <= NINE_PIO4_HIGH_WORD ) { + // Case: |x| ~<= 7π/4 + if ( ix <= SEVEN_PIO4_HIGH_WORD ) { + // Case: |x| ~= 3π/2 + if ( ix == THREE_PIO2_HIGH_WORD ) { + goto medium; + } + if ( hx > 0 ) { + z = x - ( 3 * PIO2_1 ); + *rem1 = z - ( 3 * PIO2_1T ); + *rem2 = ( z - *rem1 ) - ( 3 * PIO2_1T ); + return 3; + } else { + z = x + ( 3 * PIO2_1 ); + *rem1 = z + ( 3 * PIO2_1T ); + *rem2 = ( z - *rem1 ) + ( 3.0 * PIO2_1T ); + return -3; + } + } else { + // Case: |x| ~= 4π/2 + if ( ix == TWO_PI_HIGH_WORD ) { + goto medium; + } + if ( hx > 0 ) { + z = x - ( 4 * PIO2_1 ); + *rem1 = z - ( 4 * PIO2_1T ); + *rem2 = ( z - *rem1 ) - ( 4 * PIO2_1T ); + return 4; + } else { + z = x + ( 4 * PIO2_1 ); + *rem1 = z + ( 4 * PIO2_1T ); + *rem2 = ( z - *rem1 ) + ( 4 * PIO2_1T ); + return -4; + } + } + } + // Case: |x| ~< 2^20*π/2 (medium size) + if ( ix < MEDIUM ) { +medium: + n = stdlib_base_round( x * INVPIO2 ); + r = x - ( n * PIO2_1 ); + w = n * PIO2_1T; + + // First rounding (good to 85 bits)... + j = ix >> 20; + *rem1 = r - w; + stdlib_base_float64_get_high_word( *rem1, &uhigh ); + high = ( int32_t )uhigh; + i = j - ( ( high >> 20 ) & STDLIB_CONSTANT_FLOAT64_HIGH_WORD_EXPONENT_MASK ); + + // Check if a second iteration is needed (good to 118 bits)... + if ( i > 16 ) { + t = r; + w = n * PIO2_2; + r = t - w; + w = ( n * PIO2_2T ) - ( ( t - r ) - w ); + *rem1 = r - w; + stdlib_base_float64_get_high_word( *rem1, &uhigh ); + high = ( int32_t )uhigh; + i = j - ( ( high >> 20 ) & STDLIB_CONSTANT_FLOAT64_HIGH_WORD_EXPONENT_MASK ); + + // Check if a third iteration is needed (151 bits accumulated)... + if ( i > 49 ) { + // will cover all possible cases + t = r; + w = n * PIO2_3; + r = t - w; + w = ( n * PIO2_3T ) - ( ( t - r ) - w ); + *rem1 = r - w; + } + } + *rem2 = ( r - *rem1 ) - w; + return n; + } + // Case: x is NaN or infinity + if( ix >= STDLIB_CONSTANT_FLOAT64_HIGH_WORD_EXPONENT_MASK ) { + *rem1 = 0.0 / 0.0; // NaN + *rem2 = 0.0 / 0.0; // NaN + return 0; + } + // Set z = scalbn(|x|, ilogb(x)-23)... + stdlib_base_float64_get_low_word( x, &ulow ); + low = ( int32_t )ulow; + e0 = ( ix >> 20 ) - 1046; // `e0 = ilogb(z) - 23` => unbiased exponent minus 23 + stdlib_base_float64_from_words( ix - ( int32_t ) ( e0 << 20 ), low, &z ); + for ( i = 0; i < 2; i++ ) { + TX[ i ] = ( double )( ( int32_t )z ); + z = ( z - TX[ i ] ) * TWO24; + } + TX[ 2 ] = z; + nx = 3; + while ( TX[ nx - 1 ] == ZERO ) { + // Skip zero term... + nx--; + } + n = rempio2Kernel( TX, TY, e0, nx); + if ( hx < 0 ) { + *rem1 = -TY[ 0 ]; + *rem2 = -TY[ 1 ]; + return -n; + } + *rem1 = TY[ 0 ]; + *rem2 = TY[ 1 ]; + return n; +} diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js b/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js new file mode 100644 index 000000000000..42592a240311 --- /dev/null +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js @@ -0,0 +1,281 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2024 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 resolve = require( 'path' ).resolve; +var tape = require( 'tape' ); +var isInteger = require( '@stdlib/assert/is-integer' ); +var isNumber = require( '@stdlib/assert/is-number' ); +var isnan = require( '@stdlib/math/base/assert/is-nan' ); +var incrspace = require( '@stdlib/array/base/incrspace' ); +var linspace = require( '@stdlib/array/base/linspace' ); +var randu = require( '@stdlib/random/base/randu' ); +var abs = require( '@stdlib/math/base/special/abs' ); +var pow = require( '@stdlib/math/base/special/pow' ); +var PINF = require( '@stdlib/constants/float64/pinf' ); +var NINF = require( '@stdlib/constants/float64/ninf' ); +var EPS = require( '@stdlib/constants/float64/eps' ); +var PIO2 = require( '@stdlib/constants/float64/half-pi' ); +var PIO4 = require( '@stdlib/constants/float64/fourth-pi' ); +var tryRequire = require( '@stdlib/utils/try-require' ); + + +// VARIABLES // + +var rempio2 = tryRequire( resolve( __dirname, './../lib/native.js' ) ); +var opts = { + 'skip': ( rempio2 instanceof Error ) +}; + + +// TESTS // + +tape( 'main export is a function', opts, function test( t ) { + t.ok( true, __filename ); + t.strictEqual( typeof rempio2, 'function', 'main export is a function' ); + t.end(); +}); + +tape( 'the function returns `0` and sets the elements of `y` to `NaN` if provided `NaN`', opts, function test( t ) { + var y = [ 0.0, 0.0 ]; + var n = rempio2( NaN, y ); + t.strictEqual( n, 0.0, 'returns 0' ); + t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); + t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); + t.end(); +}); + +tape( 'the function returns `0` and sets the elements of `y` to `NaN` if provided positive or negative infinity', opts, function test( t ) { + var y = [ 0.0, 0.0 ]; + var n = rempio2( PINF, y ); + t.strictEqual( n, 0.0, 'returns 0' ); + t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); + t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); + + y = [ 0.0, 0.0 ]; + n = rempio2( NINF, y ); + t.strictEqual( n, 0.0, 'returns 0' ); + t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); + t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); + t.end(); +}); + +tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (positive)', opts, function test( t ) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + y = [ 0.0, 0.0 ]; + for ( i = 0; i < 200; i++ ) { + x = randu()*100.0; + n = rempio2( x, y ); + + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + + z = ( PIO2*n ) + ( y[0]+y[1] ); + delta = abs( z - x ); + tol = EPS * abs( x ); + t.strictEqual( delta <= tol, true, 'can recover input value' ); + } + t.end(); +}); + +tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (tiny positive)', opts, function test( t ) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + y = [ 0.0, 0.0 ]; + for ( i = 0; i < 200; i++ ) { + x = randu()*1.0e-100; + n = rempio2( x, y ); + + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + + z = ( PIO2*n ) + ( y[0]+y[1] ); + delta = abs( z - x ); + tol = EPS * abs( x ); + t.strictEqual( delta <= tol, true, 'can recover input value' ); + } + t.end(); +}); + +tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (negative)', opts, function test( t ) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + y = [ 0.0, 0.0 ]; + for ( i = 0; i < 200; i++ ) { + x = -100.0 * randu(); + n = rempio2( x, y ); + + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + + z = ( PIO2*n ) + ( y[0]+y[1] ); + delta = abs( z - x ); + tol = EPS * abs( x ); + t.strictEqual( delta <= tol, true, 'can recover input value' ); + } + t.end(); +}); + +tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (tiny negative)', opts, function test( t ) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + y = [ 0.0, 0.0 ]; + for ( i = 0; i < 200; i++ ) { + x = -1.0e-100 * randu(); + n = rempio2( x, y ); + + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + + z = ( PIO2*n ) + ( y[0]+y[1] ); + delta = abs( z - x ); + tol = EPS * abs( x ); + t.strictEqual( delta <= tol, true, 'can recover input value' ); + } + t.end(); +}); + +tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (multiples of π/4)', opts, function test( t ) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + x = incrspace( -10.0*PIO4, 10.0*PIO4, PIO4 ); + y = [ 0.0, 0.0 ]; + for ( i = 0; i < x.length; i++ ) { + n = rempio2( x[i], y ); + + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + + z = ( PIO2*n ) + ( y[0]+y[1] ); + delta = abs( z - x[i] ); + tol = EPS * abs( x[i] ); + t.strictEqual( delta <= tol, true, 'can recover input value' ); + } + t.end(); +}); + +tape( 'for large positive input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`', opts, function test( t ) { + var x; + var y; + var n; + var i; + + x = linspace( pow( 2.0, 20 )*PIO2, pow( 2.0, 60 )*PIO2, 4000 ); + y = [ 0.0, 0.0 ]; + for ( i = 0; i < x.length; i++ ) { + n = rempio2( x[i], y ); + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + } + t.end(); +}); + +tape( 'for large negative input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`', opts, function test( t ) { + var x; + var y; + var n; + var i; + + x = linspace( -pow( 2.0, 20 )*PIO2, -pow( 2.0, 60 )*PIO2, 4000 ); + y = [ 0.0, 0.0 ]; + for ( i = 0; i < x.length; i++ ) { + n = rempio2( x[i], y ); + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + } + t.end(); +}); + +tape( 'for huge positive input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`', opts, function test( t ) { + var x; + var y; + var n; + var i; + + x = linspace( pow( 2.0, 60 )*PIO2, pow( 2.0, 1020 )*PIO2, 6000 ); + y = [ 0.0, 0.0 ]; + for ( i = 0; i < x.length; i++ ) { + n = rempio2( x[i], y ); + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + } + t.end(); +}); + +tape( 'for huge negative input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`', opts, function test( t ) { + var x; + var y; + var n; + var i; + + x = linspace( -pow( 2.0, 60 )*PIO2, -pow( 2.0, 1020 )*PIO2, 6000 ); + y = [ 0.0, 0.0 ]; + for ( i = 0; i < x.length; i++ ) { + n = rempio2( x[i], y ); + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + } + t.end(); +}); From d216eea14168734ea9ced83c7401c1238b47fef9 Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 17:52:19 +0530 Subject: [PATCH 02/28] feat: draft commit --- .../math/base/special/modf/examples/c/example.c | 15 --------------- .../include/stdlib/math/base/special/rempio2.h | 2 ++ 2 files changed, 2 insertions(+), 15 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/modf/examples/c/example.c b/lib/node_modules/@stdlib/math/base/special/modf/examples/c/example.c index 536f29780128..e974a8225c5c 100644 --- a/lib/node_modules/@stdlib/math/base/special/modf/examples/c/example.c +++ b/lib/node_modules/@stdlib/math/base/special/modf/examples/c/example.c @@ -15,18 +15,3 @@ * See the License for the specific language governing permissions and * limitations under the License. */ - -#include "stdlib/math/base/special/modf.h" -#include - -int main( void ) { - const double x[] = { 4.0, 0.0, -0.0, 1.0, -1.0, 3.14, -3.14, 1.0e308, -1.0e308, 1.0/0.0, -1.0/0.0, 0.0/0.0 }; - - double integral; - double frac; - int i; - for ( i = 0; i < 12; i++ ) { - stdlib_base_modf( x[ i ], &integral, &frac ); - printf( "x: %lf => integral: %lf, frac: %lf\n", x[ i ], integral, frac ); - } -} diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/include/stdlib/math/base/special/rempio2.h b/lib/node_modules/@stdlib/math/base/special/rempio2/include/stdlib/math/base/special/rempio2.h index ea24a327d940..26badf607591 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/include/stdlib/math/base/special/rempio2.h +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/include/stdlib/math/base/special/rempio2.h @@ -22,6 +22,8 @@ #ifndef STDLIB_MATH_BASE_SPECIAL_REMPIO2_H #define STDLIB_MATH_BASE_SPECIAL_REMPIO2_H +#include + /* * If C++, prevent name mangling so that the compiler emits a binary file having undecorated names, thus mirroring the behavior of a C compiler. */ From 1a95200c2d9c51bb80ef500f0c915f3f79995a74 Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 17:54:01 +0530 Subject: [PATCH 03/28] feat: draft commit --- .../base/special/rempio2/examples/c/example.c | 16 ---------------- 1 file changed, 16 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c index 8286358afd0a..e974a8225c5c 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c @@ -15,19 +15,3 @@ * See the License for the specific language governing permissions and * limitations under the License. */ - -#include "stdlib/math/base/special/rempio2.h" -#include - -int main( void ) { - const double x[] = { 0.0, 1.0, 4.0, 128.0 }; - - double rem1; - double rem2; - double n; - int i; - for ( i = 0; i < 4; i++ ) { - n = stdlib_base_rempio2( x[ i ], &rem1, &rem2 ); - printf( "%lf - %lfπ/2 = %lf + %lf", x[ i ], n, rem1, rem2 ); - } -} From 26cd2a16f84d98ba50bf232c0e873e91ce07eb2c Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 17:55:42 +0530 Subject: [PATCH 04/28] feat: draft commit --- .../base/special/rempio2/examples/c/example.c | 15 +++++++++++++++ 1 file changed, 15 insertions(+) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c index e974a8225c5c..ce8efb9ea9da 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c @@ -15,3 +15,18 @@ * See the License for the specific language governing permissions and * limitations under the License. */ + +#include + +int main( void ) { + const double x[] = { 0.0, 1.0, 4.0, 128.0 }; + + double rem1 = 1; + double rem2 = 2; + double n; + int i; + for ( i = 0; i < 4; i++ ) { + n = 0; + printf( "%lf - %lfπ/2 = %lf + %lf", x[ i ], n, rem1, rem2 ); + } +} From 207eff47fd8b2e926594e411c2866bd86540a9a5 Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 17:56:20 +0530 Subject: [PATCH 05/28] feat: draft commit --- .../@stdlib/math/base/special/rempio2/examples/c/example.c | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c index ce8efb9ea9da..8286358afd0a 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c @@ -16,17 +16,18 @@ * limitations under the License. */ +#include "stdlib/math/base/special/rempio2.h" #include int main( void ) { const double x[] = { 0.0, 1.0, 4.0, 128.0 }; - double rem1 = 1; - double rem2 = 2; + double rem1; + double rem2; double n; int i; for ( i = 0; i < 4; i++ ) { - n = 0; + n = stdlib_base_rempio2( x[ i ], &rem1, &rem2 ); printf( "%lf - %lfπ/2 = %lf + %lf", x[ i ], n, rem1, rem2 ); } } From b1f96c83b3aa84f6ff36814aeecdbe56bc04b49f Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 17:58:14 +0530 Subject: [PATCH 06/28] feat: draft commit --- .../math/base/special/modf/examples/c/example.c | 15 +++++++++++++++ 1 file changed, 15 insertions(+) diff --git a/lib/node_modules/@stdlib/math/base/special/modf/examples/c/example.c b/lib/node_modules/@stdlib/math/base/special/modf/examples/c/example.c index e974a8225c5c..536f29780128 100644 --- a/lib/node_modules/@stdlib/math/base/special/modf/examples/c/example.c +++ b/lib/node_modules/@stdlib/math/base/special/modf/examples/c/example.c @@ -15,3 +15,18 @@ * See the License for the specific language governing permissions and * limitations under the License. */ + +#include "stdlib/math/base/special/modf.h" +#include + +int main( void ) { + const double x[] = { 4.0, 0.0, -0.0, 1.0, -1.0, 3.14, -3.14, 1.0e308, -1.0e308, 1.0/0.0, -1.0/0.0, 0.0/0.0 }; + + double integral; + double frac; + int i; + for ( i = 0; i < 12; i++ ) { + stdlib_base_modf( x[ i ], &integral, &frac ); + printf( "x: %lf => integral: %lf, frac: %lf\n", x[ i ], integral, frac ); + } +} From bc0edcaeeb1d01e3f436d8290d18e8caae94d6ca Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 17:58:54 +0530 Subject: [PATCH 07/28] feat: draft commit --- .../@stdlib/math/base/special/rempio2/examples/c/example.c | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c index 8286358afd0a..9e3abdc4e005 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c @@ -16,18 +16,17 @@ * limitations under the License. */ -#include "stdlib/math/base/special/rempio2.h" #include int main( void ) { const double x[] = { 0.0, 1.0, 4.0, 128.0 }; - double rem1; - double rem2; + double rem1 = 1; + double rem2 = 2; double n; int i; for ( i = 0; i < 4; i++ ) { - n = stdlib_base_rempio2( x[ i ], &rem1, &rem2 ); + n = 1; printf( "%lf - %lfπ/2 = %lf + %lf", x[ i ], n, rem1, rem2 ); } } From 78fd73bd5cdcf159da7fde7b49a0bcbe25e9ac7d Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 18:01:40 +0530 Subject: [PATCH 08/28] feat: draft commit --- .../rempio2/benchmark/c/native/benchmark.c | 138 ------------------ .../base/special/rempio2/examples/c/example.c | 32 ---- 2 files changed, 170 deletions(-) delete mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/benchmark.c delete mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/benchmark.c b/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/benchmark.c deleted file mode 100644 index 78439da428ff..000000000000 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/benchmark.c +++ /dev/null @@ -1,138 +0,0 @@ -/** -* @license Apache-2.0 -* -* Copyright (c) 2024 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. -*/ - -/** -* Benchmark `rempio2`. -*/ -#include "stdlib/math/base/special/rempio2.h" -#include -#include -#include -#include -#include - -#define NAME "rempio2" -#define ITERATIONS 1000000 -#define REPEATS 3 - -/** -* Prints the TAP version. -*/ -void print_version() { - printf( "TAP version 13\n" ); -} - -/** -* Prints the TAP summary. -* -* @param total total number of tests -* @param passing total number of passing tests -*/ -void print_summary( int total, int passing ) { - printf( "#\n" ); - printf( "1..%d\n", total ); // TAP plan - printf( "# total %d\n", total ); - printf( "# pass %d\n", passing ); - printf( "#\n" ); - printf( "# ok\n" ); -} - -/** -* Prints benchmarks results. -* -* @param elapsed elapsed time in seconds -*/ -void print_results( double elapsed ) { - double rate = (double)ITERATIONS / elapsed; - printf( " ---\n" ); - printf( " iterations: %d\n", ITERATIONS ); - printf( " elapsed: %0.9f\n", elapsed ); - printf( " rate: %0.9f\n", rate ); - printf( " ...\n" ); -} - -/** -* Returns a clock time. -* -* @return clock time -*/ -double tic() { - struct timeval now; - gettimeofday( &now, NULL ); - return (double)now.tv_sec + (double)now.tv_usec/1.0e6; -} - -/** -* Generates a random number on the interval [0,1]. -* -* @return random number -*/ -double rand_double() { - int r = rand(); - return (double)r / ( (double)RAND_MAX + 1.0 ); -} - -/** -* Runs a benchmark. -* -* @return elapsed time in seconds -*/ -double benchmark() { - double elapsed; - double rem1; - double rem2; - double x; - double y; - double t; - int i; - - t = tic(); - for ( i = 0; i < ITERATIONS; i++ ) { - x = ( rand_double() * 200.0 ) - 100.0; - y = stdlib_base_rempio2( x, &rem1, &rem2 ); - if ( y != y ) { - printf( "should not return NaN\n" ); - break; - } - } - elapsed = tic() - t; - if ( y != y ) { - printf( "should not return NaN\n" ); - } - return elapsed; -} - -/** -* Main execution sequence. -*/ -int main( void ) { - double elapsed; - int i; - - // Use the current time to seed the random number generator: - srand( time( NULL ) ); - - print_version(); - for ( i = 0; i < REPEATS; i++ ) { - printf( "# c::native::%s\n", NAME ); - elapsed = benchmark(); - print_results( elapsed ); - printf( "ok %d benchmark finished\n", i+1 ); - } - print_summary( REPEATS, REPEATS ); -} diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c deleted file mode 100644 index 9e3abdc4e005..000000000000 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c +++ /dev/null @@ -1,32 +0,0 @@ -/** -* @license Apache-2.0 -* -* Copyright (c) 2024 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. -*/ - -#include - -int main( void ) { - const double x[] = { 0.0, 1.0, 4.0, 128.0 }; - - double rem1 = 1; - double rem2 = 2; - double n; - int i; - for ( i = 0; i < 4; i++ ) { - n = 1; - printf( "%lf - %lfπ/2 = %lf + %lf", x[ i ], n, rem1, rem2 ); - } -} From 287b4b416870f6659070e19b6060e22c536d1bdc Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 18:02:43 +0530 Subject: [PATCH 09/28] feat: draft commit --- .../base/special/rempio2/test/test.native.js | 281 ------------------ 1 file changed, 281 deletions(-) delete mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js b/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js deleted file mode 100644 index 42592a240311..000000000000 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js +++ /dev/null @@ -1,281 +0,0 @@ -/** -* @license Apache-2.0 -* -* Copyright (c) 2024 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 resolve = require( 'path' ).resolve; -var tape = require( 'tape' ); -var isInteger = require( '@stdlib/assert/is-integer' ); -var isNumber = require( '@stdlib/assert/is-number' ); -var isnan = require( '@stdlib/math/base/assert/is-nan' ); -var incrspace = require( '@stdlib/array/base/incrspace' ); -var linspace = require( '@stdlib/array/base/linspace' ); -var randu = require( '@stdlib/random/base/randu' ); -var abs = require( '@stdlib/math/base/special/abs' ); -var pow = require( '@stdlib/math/base/special/pow' ); -var PINF = require( '@stdlib/constants/float64/pinf' ); -var NINF = require( '@stdlib/constants/float64/ninf' ); -var EPS = require( '@stdlib/constants/float64/eps' ); -var PIO2 = require( '@stdlib/constants/float64/half-pi' ); -var PIO4 = require( '@stdlib/constants/float64/fourth-pi' ); -var tryRequire = require( '@stdlib/utils/try-require' ); - - -// VARIABLES // - -var rempio2 = tryRequire( resolve( __dirname, './../lib/native.js' ) ); -var opts = { - 'skip': ( rempio2 instanceof Error ) -}; - - -// TESTS // - -tape( 'main export is a function', opts, function test( t ) { - t.ok( true, __filename ); - t.strictEqual( typeof rempio2, 'function', 'main export is a function' ); - t.end(); -}); - -tape( 'the function returns `0` and sets the elements of `y` to `NaN` if provided `NaN`', opts, function test( t ) { - var y = [ 0.0, 0.0 ]; - var n = rempio2( NaN, y ); - t.strictEqual( n, 0.0, 'returns 0' ); - t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); - t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); - t.end(); -}); - -tape( 'the function returns `0` and sets the elements of `y` to `NaN` if provided positive or negative infinity', opts, function test( t ) { - var y = [ 0.0, 0.0 ]; - var n = rempio2( PINF, y ); - t.strictEqual( n, 0.0, 'returns 0' ); - t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); - t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); - - y = [ 0.0, 0.0 ]; - n = rempio2( NINF, y ); - t.strictEqual( n, 0.0, 'returns 0' ); - t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); - t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); - t.end(); -}); - -tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (positive)', opts, function test( t ) { - var delta; - var tol; - var x; - var y; - var z; - var n; - var i; - - y = [ 0.0, 0.0 ]; - for ( i = 0; i < 200; i++ ) { - x = randu()*100.0; - n = rempio2( x, y ); - - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - - z = ( PIO2*n ) + ( y[0]+y[1] ); - delta = abs( z - x ); - tol = EPS * abs( x ); - t.strictEqual( delta <= tol, true, 'can recover input value' ); - } - t.end(); -}); - -tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (tiny positive)', opts, function test( t ) { - var delta; - var tol; - var x; - var y; - var z; - var n; - var i; - - y = [ 0.0, 0.0 ]; - for ( i = 0; i < 200; i++ ) { - x = randu()*1.0e-100; - n = rempio2( x, y ); - - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - - z = ( PIO2*n ) + ( y[0]+y[1] ); - delta = abs( z - x ); - tol = EPS * abs( x ); - t.strictEqual( delta <= tol, true, 'can recover input value' ); - } - t.end(); -}); - -tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (negative)', opts, function test( t ) { - var delta; - var tol; - var x; - var y; - var z; - var n; - var i; - - y = [ 0.0, 0.0 ]; - for ( i = 0; i < 200; i++ ) { - x = -100.0 * randu(); - n = rempio2( x, y ); - - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - - z = ( PIO2*n ) + ( y[0]+y[1] ); - delta = abs( z - x ); - tol = EPS * abs( x ); - t.strictEqual( delta <= tol, true, 'can recover input value' ); - } - t.end(); -}); - -tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (tiny negative)', opts, function test( t ) { - var delta; - var tol; - var x; - var y; - var z; - var n; - var i; - - y = [ 0.0, 0.0 ]; - for ( i = 0; i < 200; i++ ) { - x = -1.0e-100 * randu(); - n = rempio2( x, y ); - - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - - z = ( PIO2*n ) + ( y[0]+y[1] ); - delta = abs( z - x ); - tol = EPS * abs( x ); - t.strictEqual( delta <= tol, true, 'can recover input value' ); - } - t.end(); -}); - -tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (multiples of π/4)', opts, function test( t ) { - var delta; - var tol; - var x; - var y; - var z; - var n; - var i; - - x = incrspace( -10.0*PIO4, 10.0*PIO4, PIO4 ); - y = [ 0.0, 0.0 ]; - for ( i = 0; i < x.length; i++ ) { - n = rempio2( x[i], y ); - - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - - z = ( PIO2*n ) + ( y[0]+y[1] ); - delta = abs( z - x[i] ); - tol = EPS * abs( x[i] ); - t.strictEqual( delta <= tol, true, 'can recover input value' ); - } - t.end(); -}); - -tape( 'for large positive input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`', opts, function test( t ) { - var x; - var y; - var n; - var i; - - x = linspace( pow( 2.0, 20 )*PIO2, pow( 2.0, 60 )*PIO2, 4000 ); - y = [ 0.0, 0.0 ]; - for ( i = 0; i < x.length; i++ ) { - n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); - } - t.end(); -}); - -tape( 'for large negative input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`', opts, function test( t ) { - var x; - var y; - var n; - var i; - - x = linspace( -pow( 2.0, 20 )*PIO2, -pow( 2.0, 60 )*PIO2, 4000 ); - y = [ 0.0, 0.0 ]; - for ( i = 0; i < x.length; i++ ) { - n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); - } - t.end(); -}); - -tape( 'for huge positive input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`', opts, function test( t ) { - var x; - var y; - var n; - var i; - - x = linspace( pow( 2.0, 60 )*PIO2, pow( 2.0, 1020 )*PIO2, 6000 ); - y = [ 0.0, 0.0 ]; - for ( i = 0; i < x.length; i++ ) { - n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); - } - t.end(); -}); - -tape( 'for huge negative input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`', opts, function test( t ) { - var x; - var y; - var n; - var i; - - x = linspace( -pow( 2.0, 60 )*PIO2, -pow( 2.0, 1020 )*PIO2, 6000 ); - y = [ 0.0, 0.0 ]; - for ( i = 0; i < x.length; i++ ) { - n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); - } - t.end(); -}); From 62dd882b31626754d8ec8b4597ef2214817b7433 Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 18:07:32 +0530 Subject: [PATCH 10/28] Create test.native.js Signed-off-by: GUNJ JOSHI --- .../base/special/rempio2/test/test.native.js | 323 ++++++++++++++++++ 1 file changed, 323 insertions(+) create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js b/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js new file mode 100644 index 000000000000..19180df6084c --- /dev/null +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js @@ -0,0 +1,323 @@ +/** + * @license Apache-2.0 + * + * Copyright (c) 2024 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 resolve = require("path").resolve; +var tape = require("tape"); +var isInteger = require("@stdlib/assert/is-integer"); +var isNumber = require("@stdlib/assert/is-number"); +var isnan = require("@stdlib/math/base/assert/is-nan"); +var incrspace = require("@stdlib/array/base/incrspace"); +var linspace = require("@stdlib/array/base/linspace"); +var randu = require("@stdlib/random/base/randu"); +var abs = require("@stdlib/math/base/special/abs"); +var pow = require("@stdlib/math/base/special/pow"); +var PINF = require("@stdlib/constants/float64/pinf"); +var NINF = require("@stdlib/constants/float64/ninf"); +var EPS = require("@stdlib/constants/float64/eps"); +var PIO2 = require("@stdlib/constants/float64/half-pi"); +var PIO4 = require("@stdlib/constants/float64/fourth-pi"); +var tryRequire = require("@stdlib/utils/try-require"); + +// VARIABLES // + +var rempio2 = tryRequire(resolve(__dirname, "./../lib/native.js")); +var opts = { + skip: rempio2 instanceof Error, +}; + +// TESTS // + +tape("main export is a function", opts, function test(t) { + t.ok(true, __filename); + t.strictEqual(typeof rempio2, "function", "main export is a function"); + t.end(); +}); + +tape( + "the function returns `0` and sets the elements of `y` to `NaN` if provided `NaN`", + opts, + function test(t) { + var y = [0.0, 0.0]; + var n = rempio2(NaN, y); + t.strictEqual(n, 0.0, "returns 0"); + t.strictEqual(isnan(y[0]), true, "sets y[0] to NaN"); + t.strictEqual(isnan(y[1]), true, "sets y[1] to NaN"); + t.end(); + } +); + +tape( + "the function returns `0` and sets the elements of `y` to `NaN` if provided positive or negative infinity", + opts, + function test(t) { + var y = [0.0, 0.0]; + var n = rempio2(PINF, y); + t.strictEqual(n, 0.0, "returns 0"); + t.strictEqual(isnan(y[0]), true, "sets y[0] to NaN"); + t.strictEqual(isnan(y[1]), true, "sets y[1] to NaN"); + + y = [0.0, 0.0]; + n = rempio2(NINF, y); + t.strictEqual(n, 0.0, "returns 0"); + t.strictEqual(isnan(y[0]), true, "sets y[0] to NaN"); + t.strictEqual(isnan(y[1]), true, "sets y[1] to NaN"); + t.end(); + } +); + +tape( + "the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (positive)", + opts, + function test(t) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + y = [0.0, 0.0]; + for (i = 0; i < 200; i++) { + x = randu() * 100.0; + n = rempio2(x, y); + + t.strictEqual(isInteger(n), true, "returns integer"); + t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); + t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); + + z = PIO2 * n + (y[0] + y[1]); + delta = abs(z - x); + tol = EPS * abs(x); + t.strictEqual(delta <= tol, true, "can recover input value"); + } + t.end(); + } +); + +tape( + "the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (tiny positive)", + opts, + function test(t) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + y = [0.0, 0.0]; + for (i = 0; i < 200; i++) { + x = randu() * 1.0e-100; + n = rempio2(x, y); + + t.strictEqual(isInteger(n), true, "returns integer"); + t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); + t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); + + z = PIO2 * n + (y[0] + y[1]); + delta = abs(z - x); + tol = EPS * abs(x); + t.strictEqual(delta <= tol, true, "can recover input value"); + } + t.end(); + } +); + +tape( + "the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (negative)", + opts, + function test(t) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + y = [0.0, 0.0]; + for (i = 0; i < 200; i++) { + x = -100.0 * randu(); + n = rempio2(x, y); + + t.strictEqual(isInteger(n), true, "returns integer"); + t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); + t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); + + z = PIO2 * n + (y[0] + y[1]); + delta = abs(z - x); + tol = EPS * abs(x); + t.strictEqual(delta <= tol, true, "can recover input value"); + } + t.end(); + } +); + +tape( + "the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (tiny negative)", + opts, + function test(t) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + y = [0.0, 0.0]; + for (i = 0; i < 200; i++) { + x = -1.0e-100 * randu(); + n = rempio2(x, y); + + t.strictEqual(isInteger(n), true, "returns integer"); + t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); + t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); + + z = PIO2 * n + (y[0] + y[1]); + delta = abs(z - x); + tol = EPS * abs(x); + t.strictEqual(delta <= tol, true, "can recover input value"); + } + t.end(); + } +); + +tape( + "the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (multiples of π/4)", + opts, + function test(t) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + x = incrspace(-10.0 * PIO4, 10.0 * PIO4, PIO4); + y = [0.0, 0.0]; + for (i = 0; i < x.length; i++) { + n = rempio2(x[i], y); + + t.strictEqual(isInteger(n), true, "returns integer"); + t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); + t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); + + z = PIO2 * n + (y[0] + y[1]); + delta = abs(z - x[i]); + tol = EPS * abs(x[i]); + t.strictEqual(delta <= tol, true, "can recover input value"); + } + t.end(); + } +); + +tape( + "for large positive input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`", + opts, + function test(t) { + var x; + var y; + var n; + var i; + + x = linspace(pow(2.0, 20) * PIO2, pow(2.0, 60) * PIO2, 4000); + y = [0.0, 0.0]; + for (i = 0; i < x.length; i++) { + n = rempio2(x[i], y); + t.strictEqual(isInteger(n), true, "returns integer"); + t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); + t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); + t.strictEqual(n < 8, true, "returns a value smaller than 8"); + } + t.end(); + } +); + +tape( + "for large negative input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`", + opts, + function test(t) { + var x; + var y; + var n; + var i; + + x = linspace(-pow(2.0, 20) * PIO2, -pow(2.0, 60) * PIO2, 4000); + y = [0.0, 0.0]; + for (i = 0; i < x.length; i++) { + n = rempio2(x[i], y); + t.strictEqual(isInteger(n), true, "returns integer"); + t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); + t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); + t.strictEqual(n < 8, true, "returns a value smaller than 8"); + } + t.end(); + } +); + +tape( + "for huge positive input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`", + opts, + function test(t) { + var x; + var y; + var n; + var i; + + x = linspace(pow(2.0, 60) * PIO2, pow(2.0, 1020) * PIO2, 6000); + y = [0.0, 0.0]; + for (i = 0; i < x.length; i++) { + n = rempio2(x[i], y); + t.strictEqual(isInteger(n), true, "returns integer"); + t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); + t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); + t.strictEqual(n < 8, true, "returns a value smaller than 8"); + } + t.end(); + } +); + +tape( + "for huge negative input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`", + opts, + function test(t) { + var x; + var y; + var n; + var i; + + x = linspace(-pow(2.0, 60) * PIO2, -pow(2.0, 1020) * PIO2, 6000); + y = [0.0, 0.0]; + for (i = 0; i < x.length; i++) { + n = rempio2(x[i], y); + t.strictEqual(isInteger(n), true, "returns integer"); + t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); + t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); + t.strictEqual(n < 8, true, "returns a value smaller than 8"); + } + t.end(); + } +); From 8a18ebac42265e730798f1b0faf5cf6a305b7b9d Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 18:08:12 +0530 Subject: [PATCH 11/28] Create benchmark.c Signed-off-by: GUNJ JOSHI --- .../rempio2/benchmark/c/native/benchmark.c | 138 ++++++++++++++++++ 1 file changed, 138 insertions(+) create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/benchmark.c diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/benchmark.c b/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/benchmark.c new file mode 100644 index 000000000000..78439da428ff --- /dev/null +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/benchmark/c/native/benchmark.c @@ -0,0 +1,138 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2024 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. +*/ + +/** +* Benchmark `rempio2`. +*/ +#include "stdlib/math/base/special/rempio2.h" +#include +#include +#include +#include +#include + +#define NAME "rempio2" +#define ITERATIONS 1000000 +#define REPEATS 3 + +/** +* Prints the TAP version. +*/ +void print_version() { + printf( "TAP version 13\n" ); +} + +/** +* Prints the TAP summary. +* +* @param total total number of tests +* @param passing total number of passing tests +*/ +void print_summary( int total, int passing ) { + printf( "#\n" ); + printf( "1..%d\n", total ); // TAP plan + printf( "# total %d\n", total ); + printf( "# pass %d\n", passing ); + printf( "#\n" ); + printf( "# ok\n" ); +} + +/** +* Prints benchmarks results. +* +* @param elapsed elapsed time in seconds +*/ +void print_results( double elapsed ) { + double rate = (double)ITERATIONS / elapsed; + printf( " ---\n" ); + printf( " iterations: %d\n", ITERATIONS ); + printf( " elapsed: %0.9f\n", elapsed ); + printf( " rate: %0.9f\n", rate ); + printf( " ...\n" ); +} + +/** +* Returns a clock time. +* +* @return clock time +*/ +double tic() { + struct timeval now; + gettimeofday( &now, NULL ); + return (double)now.tv_sec + (double)now.tv_usec/1.0e6; +} + +/** +* Generates a random number on the interval [0,1]. +* +* @return random number +*/ +double rand_double() { + int r = rand(); + return (double)r / ( (double)RAND_MAX + 1.0 ); +} + +/** +* Runs a benchmark. +* +* @return elapsed time in seconds +*/ +double benchmark() { + double elapsed; + double rem1; + double rem2; + double x; + double y; + double t; + int i; + + t = tic(); + for ( i = 0; i < ITERATIONS; i++ ) { + x = ( rand_double() * 200.0 ) - 100.0; + y = stdlib_base_rempio2( x, &rem1, &rem2 ); + if ( y != y ) { + printf( "should not return NaN\n" ); + break; + } + } + elapsed = tic() - t; + if ( y != y ) { + printf( "should not return NaN\n" ); + } + return elapsed; +} + +/** +* Main execution sequence. +*/ +int main( void ) { + double elapsed; + int i; + + // Use the current time to seed the random number generator: + srand( time( NULL ) ); + + print_version(); + for ( i = 0; i < REPEATS; i++ ) { + printf( "# c::native::%s\n", NAME ); + elapsed = benchmark(); + print_results( elapsed ); + printf( "ok %d benchmark finished\n", i+1 ); + } + print_summary( REPEATS, REPEATS ); +} From bb235e53afe74c181f97842b15859fb0cb4708ec Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 18:08:44 +0530 Subject: [PATCH 12/28] Create example.c Signed-off-by: GUNJ JOSHI --- .../base/special/rempio2/examples/c/example.c | 33 +++++++++++++++++++ 1 file changed, 33 insertions(+) create mode 100644 lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c new file mode 100644 index 000000000000..8286358afd0a --- /dev/null +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c @@ -0,0 +1,33 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2024 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. +*/ + +#include "stdlib/math/base/special/rempio2.h" +#include + +int main( void ) { + const double x[] = { 0.0, 1.0, 4.0, 128.0 }; + + double rem1; + double rem2; + double n; + int i; + for ( i = 0; i < 4; i++ ) { + n = stdlib_base_rempio2( x[ i ], &rem1, &rem2 ); + printf( "%lf - %lfπ/2 = %lf + %lf", x[ i ], n, rem1, rem2 ); + } +} From 12b675f15c233bfbf76b2ea6c4eea3987dc1d88b Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 18:21:37 +0530 Subject: [PATCH 13/28] Update main.c Signed-off-by: GUNJ JOSHI --- lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c index 9ca6f82ad658..1c08be6388e2 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c @@ -281,7 +281,7 @@ double rempio2Kernel( const double* x, const double* y, int32_t e0, int32_t nx ) } // Chop off zero terms... - if ( z == 1.0 ) { + if ( z == 0.0 ) { jz -= 1; q0 -= 24; while ( IQ[ jz ] == 0 ) { From b20b442ed186ef413e92d893a0ad1ae63a927078 Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 18:22:23 +0530 Subject: [PATCH 14/28] Update native.js Signed-off-by: GUNJ JOSHI --- .../@stdlib/math/base/special/rempio2/lib/native.js | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/lib/native.js b/lib/node_modules/@stdlib/math/base/special/rempio2/lib/native.js index 8cd9937c4e7f..8f18bc0ba69b 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/lib/native.js +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/lib/native.js @@ -40,10 +40,10 @@ var addon = require( './../src/addon.node' ); * // returns 81 * * var y1 = y[ 0 ]; -* // returns 0 +* // returns ~0.765 * * var y2 = y[ 1 ]; -* // returns 0 +* // returns ~3.618e-17 */ function rempio2( x ) { var out = new Float64Array( 2 ); From 9a7e62bd03a62687229ff18ebee649a3def6c1f3 Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 19:28:35 +0530 Subject: [PATCH 15/28] removed const from rempio2Kernel parameters Signed-off-by: GUNJ JOSHI --- lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c index 1c08be6388e2..010c54ce10e0 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c @@ -140,7 +140,7 @@ int32_t IQ[ 20 ]; * @param e0 the exponent of `x[0]` (must be <= 16360) * @param nx dimension of `x[]` */ -double rempio2Kernel( const double* x, const double* y, int32_t e0, int32_t nx ) { +double rempio2Kernel( double* x, double* y, int32_t e0, int32_t nx ) { int32_t carry; int32_t jz; int32_t jx; From d2584762a846fc6ce23d247baf3d28b1dafcb80e Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 20:03:44 +0530 Subject: [PATCH 16/28] Update lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json Co-authored-by: Philipp Burckhardt Signed-off-by: GUNJ JOSHI --- .../@stdlib/math/base/special/rempio2/manifest.json | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json b/lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json index a8adc8eb7317..fa8a8c4fb159 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json @@ -58,9 +58,9 @@ "@stdlib/math/base/special/round", "@stdlib/math/base/special/floor", "@stdlib/math/base/special/ldexp", - "@stdlib/number/float64/base/get_high_word", - "@stdlib/number/float64/base/get_low_word", - "@stdlib/number/float64/base/from_words" + "@stdlib/number/float64/base/get-high-word", + "@stdlib/number/float64/base/get-low-word", + "@stdlib/number/float64/base/from-words" ] }, { From 0772434a31b4a40f220f15075ac8cd3f266823ce Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Fri, 31 May 2024 20:03:51 +0530 Subject: [PATCH 17/28] Update lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json Co-authored-by: Philipp Burckhardt Signed-off-by: GUNJ JOSHI --- .../@stdlib/math/base/special/rempio2/manifest.json | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json b/lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json index fa8a8c4fb159..3d071ae81bb3 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json @@ -77,9 +77,9 @@ "@stdlib/math/base/special/round", "@stdlib/math/base/special/floor", "@stdlib/math/base/special/ldexp", - "@stdlib/number/float64/base/get_high_word", - "@stdlib/number/float64/base/get_low_word", - "@stdlib/number/float64/base/from_words" + "@stdlib/number/float64/base/get-high-word", + "@stdlib/number/float64/base/get-low-word", + "@stdlib/number/float64/base/from-words" ] } ] From a6bf83f46069f9062e434e46afe0ca37d82759b5 Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Sat, 1 Jun 2024 09:21:48 +0530 Subject: [PATCH 18/28] Update test.native.js Signed-off-by: GUNJ JOSHI --- .../base/special/rempio2/test/test.native.js | 576 ++++++++---------- 1 file changed, 267 insertions(+), 309 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js b/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js index 19180df6084c..42592a240311 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js @@ -1,323 +1,281 @@ /** - * @license Apache-2.0 - * - * Copyright (c) 2024 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"; +* @license Apache-2.0 +* +* Copyright (c) 2024 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 resolve = require("path").resolve; -var tape = require("tape"); -var isInteger = require("@stdlib/assert/is-integer"); -var isNumber = require("@stdlib/assert/is-number"); -var isnan = require("@stdlib/math/base/assert/is-nan"); -var incrspace = require("@stdlib/array/base/incrspace"); -var linspace = require("@stdlib/array/base/linspace"); -var randu = require("@stdlib/random/base/randu"); -var abs = require("@stdlib/math/base/special/abs"); -var pow = require("@stdlib/math/base/special/pow"); -var PINF = require("@stdlib/constants/float64/pinf"); -var NINF = require("@stdlib/constants/float64/ninf"); -var EPS = require("@stdlib/constants/float64/eps"); -var PIO2 = require("@stdlib/constants/float64/half-pi"); -var PIO4 = require("@stdlib/constants/float64/fourth-pi"); -var tryRequire = require("@stdlib/utils/try-require"); +var resolve = require( 'path' ).resolve; +var tape = require( 'tape' ); +var isInteger = require( '@stdlib/assert/is-integer' ); +var isNumber = require( '@stdlib/assert/is-number' ); +var isnan = require( '@stdlib/math/base/assert/is-nan' ); +var incrspace = require( '@stdlib/array/base/incrspace' ); +var linspace = require( '@stdlib/array/base/linspace' ); +var randu = require( '@stdlib/random/base/randu' ); +var abs = require( '@stdlib/math/base/special/abs' ); +var pow = require( '@stdlib/math/base/special/pow' ); +var PINF = require( '@stdlib/constants/float64/pinf' ); +var NINF = require( '@stdlib/constants/float64/ninf' ); +var EPS = require( '@stdlib/constants/float64/eps' ); +var PIO2 = require( '@stdlib/constants/float64/half-pi' ); +var PIO4 = require( '@stdlib/constants/float64/fourth-pi' ); +var tryRequire = require( '@stdlib/utils/try-require' ); + // VARIABLES // -var rempio2 = tryRequire(resolve(__dirname, "./../lib/native.js")); +var rempio2 = tryRequire( resolve( __dirname, './../lib/native.js' ) ); var opts = { - skip: rempio2 instanceof Error, + 'skip': ( rempio2 instanceof Error ) }; + // TESTS // -tape("main export is a function", opts, function test(t) { - t.ok(true, __filename); - t.strictEqual(typeof rempio2, "function", "main export is a function"); - t.end(); +tape( 'main export is a function', opts, function test( t ) { + t.ok( true, __filename ); + t.strictEqual( typeof rempio2, 'function', 'main export is a function' ); + t.end(); +}); + +tape( 'the function returns `0` and sets the elements of `y` to `NaN` if provided `NaN`', opts, function test( t ) { + var y = [ 0.0, 0.0 ]; + var n = rempio2( NaN, y ); + t.strictEqual( n, 0.0, 'returns 0' ); + t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); + t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); + t.end(); +}); + +tape( 'the function returns `0` and sets the elements of `y` to `NaN` if provided positive or negative infinity', opts, function test( t ) { + var y = [ 0.0, 0.0 ]; + var n = rempio2( PINF, y ); + t.strictEqual( n, 0.0, 'returns 0' ); + t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); + t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); + + y = [ 0.0, 0.0 ]; + n = rempio2( NINF, y ); + t.strictEqual( n, 0.0, 'returns 0' ); + t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); + t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); + t.end(); +}); + +tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (positive)', opts, function test( t ) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + y = [ 0.0, 0.0 ]; + for ( i = 0; i < 200; i++ ) { + x = randu()*100.0; + n = rempio2( x, y ); + + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + + z = ( PIO2*n ) + ( y[0]+y[1] ); + delta = abs( z - x ); + tol = EPS * abs( x ); + t.strictEqual( delta <= tol, true, 'can recover input value' ); + } + t.end(); +}); + +tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (tiny positive)', opts, function test( t ) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + y = [ 0.0, 0.0 ]; + for ( i = 0; i < 200; i++ ) { + x = randu()*1.0e-100; + n = rempio2( x, y ); + + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + + z = ( PIO2*n ) + ( y[0]+y[1] ); + delta = abs( z - x ); + tol = EPS * abs( x ); + t.strictEqual( delta <= tol, true, 'can recover input value' ); + } + t.end(); +}); + +tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (negative)', opts, function test( t ) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + y = [ 0.0, 0.0 ]; + for ( i = 0; i < 200; i++ ) { + x = -100.0 * randu(); + n = rempio2( x, y ); + + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + + z = ( PIO2*n ) + ( y[0]+y[1] ); + delta = abs( z - x ); + tol = EPS * abs( x ); + t.strictEqual( delta <= tol, true, 'can recover input value' ); + } + t.end(); }); -tape( - "the function returns `0` and sets the elements of `y` to `NaN` if provided `NaN`", - opts, - function test(t) { - var y = [0.0, 0.0]; - var n = rempio2(NaN, y); - t.strictEqual(n, 0.0, "returns 0"); - t.strictEqual(isnan(y[0]), true, "sets y[0] to NaN"); - t.strictEqual(isnan(y[1]), true, "sets y[1] to NaN"); - t.end(); - } -); - -tape( - "the function returns `0` and sets the elements of `y` to `NaN` if provided positive or negative infinity", - opts, - function test(t) { - var y = [0.0, 0.0]; - var n = rempio2(PINF, y); - t.strictEqual(n, 0.0, "returns 0"); - t.strictEqual(isnan(y[0]), true, "sets y[0] to NaN"); - t.strictEqual(isnan(y[1]), true, "sets y[1] to NaN"); - - y = [0.0, 0.0]; - n = rempio2(NINF, y); - t.strictEqual(n, 0.0, "returns 0"); - t.strictEqual(isnan(y[0]), true, "sets y[0] to NaN"); - t.strictEqual(isnan(y[1]), true, "sets y[1] to NaN"); - t.end(); - } -); - -tape( - "the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (positive)", - opts, - function test(t) { - var delta; - var tol; - var x; - var y; - var z; - var n; - var i; - - y = [0.0, 0.0]; - for (i = 0; i < 200; i++) { - x = randu() * 100.0; - n = rempio2(x, y); - - t.strictEqual(isInteger(n), true, "returns integer"); - t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); - t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); - - z = PIO2 * n + (y[0] + y[1]); - delta = abs(z - x); - tol = EPS * abs(x); - t.strictEqual(delta <= tol, true, "can recover input value"); - } - t.end(); - } -); - -tape( - "the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (tiny positive)", - opts, - function test(t) { - var delta; - var tol; - var x; - var y; - var z; - var n; - var i; - - y = [0.0, 0.0]; - for (i = 0; i < 200; i++) { - x = randu() * 1.0e-100; - n = rempio2(x, y); - - t.strictEqual(isInteger(n), true, "returns integer"); - t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); - t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); - - z = PIO2 * n + (y[0] + y[1]); - delta = abs(z - x); - tol = EPS * abs(x); - t.strictEqual(delta <= tol, true, "can recover input value"); - } - t.end(); - } -); - -tape( - "the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (negative)", - opts, - function test(t) { - var delta; - var tol; - var x; - var y; - var z; - var n; - var i; - - y = [0.0, 0.0]; - for (i = 0; i < 200; i++) { - x = -100.0 * randu(); - n = rempio2(x, y); - - t.strictEqual(isInteger(n), true, "returns integer"); - t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); - t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); - - z = PIO2 * n + (y[0] + y[1]); - delta = abs(z - x); - tol = EPS * abs(x); - t.strictEqual(delta <= tol, true, "can recover input value"); - } - t.end(); - } -); - -tape( - "the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (tiny negative)", - opts, - function test(t) { - var delta; - var tol; - var x; - var y; - var z; - var n; - var i; - - y = [0.0, 0.0]; - for (i = 0; i < 200; i++) { - x = -1.0e-100 * randu(); - n = rempio2(x, y); - - t.strictEqual(isInteger(n), true, "returns integer"); - t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); - t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); - - z = PIO2 * n + (y[0] + y[1]); - delta = abs(z - x); - tol = EPS * abs(x); - t.strictEqual(delta <= tol, true, "can recover input value"); - } - t.end(); - } -); - -tape( - "the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (multiples of π/4)", - opts, - function test(t) { - var delta; - var tol; - var x; - var y; - var z; - var n; - var i; - - x = incrspace(-10.0 * PIO4, 10.0 * PIO4, PIO4); - y = [0.0, 0.0]; - for (i = 0; i < x.length; i++) { - n = rempio2(x[i], y); - - t.strictEqual(isInteger(n), true, "returns integer"); - t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); - t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); - - z = PIO2 * n + (y[0] + y[1]); - delta = abs(z - x[i]); - tol = EPS * abs(x[i]); - t.strictEqual(delta <= tol, true, "can recover input value"); - } - t.end(); - } -); - -tape( - "for large positive input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`", - opts, - function test(t) { - var x; - var y; - var n; - var i; - - x = linspace(pow(2.0, 20) * PIO2, pow(2.0, 60) * PIO2, 4000); - y = [0.0, 0.0]; - for (i = 0; i < x.length; i++) { - n = rempio2(x[i], y); - t.strictEqual(isInteger(n), true, "returns integer"); - t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); - t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); - t.strictEqual(n < 8, true, "returns a value smaller than 8"); - } - t.end(); - } -); - -tape( - "for large negative input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`", - opts, - function test(t) { - var x; - var y; - var n; - var i; - - x = linspace(-pow(2.0, 20) * PIO2, -pow(2.0, 60) * PIO2, 4000); - y = [0.0, 0.0]; - for (i = 0; i < x.length; i++) { - n = rempio2(x[i], y); - t.strictEqual(isInteger(n), true, "returns integer"); - t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); - t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); - t.strictEqual(n < 8, true, "returns a value smaller than 8"); - } - t.end(); - } -); - -tape( - "for huge positive input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`", - opts, - function test(t) { - var x; - var y; - var n; - var i; - - x = linspace(pow(2.0, 60) * PIO2, pow(2.0, 1020) * PIO2, 6000); - y = [0.0, 0.0]; - for (i = 0; i < x.length; i++) { - n = rempio2(x[i], y); - t.strictEqual(isInteger(n), true, "returns integer"); - t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); - t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); - t.strictEqual(n < 8, true, "returns a value smaller than 8"); - } - t.end(); - } -); - -tape( - "for huge negative input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`", - opts, - function test(t) { - var x; - var y; - var n; - var i; - - x = linspace(-pow(2.0, 60) * PIO2, -pow(2.0, 1020) * PIO2, 6000); - y = [0.0, 0.0]; - for (i = 0; i < x.length; i++) { - n = rempio2(x[i], y); - t.strictEqual(isInteger(n), true, "returns integer"); - t.strictEqual(isNumber(y[0]), true, "assigns a number value to y[0]"); - t.strictEqual(isNumber(y[1]), true, "assigns a number value to y[1]"); - t.strictEqual(n < 8, true, "returns a value smaller than 8"); - } - t.end(); - } -); +tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (tiny negative)', opts, function test( t ) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + y = [ 0.0, 0.0 ]; + for ( i = 0; i < 200; i++ ) { + x = -1.0e-100 * randu(); + n = rempio2( x, y ); + + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + + z = ( PIO2*n ) + ( y[0]+y[1] ); + delta = abs( z - x ); + tol = EPS * abs( x ); + t.strictEqual( delta <= tol, true, 'can recover input value' ); + } + t.end(); +}); + +tape( 'the function returns `n` and stores `r` as two double-precision floating points numbers in `y` such that `x - nπ/2 = r` (multiples of π/4)', opts, function test( t ) { + var delta; + var tol; + var x; + var y; + var z; + var n; + var i; + + x = incrspace( -10.0*PIO4, 10.0*PIO4, PIO4 ); + y = [ 0.0, 0.0 ]; + for ( i = 0; i < x.length; i++ ) { + n = rempio2( x[i], y ); + + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + + z = ( PIO2*n ) + ( y[0]+y[1] ); + delta = abs( z - x[i] ); + tol = EPS * abs( x[i] ); + t.strictEqual( delta <= tol, true, 'can recover input value' ); + } + t.end(); +}); + +tape( 'for large positive input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`', opts, function test( t ) { + var x; + var y; + var n; + var i; + + x = linspace( pow( 2.0, 20 )*PIO2, pow( 2.0, 60 )*PIO2, 4000 ); + y = [ 0.0, 0.0 ]; + for ( i = 0; i < x.length; i++ ) { + n = rempio2( x[i], y ); + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + } + t.end(); +}); + +tape( 'for large negative input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`', opts, function test( t ) { + var x; + var y; + var n; + var i; + + x = linspace( -pow( 2.0, 20 )*PIO2, -pow( 2.0, 60 )*PIO2, 4000 ); + y = [ 0.0, 0.0 ]; + for ( i = 0; i < x.length; i++ ) { + n = rempio2( x[i], y ); + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + } + t.end(); +}); + +tape( 'for huge positive input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`', opts, function test( t ) { + var x; + var y; + var n; + var i; + + x = linspace( pow( 2.0, 60 )*PIO2, pow( 2.0, 1020 )*PIO2, 6000 ); + y = [ 0.0, 0.0 ]; + for ( i = 0; i < x.length; i++ ) { + n = rempio2( x[i], y ); + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + } + t.end(); +}); + +tape( 'for huge negative input values, the function returns the last three binary digits of `n` and stores `r` as two double-precision floating points numbers in `y`', opts, function test( t ) { + var x; + var y; + var n; + var i; + + x = linspace( -pow( 2.0, 60 )*PIO2, -pow( 2.0, 1020 )*PIO2, 6000 ); + y = [ 0.0, 0.0 ]; + for ( i = 0; i < x.length; i++ ) { + n = rempio2( x[i], y ); + t.strictEqual( isInteger( n ), true, 'returns integer' ); + t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); + t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + } + t.end(); +}); From f47d9e35361ffdc5f4cd9d1d450be7a9ce3df878 Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Sat, 1 Jun 2024 09:47:39 +0530 Subject: [PATCH 19/28] Update main.c Signed-off-by: GUNJ JOSHI --- lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c index 010c54ce10e0..227fb894ae16 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c @@ -140,7 +140,7 @@ int32_t IQ[ 20 ]; * @param e0 the exponent of `x[0]` (must be <= 16360) * @param nx dimension of `x[]` */ -double rempio2Kernel( double* x, double* y, int32_t e0, int32_t nx ) { +double rempio2Kernel( const double* x, double* y, int32_t e0, int32_t nx ) { int32_t carry; int32_t jz; int32_t jx; From 0df78ac6dd54e79bdb066f7196fe127a92e0ec29 Mon Sep 17 00:00:00 2001 From: Athan Reines Date: Fri, 7 Jun 2024 02:37:06 -0700 Subject: [PATCH 20/28] fix: refactor addon, update docs, and update tests --- .../math/base/special/rempio2/README.md | 29 ++- .../base/special/rempio2/examples/c/example.c | 6 +- .../stdlib/math/base/special/rempio2.h | 2 +- .../math/base/special/rempio2/lib/native.js | 12 +- .../math/base/special/rempio2/manifest.json | 182 ++++++++++-------- .../math/base/special/rempio2/src/addon.c | 98 ++-------- .../math/base/special/rempio2/src/main.c | 14 +- .../math/base/special/rempio2/test/test.js | 90 ++++----- .../base/special/rempio2/test/test.native.js | 105 +++++----- 9 files changed, 245 insertions(+), 293 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/README.md b/lib/node_modules/@stdlib/math/base/special/rempio2/README.md index 8c32bb0e53c5..978033ba6f65 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/README.md +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/README.md @@ -32,7 +32,7 @@ var rempio2 = require( '@stdlib/math/base/special/rempio2' ); #### rempio2( x, y ) -Computes `x - nπ/2 = r`. The function returns `n` and stores the remainder `r` as two numbers in `y`, such that `y[0]+y[1] = r`. +Computes `x - nπ/2 = r`. ```javascript var y = [ 0.0, 0.0 ]; @@ -46,7 +46,7 @@ var y2 = y[ 1 ]; // returns ~3.618e-17 ``` -When `x` is `NaN` or infinite, the function returns zero and sets the elements of `y` to `NaN`. +When `x` is `NaN` or infinite, the function returns `0` and sets the elements of `y` to `NaN`. ```javascript var y = [ 0.0, 0.0 ]; @@ -80,6 +80,7 @@ y2 = y[ 1 ]; ## Notes +- The function returns `n` and stores the remainder `r` as two numbers in `y`, such that `y[0]+y[1] = r`. - For input values larger than `2^20*π/2` in magnitude, the function **only** returns the last three binary digits of `n` and not the full result.
@@ -139,23 +140,25 @@ for ( i = 0; i < x.length; i++ ) { #### stdlib_base_rempio2( x, &rem1, &rem2 ) -Computes `x - nπ/2 = r`. The function returns `n` and stores the remainder `r` as two numbers in `y`, such that `y[0]+y[1] = r`. +Computes `x - nπ/2 = r`. ```c +#include + double rem1; double rem2; -stdlib_base_rempio2( 4.0, &rem1, &rem2 ); +int32_t n = stdlib_base_rempio2( 4.0, &rem1, &rem2 ); ``` The function accepts the following arguments: - **x**: `[in] double` input value. -- **rem1**: `[out] double*` destination for rem1. -- **rem2**: `[out] double*` destination for rem2. +- **rem1**: `[out] double*` destination for first remainder number. +- **rem2**: `[out] double*` destination for second remainder number. ```c -void stdlib_base_rempio2( const double x, double *rem1, double *rem2 ); +int32_t stdlib_base_rempio2( const double x, double *rem1, double *rem2 ); ``` @@ -166,6 +169,10 @@ void stdlib_base_rempio2( const double x, double *rem1, double *rem2 );
+### Notes + +- The function returns `n` and stores the remainder `r` as two numbers in `rem1` and `rem2`, respectively, such that `rem1+rem2 = r`. +
@@ -179,17 +186,19 @@ void stdlib_base_rempio2( const double x, double *rem1, double *rem2 ); ```c #include "stdlib/math/base/special/rempio2.h" #include +#include +#include int main( void ) { const double x[] = { 0.0, 1.0, 4.0, 128.0 }; double rem1; double rem2; - double n; + int32_t n; int i; for ( i = 0; i < 4; i++ ) { - stdlib_base_rempio2( x[ i ], &rem1, &rem2 ); - printf( "%lf - %lfπ/2 = %lf + %lf", x[ i ], n, rem1, rem2 ); + n = stdlib_base_rempio2( x[ i ], &rem1, &rem2 ); + printf( "%lf - %"PRId32"π/2 = %lf + %lf\n", x[ i ], n, rem1, rem2 ); } } ``` diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c index 8286358afd0a..4b314d364054 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/examples/c/example.c @@ -18,16 +18,18 @@ #include "stdlib/math/base/special/rempio2.h" #include +#include +#include int main( void ) { const double x[] = { 0.0, 1.0, 4.0, 128.0 }; double rem1; double rem2; - double n; + int32_t n; int i; for ( i = 0; i < 4; i++ ) { n = stdlib_base_rempio2( x[ i ], &rem1, &rem2 ); - printf( "%lf - %lfπ/2 = %lf + %lf", x[ i ], n, rem1, rem2 ); + printf( "%lf - %"PRId32"π/2 = %lf + %lf\n", x[ i ], n, rem1, rem2 ); } } diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/include/stdlib/math/base/special/rempio2.h b/lib/node_modules/@stdlib/math/base/special/rempio2/include/stdlib/math/base/special/rempio2.h index 26badf607591..dcba0ab4172a 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/include/stdlib/math/base/special/rempio2.h +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/include/stdlib/math/base/special/rempio2.h @@ -34,7 +34,7 @@ extern "C" { /** * Computes `x - nπ/2 = r`. */ -double stdlib_base_rempio2( const double x, double* rem1, double* rem2 ); +int32_t stdlib_base_rempio2( const double x, double *rem1, double *rem2 ); #ifdef __cplusplus } diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/lib/native.js b/lib/node_modules/@stdlib/math/base/special/rempio2/lib/native.js index 8f18bc0ba69b..ab9ea199ea87 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/lib/native.js +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/lib/native.js @@ -20,7 +20,6 @@ // MODULES // -var Float64Array = require( '@stdlib/array/float64' ); var addon = require( './../src/addon.node' ); @@ -31,11 +30,13 @@ var addon = require( './../src/addon.node' ); * * @private * @param {number} x - input value -* @param {(Array|TypedArray|Object)} y - remainder elements +* @param {Float64Array} y - remainder elements * @returns {integer} factor of `π/2` * * @example -* var y = [ 0.0, 0.0 ]; +* var Float64Array = require( '@stdlib/array/float64' ); +* +* var y = new Float64Array( [ 0.0, 0.0 ] ); * var n = rempio2( 128.0, y ); * // returns 81 * @@ -45,9 +46,8 @@ var addon = require( './../src/addon.node' ); * var y2 = y[ 1 ]; * // returns ~3.618e-17 */ -function rempio2( x ) { - var out = new Float64Array( 2 ); - return addon( out, x ); +function rempio2( x, y ) { + return addon( x, y ); } diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json b/lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json index 3d071ae81bb3..11524178fb99 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/manifest.json @@ -1,86 +1,100 @@ { - "options": { - "task": "build" - }, - "fields": [ - { - "field": "src", - "resolve": true, - "relative": true - }, - { - "field": "include", - "resolve": true, - "relative": true - }, - { - "field": "libraries", - "resolve": false, - "relative": false - }, - { - "field": "libpath", - "resolve": true, - "relative": false - } - ], - "confs": [ - { - "task": "build", - "src": [ - "./src/main.c" - ], - "include": [ - "./include" - ], - "libraries": [], - "libpath": [], - "dependencies": [ - "@stdlib/math/base/special/round", - "@stdlib/math/base/special/floor", - "@stdlib/math/base/special/ldexp", - "@stdlib/number/float64/base/get-high-word", - "@stdlib/number/float64/base/get-low-word", - "@stdlib/number/float64/base/from-words" - ] - }, - { - "task": "benchmark", - "src": [ - "./src/main.c" - ], - "include": [ - "./include" - ], - "libraries": [], - "libpath": [], - "dependencies": [ - "@stdlib/math/base/special/round", - "@stdlib/math/base/special/floor", - "@stdlib/math/base/special/ldexp", - "@stdlib/number/float64/base/get-high-word", - "@stdlib/number/float64/base/get-low-word", - "@stdlib/number/float64/base/from-words" - ] - }, - { - "task": "examples", - "src": [ - "./src/main.c" - ], - "include": [ - "./include" - ], - "libraries": [], - "libpath": [], - "dependencies": [ - "@stdlib/math/base/special/round", - "@stdlib/math/base/special/floor", - "@stdlib/math/base/special/ldexp", - "@stdlib/number/float64/base/get-high-word", - "@stdlib/number/float64/base/get-low-word", - "@stdlib/number/float64/base/from-words" - ] - } - ] + "options": { + "task": "build" + }, + "fields": [ + { + "field": "src", + "resolve": true, + "relative": true + }, + { + "field": "include", + "resolve": true, + "relative": true + }, + { + "field": "libraries", + "resolve": false, + "relative": false + }, + { + "field": "libpath", + "resolve": true, + "relative": false + } + ], + "confs": [ + { + "task": "build", + "src": [ + "./src/main.c" + ], + "include": [ + "./include" + ], + "libraries": [], + "libpath": [], + "dependencies": [ + "@stdlib/napi/argv", + "@stdlib/napi/argv-double", + "@stdlib/napi/argv-float64array", + "@stdlib/napi/create-double", + "@stdlib/napi/export", + "@stdlib/math/base/special/round", + "@stdlib/math/base/special/floor", + "@stdlib/math/base/special/ldexp", + "@stdlib/number/float64/base/get-high-word", + "@stdlib/number/float64/base/get-low-word", + "@stdlib/number/float64/base/from-words", + "@stdlib/constants/float64/high-word-abs-mask", + "@stdlib/constants/float64/high-word-exponent-mask", + "@stdlib/constants/float64/high-word-sign-mask" + ] + }, + { + "task": "benchmark", + "src": [ + "./src/main.c" + ], + "include": [ + "./include" + ], + "libraries": [], + "libpath": [], + "dependencies": [ + "@stdlib/math/base/special/round", + "@stdlib/math/base/special/floor", + "@stdlib/math/base/special/ldexp", + "@stdlib/number/float64/base/get-high-word", + "@stdlib/number/float64/base/get-low-word", + "@stdlib/number/float64/base/from-words", + "@stdlib/constants/float64/high-word-abs-mask", + "@stdlib/constants/float64/high-word-exponent-mask", + "@stdlib/constants/float64/high-word-sign-mask" + ] + }, + { + "task": "examples", + "src": [ + "./src/main.c" + ], + "include": [ + "./include" + ], + "libraries": [], + "libpath": [], + "dependencies": [ + "@stdlib/math/base/special/round", + "@stdlib/math/base/special/floor", + "@stdlib/math/base/special/ldexp", + "@stdlib/number/float64/base/get-high-word", + "@stdlib/number/float64/base/get-low-word", + "@stdlib/number/float64/base/from-words", + "@stdlib/constants/float64/high-word-abs-mask", + "@stdlib/constants/float64/high-word-exponent-mask", + "@stdlib/constants/float64/high-word-sign-mask" + ] + } + ] } diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/src/addon.c b/lib/node_modules/@stdlib/math/base/special/rempio2/src/addon.c index 6787d85bdca8..daf621386400 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/src/addon.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/src/addon.c @@ -17,9 +17,12 @@ */ #include "stdlib/math/base/special/rempio2.h" +#include "stdlib/napi/argv.h" +#include "stdlib/napi/argv_double.h" +#include "stdlib/napi/argv_float64array.h" +#include "stdlib/napi/create_double.h" +#include "stdlib/napi/export.h" #include -#include -#include /** * Receives JavaScript callback invocation data. @@ -30,90 +33,11 @@ * @return Node-API value */ static napi_value addon( napi_env env, napi_callback_info info ) { - napi_status status; - - // Get callback arguments: - size_t argc = 2; - napi_value argv[ 2 ]; - status = napi_get_cb_info( env, info, &argc, argv, NULL, NULL ); - assert( status == napi_ok ); - - // Check whether we were provided the correct number of arguments: - if ( argc < 2 ) { - status = napi_throw_error( env, NULL, "invalid invocation. Insufficient arguments." ); - assert( status == napi_ok ); - return NULL; - } - if ( argc > 2 ) { - status = napi_throw_error( env, NULL, "invalid invocation. Too many arguments." ); - assert( status == napi_ok ); - return NULL; - } - - bool res; - status = napi_is_typedarray( env, argv[ 0 ], &res ); - assert( status == napi_ok ); - if ( res == false ) { - status = napi_throw_type_error( env, NULL, "invalid argument. First argument must be a Float64Array." ); - assert( status == napi_ok ); - return NULL; - } - - napi_valuetype vtype1; - status = napi_typeof( env, argv[ 1 ], &vtype1 ); - assert( status == napi_ok ); - if ( vtype1 != napi_number ) { - status = napi_throw_type_error( env, NULL, "invalid argument. Second argument must be a number." ); - assert( status == napi_ok ); - return NULL; - } - - napi_typedarray_type vtype0; - size_t len; - void *Out; - status = napi_get_typedarray_info( env, argv[ 0 ], &vtype0, &len, &Out, NULL, NULL ); - assert( status == napi_ok ); - if ( vtype0 != napi_float64_array ) { - status = napi_throw_type_error( env, NULL, "invalid argument. First argument must be a Float64Array." ); - assert( status == napi_ok ); - return NULL; - } - if ( len != 2 ) { - status = napi_throw_range_error( env, NULL, "invalid argument. First argument must have 2 elements." ); - assert( status == napi_ok ); - return NULL; - } - - double value; - status = napi_get_value_double( env, argv[ 1 ], &value ); - assert( status == napi_ok ); - - double rem1; - double rem2; - napi_value v; - status = napi_create_double( env, stdlib_base_rempio2( value, &rem1, &rem2 ), &v ); - assert( status == napi_ok ); - - double *op = (double *)Out; - op[ 0 ] = rem1; - op[ 1 ] = rem2; - - return v; -} - -/** -* Initializes a Node-API module. -* -* @private -* @param env environment under which the function is invoked -* @param exports exports object -* @return main export -*/ -static napi_value init( napi_env env, napi_value exports ) { - napi_value fcn; - napi_status status = napi_create_function( env, "exports", NAPI_AUTO_LENGTH, addon, NULL, &fcn ); - assert( status == napi_ok ); - return fcn; + STDLIB_NAPI_ARGV( env, info, argv, argc, 2 ); + STDLIB_NAPI_ARGV_DOUBLE( env, x, argv, 0 ); + STDLIB_NAPI_ARGV_FLOAT64ARRAY( env, y, ylen, argv, 1 ); + STDLIB_NAPI_CREATE_DOUBLE( env, stdlib_base_rempio2( x, &y[0], &y[1] ), n ); + return n; } -NAPI_MODULE( NODE_GYP_MODULE_NAME, init ) +STDLIB_NAPI_MODULE_EXPORT_FCN( addon ) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c index 227fb894ae16..3b938164f352 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c @@ -32,11 +32,11 @@ #include "stdlib/math/base/special/rempio2.h" #include "stdlib/math/base/special/round.h" +#include "stdlib/math/base/special/floor.h" +#include "stdlib/math/base/special/ldexp.h" #include "stdlib/number/float64/base/get_high_word.h" #include "stdlib/number/float64/base/get_low_word.h" #include "stdlib/number/float64/base/from_words.h" -#include "stdlib/math/base/special/floor.h" -#include "stdlib/math/base/special/ldexp.h" #include "stdlib/constants/float64/high_word_abs_mask.h" #include "stdlib/constants/float64/high_word_exponent_mask.h" #include "stdlib/constants/float64/high_word_sign_mask.h" @@ -139,8 +139,9 @@ int32_t IQ[ 20 ]; * @param rem2 remainder element 2 * @param e0 the exponent of `x[0]` (must be <= 16360) * @param nx dimension of `x[]` +* @return last 3 binary digits */ -double rempio2Kernel( const double* x, double* y, int32_t e0, int32_t nx ) { +int32_t rempio2Kernel( const double* x, double *y, int32_t e0, int32_t nx ) { int32_t carry; int32_t jz; int32_t jx; @@ -350,15 +351,16 @@ double rempio2Kernel( const double* x, double* y, int32_t e0, int32_t nx ) { * @param x input value * @param rem1 remainder element 1 * @param rem2 remainder element 2 +* @return factor of `π/2` * * @example -* double x = 128.0; +* #include * * double rem1; * double rem2; -* stdlib_base_rempio2( x, &rem1, &rem2 ); +* int32_t n = stdlib_base_rempio2( 128.0, &rem1, &rem2 ); */ -double stdlib_base_rempio2( const double x, double* rem1, double* rem2 ) { +int32_t stdlib_base_rempio2( const double x, double *rem1, double *rem2 ) { uint32_t uhigh; uint32_t ulow; uint32_t uhx; diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.js b/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.js index 5ff347e0d64b..8db709e784e3 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.js +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.js @@ -48,24 +48,24 @@ tape( 'main export is a function', function test( t ) { tape( 'the function returns `0` and sets the elements of `y` to `NaN` if provided `NaN`', function test( t ) { var y = [ 0.0, 0.0 ]; var n = rempio2( NaN, y ); - t.strictEqual( n, 0.0, 'returns 0' ); - t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); - t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); + t.strictEqual( n, 0.0, 'returns expected value' ); + t.strictEqual( isnan( y[0] ), true, 'returns expected value' ); + t.strictEqual( isnan( y[1] ), true, 'returns expected value' ); t.end(); }); tape( 'the function returns `0` and sets the elements of `y` to `NaN` if provided positive or negative infinity', function test( t ) { var y = [ 0.0, 0.0 ]; var n = rempio2( PINF, y ); - t.strictEqual( n, 0.0, 'returns 0' ); - t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); - t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); + t.strictEqual( n, 0.0, 'returns expected value' ); + t.strictEqual( isnan( y[0] ), true, 'returns expected value' ); + t.strictEqual( isnan( y[1] ), true, 'returns expected value' ); y = [ 0.0, 0.0 ]; n = rempio2( NINF, y ); - t.strictEqual( n, 0.0, 'returns 0' ); - t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); - t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); + t.strictEqual( n, 0.0, 'returns expected value' ); + t.strictEqual( isnan( y[0] ), true, 'returns expected value' ); + t.strictEqual( isnan( y[1] ), true, 'returns expected value' ); t.end(); }); @@ -83,14 +83,14 @@ tape( 'the function returns `n` and stores `r` as two double-precision floating x = randu()*100.0; n = rempio2( x, y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); z = ( PIO2*n ) + ( y[0]+y[1] ); delta = abs( z - x ); tol = EPS * abs( x ); - t.strictEqual( delta <= tol, true, 'can recover input value' ); + t.strictEqual( delta <= tol, true, 'returns expected value' ); } t.end(); }); @@ -109,14 +109,14 @@ tape( 'the function returns `n` and stores `r` as two double-precision floating x = randu()*1.0e-100; n = rempio2( x, y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); z = ( PIO2*n ) + ( y[0]+y[1] ); delta = abs( z - x ); tol = EPS * abs( x ); - t.strictEqual( delta <= tol, true, 'can recover input value' ); + t.strictEqual( delta <= tol, true, 'returns expected value' ); } t.end(); }); @@ -135,14 +135,14 @@ tape( 'the function returns `n` and stores `r` as two double-precision floating x = -100.0 * randu(); n = rempio2( x, y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); z = ( PIO2*n ) + ( y[0]+y[1] ); delta = abs( z - x ); tol = EPS * abs( x ); - t.strictEqual( delta <= tol, true, 'can recover input value' ); + t.strictEqual( delta <= tol, true, 'returns expected value' ); } t.end(); }); @@ -161,14 +161,14 @@ tape( 'the function returns `n` and stores `r` as two double-precision floating x = -1.0e-100 * randu(); n = rempio2( x, y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); z = ( PIO2*n ) + ( y[0]+y[1] ); delta = abs( z - x ); tol = EPS * abs( x ); - t.strictEqual( delta <= tol, true, 'can recover input value' ); + t.strictEqual( delta <= tol, true, 'returns expected value' ); } t.end(); }); @@ -187,14 +187,14 @@ tape( 'the function returns `n` and stores `r` as two double-precision floating for ( i = 0; i < x.length; i++ ) { n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); z = ( PIO2*n ) + ( y[0]+y[1] ); delta = abs( z - x[i] ); tol = EPS * abs( x[i] ); - t.strictEqual( delta <= tol, true, 'can recover input value' ); + t.strictEqual( delta <= tol, true, 'returns expected value' ); } t.end(); }); @@ -209,10 +209,10 @@ tape( 'for large positive input values, the function returns the last three bina y = [ 0.0, 0.0 ]; for ( i = 0; i < x.length; i++ ) { n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); + t.strictEqual( n < 8, true, 'returns expected value' ); } t.end(); }); @@ -227,10 +227,10 @@ tape( 'for large negative input values, the function returns the last three bina y = [ 0.0, 0.0 ]; for ( i = 0; i < x.length; i++ ) { n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); + t.strictEqual( n < 8, true, 'returns expected value' ); } t.end(); }); @@ -245,10 +245,10 @@ tape( 'for huge positive input values, the function returns the last three binar y = [ 0.0, 0.0 ]; for ( i = 0; i < x.length; i++ ) { n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); + t.strictEqual( n < 8, true, 'returns expected value' ); } t.end(); }); @@ -263,10 +263,10 @@ tape( 'for huge negative input values, the function returns the last three binar y = [ 0.0, 0.0 ]; for ( i = 0; i < x.length; i++ ) { n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); + t.strictEqual( n < 8, true, 'returns expected value' ); } t.end(); }); diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js b/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js index 42592a240311..3835d343de47 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/test/test.native.js @@ -35,6 +35,7 @@ var NINF = require( '@stdlib/constants/float64/ninf' ); var EPS = require( '@stdlib/constants/float64/eps' ); var PIO2 = require( '@stdlib/constants/float64/half-pi' ); var PIO4 = require( '@stdlib/constants/float64/fourth-pi' ); +var Float64Array = require( '@stdlib/array/float64' ); var tryRequire = require( '@stdlib/utils/try-require' ); @@ -55,26 +56,26 @@ tape( 'main export is a function', opts, function test( t ) { }); tape( 'the function returns `0` and sets the elements of `y` to `NaN` if provided `NaN`', opts, function test( t ) { - var y = [ 0.0, 0.0 ]; + var y = new Float64Array( [ 0.0, 0.0 ] ); var n = rempio2( NaN, y ); - t.strictEqual( n, 0.0, 'returns 0' ); - t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); - t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); + t.strictEqual( n, 0.0, 'returns expected value' ); + t.strictEqual( isnan( y[0] ), true, 'returns expected value' ); + t.strictEqual( isnan( y[1] ), true, 'returns expected value' ); t.end(); }); tape( 'the function returns `0` and sets the elements of `y` to `NaN` if provided positive or negative infinity', opts, function test( t ) { - var y = [ 0.0, 0.0 ]; + var y = new Float64Array( [ 0.0, 0.0 ] ); var n = rempio2( PINF, y ); - t.strictEqual( n, 0.0, 'returns 0' ); - t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); - t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); + t.strictEqual( n, 0.0, 'returns expected value' ); + t.strictEqual( isnan( y[0] ), true, 'returns expected value' ); + t.strictEqual( isnan( y[1] ), true, 'returns expected value' ); - y = [ 0.0, 0.0 ]; + y = new Float64Array( [ 0.0, 0.0 ] ); n = rempio2( NINF, y ); - t.strictEqual( n, 0.0, 'returns 0' ); - t.strictEqual( isnan( y[0] ), true, 'sets y[0] to NaN' ); - t.strictEqual( isnan( y[1] ), true, 'sets y[1] to NaN' ); + t.strictEqual( n, 0.0, 'returns expected value' ); + t.strictEqual( isnan( y[0] ), true, 'returns expected value' ); + t.strictEqual( isnan( y[1] ), true, 'returns expected value' ); t.end(); }); @@ -87,14 +88,14 @@ tape( 'the function returns `n` and stores `r` as two double-precision floating var n; var i; - y = [ 0.0, 0.0 ]; + y = new Float64Array( [ 0.0, 0.0 ] ); for ( i = 0; i < 200; i++ ) { x = randu()*100.0; n = rempio2( x, y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); z = ( PIO2*n ) + ( y[0]+y[1] ); delta = abs( z - x ); @@ -113,14 +114,14 @@ tape( 'the function returns `n` and stores `r` as two double-precision floating var n; var i; - y = [ 0.0, 0.0 ]; + y = new Float64Array( [ 0.0, 0.0 ] ); for ( i = 0; i < 200; i++ ) { x = randu()*1.0e-100; n = rempio2( x, y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); z = ( PIO2*n ) + ( y[0]+y[1] ); delta = abs( z - x ); @@ -139,14 +140,14 @@ tape( 'the function returns `n` and stores `r` as two double-precision floating var n; var i; - y = [ 0.0, 0.0 ]; + y = new Float64Array( [ 0.0, 0.0 ] ); for ( i = 0; i < 200; i++ ) { x = -100.0 * randu(); n = rempio2( x, y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); z = ( PIO2*n ) + ( y[0]+y[1] ); delta = abs( z - x ); @@ -165,14 +166,14 @@ tape( 'the function returns `n` and stores `r` as two double-precision floating var n; var i; - y = [ 0.0, 0.0 ]; + y = new Float64Array( [ 0.0, 0.0 ] ); for ( i = 0; i < 200; i++ ) { x = -1.0e-100 * randu(); n = rempio2( x, y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); z = ( PIO2*n ) + ( y[0]+y[1] ); delta = abs( z - x ); @@ -192,13 +193,13 @@ tape( 'the function returns `n` and stores `r` as two double-precision floating var i; x = incrspace( -10.0*PIO4, 10.0*PIO4, PIO4 ); - y = [ 0.0, 0.0 ]; + y = new Float64Array( [ 0.0, 0.0 ] ); for ( i = 0; i < x.length; i++ ) { n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); z = ( PIO2*n ) + ( y[0]+y[1] ); delta = abs( z - x[i] ); @@ -215,13 +216,13 @@ tape( 'for large positive input values, the function returns the last three bina var i; x = linspace( pow( 2.0, 20 )*PIO2, pow( 2.0, 60 )*PIO2, 4000 ); - y = [ 0.0, 0.0 ]; + y = new Float64Array( [ 0.0, 0.0 ] ); for ( i = 0; i < x.length; i++ ) { n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); + t.strictEqual( n < 8, true, 'returns expected value' ); } t.end(); }); @@ -233,13 +234,13 @@ tape( 'for large negative input values, the function returns the last three bina var i; x = linspace( -pow( 2.0, 20 )*PIO2, -pow( 2.0, 60 )*PIO2, 4000 ); - y = [ 0.0, 0.0 ]; + y = new Float64Array( [ 0.0, 0.0 ] ); for ( i = 0; i < x.length; i++ ) { n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); + t.strictEqual( n < 8, true, 'returns expected value' ); } t.end(); }); @@ -251,13 +252,13 @@ tape( 'for huge positive input values, the function returns the last three binar var i; x = linspace( pow( 2.0, 60 )*PIO2, pow( 2.0, 1020 )*PIO2, 6000 ); - y = [ 0.0, 0.0 ]; + y = new Float64Array( [ 0.0, 0.0 ] ); for ( i = 0; i < x.length; i++ ) { n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); + t.strictEqual( n < 8, true, 'returns expected value' ); } t.end(); }); @@ -269,13 +270,13 @@ tape( 'for huge negative input values, the function returns the last three binar var i; x = linspace( -pow( 2.0, 60 )*PIO2, -pow( 2.0, 1020 )*PIO2, 6000 ); - y = [ 0.0, 0.0 ]; + y = new Float64Array( [ 0.0, 0.0 ] ); for ( i = 0; i < x.length; i++ ) { n = rempio2( x[i], y ); - t.strictEqual( isInteger( n ), true, 'returns integer' ); - t.strictEqual( isNumber( y[0] ), true, 'assigns a number value to y[0]' ); - t.strictEqual( isNumber( y[1] ), true, 'assigns a number value to y[1]' ); - t.strictEqual( n < 8, true, 'returns a value smaller than 8' ); + t.strictEqual( isInteger( n ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[0] ), true, 'returns expected value' ); + t.strictEqual( isNumber( y[1] ), true, 'returns expected value' ); + t.strictEqual( n < 8, true, 'returns expected value' ); } t.end(); }); From 41993572c24c87399959bbdd3917fe0c0ed80bc4 Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Sat, 8 Jun 2024 19:34:41 +0530 Subject: [PATCH 21/28] use the same exponent mask as that in javascript implementation Signed-off-by: GUNJ JOSHI --- .../@stdlib/math/base/special/rempio2/src/main.c | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c index 3b938164f352..040157e83693 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c @@ -80,6 +80,9 @@ static const int32_t NINE_PIO4_HIGH_WORD = 0x401c463b; // 2^20*π/2 = 1647099.3291652855 => 0100000100111001001000011111101101010100010001000010110100011000 => high word => 0x413921fb = 1094263291 => 01000001001110010010000111111011 static const int32_t MEDIUM = 0x413921fb; +// Exponent mask (2047 => 0x7ff): +static const int32_t EXPONENT_MASK = 0x7ff; + /* * Table of constants for `2/π` (`396` hex digits, `476` decimal). @@ -473,7 +476,7 @@ int32_t stdlib_base_rempio2( const double x, double *rem1, double *rem2 ) { *rem1 = r - w; stdlib_base_float64_get_high_word( *rem1, &uhigh ); high = ( int32_t )uhigh; - i = j - ( ( high >> 20 ) & STDLIB_CONSTANT_FLOAT64_HIGH_WORD_EXPONENT_MASK ); + i = j - ( ( high >> 20 ) & EXPONENT_MASK ); // Check if a second iteration is needed (good to 118 bits)... if ( i > 16 ) { From ff58f1ab87bcd40bd1258e7d5da1d450f338d654 Mon Sep 17 00:00:00 2001 From: Athan Date: Sat, 8 Jun 2024 12:04:20 -0700 Subject: [PATCH 22/28] Apply suggestions from code review Signed-off-by: Athan --- lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c index 040157e83693..37e69def5a24 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c @@ -487,7 +487,7 @@ int32_t stdlib_base_rempio2( const double x, double *rem1, double *rem2 ) { *rem1 = r - w; stdlib_base_float64_get_high_word( *rem1, &uhigh ); high = ( int32_t )uhigh; - i = j - ( ( high >> 20 ) & STDLIB_CONSTANT_FLOAT64_HIGH_WORD_EXPONENT_MASK ); + i = j - ( ( high >> 20 ) & EXPONENT_MASK ); // Check if a third iteration is needed (151 bits accumulated)... if ( i > 49 ) { From d0c2b0646f787da2cf9c12e88342676eeb885699 Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Sun, 9 Jun 2024 08:20:49 +0530 Subject: [PATCH 23/28] declare arrays inside corresponding functions Signed-off-by: GUNJ JOSHI --- .../@stdlib/math/base/special/rempio2/src/main.c | 14 ++++++-------- 1 file changed, 6 insertions(+), 8 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c index 37e69def5a24..c76f5ef275e9 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c @@ -121,14 +121,6 @@ static const double PIO2[] = { 2.16741683877804819444e-51 // 0x3569F31D, 0x00000000 }; -// Arrays for storing temporary values: -double TX[ 3 ] = { 0.0 }; -double TY[ 2 ] = { 0.0 }; -double F[ 20 ]; -double Q[ 20 ]; -double FQ[ 20 ]; -int32_t IQ[ 20 ]; - /** * Returns the last three binary digits of `N` with `y = x - Nπ/2` so that `|y| < π/2`. * @@ -145,6 +137,10 @@ int32_t IQ[ 20 ]; * @return last 3 binary digits */ int32_t rempio2Kernel( const double* x, double *y, int32_t e0, int32_t nx ) { + int32_t IQ[ 20 ]; + double FQ[ 20 ]; + double F[ 20 ]; + double Q[ 20 ]; int32_t carry; int32_t jz; int32_t jx; @@ -364,6 +360,8 @@ int32_t rempio2Kernel( const double* x, double *y, int32_t e0, int32_t nx ) { * int32_t n = stdlib_base_rempio2( 128.0, &rem1, &rem2 ); */ int32_t stdlib_base_rempio2( const double x, double *rem1, double *rem2 ) { + double TX[ 3 ]; + double TY[ 2 ]; uint32_t uhigh; uint32_t ulow; uint32_t uhx; From 4f816305964693e53712d7cfac09fc91ace09a6e Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Sun, 9 Jun 2024 10:41:14 +0530 Subject: [PATCH 24/28] add extended comment for kernel rempio2 Signed-off-by: GUNJ JOSHI --- .../special/rempio2/lib/kernel_rempio2.js | 100 +++++++++++++++++- 1 file changed, 99 insertions(+), 1 deletion(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/lib/kernel_rempio2.js b/lib/node_modules/@stdlib/math/base/special/rempio2/lib/kernel_rempio2.js index e6aebb0912bf..b3e3d0c49c38 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/lib/kernel_rempio2.js +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/lib/kernel_rempio2.js @@ -272,7 +272,105 @@ function compute( x, y, jz, q, q0, jk, jv, jx, f ) { * * ## Method * -* - The method is to compute the integer (`mod 8`) and fraction parts of `2x/π` without doing the full multiplication. In general, we skip the part of the product that is known to be a huge integer (more accurately, equals `0 mod 8` ). Thus, the number of operations is independent of the exponent of the input. +* - The method is to compute the integer (mod 8) and fraction parts of (2/pi)\*x without doing the full multiplication. In general we skip the part of the product that are known to be a huge integer (more accurately, = 0 mod 8 ). Thus the number of operations are independent of the exponent of the input. +* +* - (2/pi) is represented by an array of 24-bit integers in ipio2\[\]. +* +* - Input parameters: +* +* x\[\] The input value (must be positive) is broken into nx pieces of 24-bit integers in double precision format. x\[i\] will be the i-th 24 bit of x. The scaled exponent of x\[0\] is given in input parameter e0 (i.e., x\[0\]\*2^e0 match x's up to 24 bits. +* +* Example of breaking a double positive z into x\[0\]+x\[1\]+x\[2\]: +* +* ```tex +* e0 = \mathrm{ilogb}(z) - 23 +* z = \mathrm{scalbn}(z, -e0) +* ``` +* +* for i = 0,1,2 +* +* ```tex +* x[i] = \lfloor z \rfloor +* z = (z - x[i]) \times 2^{24} +* ``` +* +* y\[\] output result in an array of double precision numbers. +* +* The dimension of y\[\] is: +* 24-bit precision 1 +* 53-bit precision 2 +* 64-bit precision 2 +* 113-bit precision 3 +* +* The actual value is the sum of them. Thus for 113-bit precision, one may have to do something like: +* +* ```tex +* \mathrm{long\ double} \: t, w, r_{\text{head}}, r_{\text{tail}}; \\ +* t &= (\mathrm{long\ double}) y[2] + (\mathrm{long\ double}) y[1]; \\ +* w &= (\mathrm{long\ double}) y[0]; \\ +* r_{\text{head}} &= t + w; \\ +* r_{\text{tail}} &= w - (r_{\text{head}} - t); +* ``` +* +* e0 The exponent of x\[0\]. Must be <= 16360 or you need to expand the ipio2 table. +* +* nx dimension of x\[\] +* +* prec an integer indicating the precision: +* 0 24 bits (single) +* 1 53 bits (double) +* 2 64 bits (extended) +* 3 113 bits (quad) +* +* - External function: +* +* double scalbn(), floor(); +* +* - Here is the description of some local variables: +* +* jk jk+1 is the initial number of terms of ipio2\[\] needed in the computation. The minimum and recommended value for jk is 3,4,4,6 for single, double, extended, and quad. jk+1 must be 2 larger than you might expect so that our recomputation test works. (Up to 24 bits in the integer part (the 24 bits of it that we compute) and 23 bits in the fraction part may be lost to cancellation before we recompute.) +* +* jz local integer variable indicating the number of terms of ipio2\[\] used. +* +* jx nx - 1 +* +* jv index for pointing to the suitable ipio2\[\] for the computation. In general, we want +* +* ```tex +* \frac{{2^{e0} \cdot x[0] \cdot \mathrm{ipio2}[jv-1] \cdot 2^{-24jv}}}{{8}} +* ``` +* +* is an integer. Thus +* +* ```tex +* e0 - 3 - 24 \cdot jv \geq 0 \quad \text{or} \quad \frac{{e0 - 3}}{{24}} \geq jv +* ``` +* +* Hence +* +* ```tex +* jv = \max(0, \frac{{e0 - 3}}{{24}}) +* ``` +* +* jp jp+1 is the number of terms in PIo2\[\] needed, jp = jk. +* +* q\[\] double array with integral value, representing the 24-bits chunk of the product of x and 2/pi. +* +* q0 the corresponding exponent of q\[0\]. Note that the exponent for q\[i\] would be q0-24\*i. +* +* PIo2\[\] double precision array, obtained by cutting pi/2 into 24 bits chunks. +* +* f\[\] ipso2\[\] in floating point +* +* iq\[\] integer array by breaking up q\[\] in 24-bits chunk. +* +* fq\[\] final product of x\*(2/pi) in fq\[0\],..,fq\[jk\] +* +* ih integer. If >0 it indicates q\[\] is >= 0.5, hence it also indicates the \*sign\* of the result. +* +* - Constants: +* +* The hexadecimal values are the intended ones for the following constants. The decimal values may be used, provided that the compiler will convert from decimal to binary accurately enough to produce the hexadecimal values shown. * * @private * @param {PositiveNumber} x - input value From fd93a9de98d2623077fa99f7f05dd808f40daba0 Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Sun, 9 Jun 2024 10:43:02 +0530 Subject: [PATCH 25/28] add extended comment for kernel rempio2 Signed-off-by: GUNJ JOSHI --- .../math/base/special/rempio2/src/main.c | 100 +++++++++++++++++- 1 file changed, 99 insertions(+), 1 deletion(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c index c76f5ef275e9..aeeacd3268ab 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c @@ -126,7 +126,105 @@ static const double PIO2[] = { * * ## Method * -* - The method is to compute the integer (`mod 8`) and fraction parts of `2x/π` without doing the full multiplication. In general, we skip the part of the product that is known to be a huge integer (more accurately, equals `0 mod 8` ). Thus, the number of operations is independent of the exponent of the input. +* - The method is to compute the integer (mod 8) and fraction parts of (2/pi)\*x without doing the full multiplication. In general we skip the part of the product that are known to be a huge integer (more accurately, = 0 mod 8 ). Thus the number of operations are independent of the exponent of the input. +* +* - (2/pi) is represented by an array of 24-bit integers in ipio2\[\]. +* +* - Input parameters: +* +* x\[\] The input value (must be positive) is broken into nx pieces of 24-bit integers in double precision format. x\[i\] will be the i-th 24 bit of x. The scaled exponent of x\[0\] is given in input parameter e0 (i.e., x\[0\]\*2^e0 match x's up to 24 bits. +* +* Example of breaking a double positive z into x\[0\]+x\[1\]+x\[2\]: +* +* ```tex +* e0 = \mathrm{ilogb}(z) - 23 +* z = \mathrm{scalbn}(z, -e0) +* ``` +* +* for i = 0,1,2 +* +* ```tex +* x[i] = \lfloor z \rfloor +* z = (z - x[i]) \times 2^{24} +* ``` +* +* y\[\] output result in an array of double precision numbers. +* +* The dimension of y\[\] is: +* 24-bit precision 1 +* 53-bit precision 2 +* 64-bit precision 2 +* 113-bit precision 3 +* +* The actual value is the sum of them. Thus for 113-bit precision, one may have to do something like: +* +* ```tex +* \mathrm{long\ double} \: t, w, r_{\text{head}}, r_{\text{tail}}; \\ +* t &= (\mathrm{long\ double}) y[2] + (\mathrm{long\ double}) y[1]; \\ +* w &= (\mathrm{long\ double}) y[0]; \\ +* r_{\text{head}} &= t + w; \\ +* r_{\text{tail}} &= w - (r_{\text{head}} - t); +* ``` +* +* e0 The exponent of x\[0\]. Must be <= 16360 or you need to expand the ipio2 table. +* +* nx dimension of x\[\] +* +* prec an integer indicating the precision: +* 0 24 bits (single) +* 1 53 bits (double) +* 2 64 bits (extended) +* 3 113 bits (quad) +* +* - External function: +* +* double scalbn(), floor(); +* +* - Here is the description of some local variables: +* +* jk jk+1 is the initial number of terms of ipio2\[\] needed in the computation. The minimum and recommended value for jk is 3,4,4,6 for single, double, extended, and quad. jk+1 must be 2 larger than you might expect so that our recomputation test works. (Up to 24 bits in the integer part (the 24 bits of it that we compute) and 23 bits in the fraction part may be lost to cancellation before we recompute.) +* +* jz local integer variable indicating the number of terms of ipio2\[\] used. +* +* jx nx - 1 +* +* jv index for pointing to the suitable ipio2\[\] for the computation. In general, we want +* +* ```tex +* \frac{{2^{e0} \cdot x[0] \cdot \mathrm{ipio2}[jv-1] \cdot 2^{-24jv}}}{{8}} +* ``` +* +* is an integer. Thus +* +* ```tex +* e0 - 3 - 24 \cdot jv \geq 0 \quad \text{or} \quad \frac{{e0 - 3}}{{24}} \geq jv +* ``` +* +* Hence +* +* ```tex +* jv = \max(0, \frac{{e0 - 3}}{{24}}) +* ``` +* +* jp jp+1 is the number of terms in PIo2\[\] needed, jp = jk. +* +* q\[\] double array with integral value, representing the 24-bits chunk of the product of x and 2/pi. +* +* q0 the corresponding exponent of q\[0\]. Note that the exponent for q\[i\] would be q0-24\*i. +* +* PIo2\[\] double precision array, obtained by cutting pi/2 into 24 bits chunks. +* +* f\[\] ipso2\[\] in floating point +* +* iq\[\] integer array by breaking up q\[\] in 24-bits chunk. +* +* fq\[\] final product of x\*(2/pi) in fq\[0\],..,fq\[jk\] +* +* ih integer. If >0 it indicates q\[\] is >= 0.5, hence it also indicates the \*sign\* of the result. +* +* - Constants: +* +* The hexadecimal values are the intended ones for the following constants. The decimal values may be used, provided that the compiler will convert from decimal to binary accurately enough to produce the hexadecimal values shown. * * @private * @param x input value From c7cad3e883fa4115fb686adcb2c9575afef43b6d Mon Sep 17 00:00:00 2001 From: GUNJ JOSHI Date: Sun, 9 Jun 2024 11:59:50 +0530 Subject: [PATCH 26/28] initialize arrays F and Q at the time of initialization Signed-off-by: GUNJ JOSHI --- lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c index aeeacd3268ab..264bd78d0900 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c @@ -235,10 +235,10 @@ static const double PIO2[] = { * @return last 3 binary digits */ int32_t rempio2Kernel( const double* x, double *y, int32_t e0, int32_t nx ) { + double F[ 20 ] = { 0.0 }; + double Q[ 20 ] = { 0.0 }; int32_t IQ[ 20 ]; double FQ[ 20 ]; - double F[ 20 ]; - double Q[ 20 ]; int32_t carry; int32_t jz; int32_t jx; From 276dd5158c6c31b0e1c9788dd3a0ec3c66a7ccfa Mon Sep 17 00:00:00 2001 From: Philipp Burckhardt Date: Sun, 9 Jun 2024 18:46:45 -0400 Subject: [PATCH 27/28] refactor: avoid duplicated check for zero --- .../@stdlib/math/base/special/rempio2/src/main.c | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c index 264bd78d0900..a1610166ee10 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c @@ -376,10 +376,7 @@ int32_t rempio2Kernel( const double* x, double *y, int32_t e0, int32_t nx ) { jz += k; goto recompute; } - } - - // Chop off zero terms... - if ( z == 0.0 ) { + // Chop off zero terms... jz -= 1; q0 -= 24; while ( IQ[ jz ] == 0 ) { From ea4c65c9ab56cee2c717c3a10400ba6612a331d7 Mon Sep 17 00:00:00 2001 From: Philipp Burckhardt Date: Sun, 9 Jun 2024 19:22:21 -0400 Subject: [PATCH 28/28] docs: use bullets and update comment --- .../special/rempio2/lib/kernel_rempio2.js | 122 +++++++++--------- .../math/base/special/rempio2/src/main.c | 118 ++++++++--------- 2 files changed, 119 insertions(+), 121 deletions(-) diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/lib/kernel_rempio2.js b/lib/node_modules/@stdlib/math/base/special/rempio2/lib/kernel_rempio2.js index b3e3d0c49c38..6417f94af16d 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/lib/kernel_rempio2.js +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/lib/kernel_rempio2.js @@ -206,9 +206,7 @@ function compute( x, y, jz, q, q0, jk, jv, jx, f ) { jz += k; return compute( x, y, jz, q, q0, jk, jv, jx, f ); } - } - // Chop off zero terms... - if ( z === 0.0 ) { + // Chop off zero terms... jz -= 1; q0 -= 24; while ( IQ[ jz ] === 0 ) { @@ -272,105 +270,105 @@ function compute( x, y, jz, q, q0, jk, jv, jx, f ) { * * ## Method * -* - The method is to compute the integer (mod 8) and fraction parts of (2/pi)\*x without doing the full multiplication. In general we skip the part of the product that are known to be a huge integer (more accurately, = 0 mod 8 ). Thus the number of operations are independent of the exponent of the input. +* - The method is to compute the integer (mod 8) and fraction parts of (2/π) * x without doing the full multiplication. In general, we skip the part of the product that are known to be a huge integer (more accurately, = 0 mod 8 ). Thus the number of operations are independent of the exponent of the input. * -* - (2/pi) is represented by an array of 24-bit integers in ipio2\[\]. +* - (2/π) is represented by an array of 24-bit integers in `ipio2[]`. * * - Input parameters: * -* x\[\] The input value (must be positive) is broken into nx pieces of 24-bit integers in double precision format. x\[i\] will be the i-th 24 bit of x. The scaled exponent of x\[0\] is given in input parameter e0 (i.e., x\[0\]\*2^e0 match x's up to 24 bits. +* - `x[]` The input value (must be positive) is broken into `nx` pieces of 24-bit integers in double precision format. `x[i]` will be the i-th 24 bit of x. The scaled exponent of `x[0]` is given in input parameter `e0` (i.e., `x[0]*2^e0` match x's up to 24 bits). * -* Example of breaking a double positive z into x\[0\]+x\[1\]+x\[2\]: +* Example of breaking a double positive `z` into `x[0]+x[1]+x[2]`: * -* ```tex -* e0 = \mathrm{ilogb}(z) - 23 -* z = \mathrm{scalbn}(z, -e0) -* ``` +* ```tex +* e0 = \mathrm{ilogb}(z) - 23 +* z = \mathrm{scalbn}(z, -e0) +* ``` * -* for i = 0,1,2 +* for `i = 0,1,2` * -* ```tex -* x[i] = \lfloor z \rfloor -* z = (z - x[i]) \times 2^{24} -* ``` +* ```tex +* x[i] = \lfloor z \rfloor +* z = (z - x[i]) \times 2^{24} +* ``` * -* y\[\] output result in an array of double precision numbers. +* - `y[]` output result in an array of double precision numbers. * -* The dimension of y\[\] is: -* 24-bit precision 1 -* 53-bit precision 2 -* 64-bit precision 2 -* 113-bit precision 3 +* The dimension of `y[]` is: +* 24-bit precision 1 +* 53-bit precision 2 +* 64-bit precision 2 +* 113-bit precision 3 * -* The actual value is the sum of them. Thus for 113-bit precision, one may have to do something like: +* The actual value is the sum of them. Thus, for 113-bit precision, one may have to do something like: * -* ```tex -* \mathrm{long\ double} \: t, w, r_{\text{head}}, r_{\text{tail}}; \\ -* t &= (\mathrm{long\ double}) y[2] + (\mathrm{long\ double}) y[1]; \\ -* w &= (\mathrm{long\ double}) y[0]; \\ -* r_{\text{head}} &= t + w; \\ -* r_{\text{tail}} &= w - (r_{\text{head}} - t); -* ``` +* ```tex +* \mathrm{long\ double} \: t, w, r_{\text{head}}, r_{\text{tail}}; \\ +* t &= (\mathrm{long\ double}) y[2] + (\mathrm{long\ double}) y[1]; \\ +* w &= (\mathrm{long\ double}) y[0]; \\ +* r_{\text{head}} &= t + w; \\ +* r_{\text{tail}} &= w - (r_{\text{head}} - t); +* ``` * -* e0 The exponent of x\[0\]. Must be <= 16360 or you need to expand the ipio2 table. +* - `e0` The exponent of `x[0]`. Must be <= 16360 or you need to expand the `ipio2` table. * -* nx dimension of x\[\] +* - `nx` dimension of `x[]` * -* prec an integer indicating the precision: -* 0 24 bits (single) -* 1 53 bits (double) -* 2 64 bits (extended) -* 3 113 bits (quad) +* - `prec` an integer indicating the precision: +* 0 24 bits (single) +* 1 53 bits (double) +* 2 64 bits (extended) +* 3 113 bits (quad) * * - External function: * -* double scalbn(), floor(); +* - double `scalbn()`, `floor()`; * * - Here is the description of some local variables: * -* jk jk+1 is the initial number of terms of ipio2\[\] needed in the computation. The minimum and recommended value for jk is 3,4,4,6 for single, double, extended, and quad. jk+1 must be 2 larger than you might expect so that our recomputation test works. (Up to 24 bits in the integer part (the 24 bits of it that we compute) and 23 bits in the fraction part may be lost to cancellation before we recompute.) +* - `jk` `jk+1` is the initial number of terms of `ipio2[]` needed in the computation. The minimum and recommended value for `jk` is 3,4,4,6 for single, double, extended, and quad. `jk+1` must be 2 larger than you might expect so that our recomputation test works. (Up to 24 bits in the integer part (the 24 bits of it that we compute) and 23 bits in the fraction part may be lost to cancellation before we recompute.) * -* jz local integer variable indicating the number of terms of ipio2\[\] used. +* - `jz` local integer variable indicating the number of terms of `ipio2[]` used. * -* jx nx - 1 +* - `jx` `nx - 1` * -* jv index for pointing to the suitable ipio2\[\] for the computation. In general, we want +* - `jv` index for pointing to the suitable `ipio2[]` for the computation. In general, we want * -* ```tex -* \frac{{2^{e0} \cdot x[0] \cdot \mathrm{ipio2}[jv-1] \cdot 2^{-24jv}}}{{8}} -* ``` +* ```tex +* \frac{{2^{e0} \cdot x[0] \cdot \mathrm{ipio2}[jv-1] \cdot 2^{-24jv}}}{{8}} +* ``` * -* is an integer. Thus +* to be an integer. Thus * -* ```tex -* e0 - 3 - 24 \cdot jv \geq 0 \quad \text{or} \quad \frac{{e0 - 3}}{{24}} \geq jv -* ``` +* ```tex +* e0 - 3 - 24 \cdot jv \geq 0 \quad \text{or} \quad \frac{{e0 - 3}}{{24}} \geq jv +* ``` * -* Hence +* Hence * -* ```tex -* jv = \max(0, \frac{{e0 - 3}}{{24}}) -* ``` +* ```tex +* jv = \max(0, \frac{{e0 - 3}}{{24}}) +* ``` * -* jp jp+1 is the number of terms in PIo2\[\] needed, jp = jk. +* - `jp` `jp+1` is the number of terms in `PIo2[]` needed, `jp = jk`. * -* q\[\] double array with integral value, representing the 24-bits chunk of the product of x and 2/pi. +* - `q[]` double array with integral value, representing the 24-bits chunk of the product of `x` and `2/π`. * -* q0 the corresponding exponent of q\[0\]. Note that the exponent for q\[i\] would be q0-24\*i. +* - `q0` the corresponding exponent of `q[0]`. Note that the exponent for `q[i]` would be `q0-24*i`. * -* PIo2\[\] double precision array, obtained by cutting pi/2 into 24 bits chunks. +* - `PIo2[]` double precision array, obtained by cutting `π/2` into 24 bits chunks. * -* f\[\] ipso2\[\] in floating point +* - `f[]` `ipso2[]` in floating point * -* iq\[\] integer array by breaking up q\[\] in 24-bits chunk. +* - `iq[]` integer array by breaking up `q[]` in 24-bits chunk. * -* fq\[\] final product of x\*(2/pi) in fq\[0\],..,fq\[jk\] +* - `fq[]` final product of `x*(2/π)` in `fq[0],..,fq[jk]` * -* ih integer. If >0 it indicates q\[\] is >= 0.5, hence it also indicates the \*sign\* of the result. +* - `ih` integer. If >0 it indicates `q[]` is >= 0.5, hence it also indicates the _sign_ of the result. * * - Constants: * -* The hexadecimal values are the intended ones for the following constants. The decimal values may be used, provided that the compiler will convert from decimal to binary accurately enough to produce the hexadecimal values shown. +* - The hexadecimal values are the intended ones for the following constants. The decimal values may be used, provided that the compiler will convert from decimal to binary accurately enough to produce the hexadecimal values shown. * * @private * @param {PositiveNumber} x - input value diff --git a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c index a1610166ee10..dc7289455ef1 100644 --- a/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c +++ b/lib/node_modules/@stdlib/math/base/special/rempio2/src/main.c @@ -126,105 +126,105 @@ static const double PIO2[] = { * * ## Method * -* - The method is to compute the integer (mod 8) and fraction parts of (2/pi)\*x without doing the full multiplication. In general we skip the part of the product that are known to be a huge integer (more accurately, = 0 mod 8 ). Thus the number of operations are independent of the exponent of the input. +* - The method is to compute the integer (mod 8) and fraction parts of (2/π) * x without doing the full multiplication. In general, we skip the part of the product that are known to be a huge integer (more accurately, = 0 mod 8 ). Thus the number of operations are independent of the exponent of the input. * -* - (2/pi) is represented by an array of 24-bit integers in ipio2\[\]. +* - (2/π) is represented by an array of 24-bit integers in `ipio2[]`. * * - Input parameters: * -* x\[\] The input value (must be positive) is broken into nx pieces of 24-bit integers in double precision format. x\[i\] will be the i-th 24 bit of x. The scaled exponent of x\[0\] is given in input parameter e0 (i.e., x\[0\]\*2^e0 match x's up to 24 bits. +* - `x[]` The input value (must be positive) is broken into `nx` pieces of 24-bit integers in double precision format. `x[i]` will be the i-th 24 bit of x. The scaled exponent of `x[0]` is given in input parameter `e0` (i.e., `x[0]*2^e0` match x's up to 24 bits). * -* Example of breaking a double positive z into x\[0\]+x\[1\]+x\[2\]: +* Example of breaking a double positive `z` into `x[0]+x[1]+x[2]`: * -* ```tex -* e0 = \mathrm{ilogb}(z) - 23 -* z = \mathrm{scalbn}(z, -e0) -* ``` +* ```tex +* e0 = \mathrm{ilogb}(z) - 23 +* z = \mathrm{scalbn}(z, -e0) +* ``` * -* for i = 0,1,2 +* for `i = 0,1,2` * -* ```tex -* x[i] = \lfloor z \rfloor -* z = (z - x[i]) \times 2^{24} -* ``` +* ```tex +* x[i] = \lfloor z \rfloor +* z = (z - x[i]) \times 2^{24} +* ``` * -* y\[\] output result in an array of double precision numbers. +* - `y[]` output result in an array of double precision numbers. * -* The dimension of y\[\] is: -* 24-bit precision 1 -* 53-bit precision 2 -* 64-bit precision 2 -* 113-bit precision 3 +* The dimension of `y[]` is: +* 24-bit precision 1 +* 53-bit precision 2 +* 64-bit precision 2 +* 113-bit precision 3 * -* The actual value is the sum of them. Thus for 113-bit precision, one may have to do something like: +* The actual value is the sum of them. Thus, for 113-bit precision, one may have to do something like: * -* ```tex -* \mathrm{long\ double} \: t, w, r_{\text{head}}, r_{\text{tail}}; \\ -* t &= (\mathrm{long\ double}) y[2] + (\mathrm{long\ double}) y[1]; \\ -* w &= (\mathrm{long\ double}) y[0]; \\ -* r_{\text{head}} &= t + w; \\ -* r_{\text{tail}} &= w - (r_{\text{head}} - t); -* ``` +* ```tex +* \mathrm{long\ double} \: t, w, r_{\text{head}}, r_{\text{tail}}; \\ +* t &= (\mathrm{long\ double}) y[2] + (\mathrm{long\ double}) y[1]; \\ +* w &= (\mathrm{long\ double}) y[0]; \\ +* r_{\text{head}} &= t + w; \\ +* r_{\text{tail}} &= w - (r_{\text{head}} - t); +* ``` * -* e0 The exponent of x\[0\]. Must be <= 16360 or you need to expand the ipio2 table. +* - `e0` The exponent of `x[0]`. Must be <= 16360 or you need to expand the `ipio2` table. * -* nx dimension of x\[\] +* - `nx` dimension of `x[]` * -* prec an integer indicating the precision: -* 0 24 bits (single) -* 1 53 bits (double) -* 2 64 bits (extended) -* 3 113 bits (quad) +* - `prec` an integer indicating the precision: +* 0 24 bits (single) +* 1 53 bits (double) +* 2 64 bits (extended) +* 3 113 bits (quad) * * - External function: * -* double scalbn(), floor(); +* - double `scalbn()`, `floor()`; * * - Here is the description of some local variables: * -* jk jk+1 is the initial number of terms of ipio2\[\] needed in the computation. The minimum and recommended value for jk is 3,4,4,6 for single, double, extended, and quad. jk+1 must be 2 larger than you might expect so that our recomputation test works. (Up to 24 bits in the integer part (the 24 bits of it that we compute) and 23 bits in the fraction part may be lost to cancellation before we recompute.) +* - `jk` `jk+1` is the initial number of terms of `ipio2[]` needed in the computation. The minimum and recommended value for `jk` is 3,4,4,6 for single, double, extended, and quad. `jk+1` must be 2 larger than you might expect so that our recomputation test works. (Up to 24 bits in the integer part (the 24 bits of it that we compute) and 23 bits in the fraction part may be lost to cancellation before we recompute.) * -* jz local integer variable indicating the number of terms of ipio2\[\] used. +* - `jz` local integer variable indicating the number of terms of `ipio2[]` used. * -* jx nx - 1 +* - `jx` `nx - 1` * -* jv index for pointing to the suitable ipio2\[\] for the computation. In general, we want +* - `jv` index for pointing to the suitable `ipio2[]` for the computation. In general, we want * -* ```tex -* \frac{{2^{e0} \cdot x[0] \cdot \mathrm{ipio2}[jv-1] \cdot 2^{-24jv}}}{{8}} -* ``` +* ```tex +* \frac{{2^{e0} \cdot x[0] \cdot \mathrm{ipio2}[jv-1] \cdot 2^{-24jv}}}{{8}} +* ``` * -* is an integer. Thus +* to be an integer. Thus * -* ```tex -* e0 - 3 - 24 \cdot jv \geq 0 \quad \text{or} \quad \frac{{e0 - 3}}{{24}} \geq jv -* ``` +* ```tex +* e0 - 3 - 24 \cdot jv \geq 0 \quad \text{or} \quad \frac{{e0 - 3}}{{24}} \geq jv +* ``` * -* Hence +* Hence * -* ```tex -* jv = \max(0, \frac{{e0 - 3}}{{24}}) -* ``` +* ```tex +* jv = \max(0, \frac{{e0 - 3}}{{24}}) +* ``` * -* jp jp+1 is the number of terms in PIo2\[\] needed, jp = jk. +* - `jp` `jp+1` is the number of terms in `PIo2[]` needed, `jp = jk`. * -* q\[\] double array with integral value, representing the 24-bits chunk of the product of x and 2/pi. +* - `q[]` double array with integral value, representing the 24-bits chunk of the product of `x` and `2/π`. * -* q0 the corresponding exponent of q\[0\]. Note that the exponent for q\[i\] would be q0-24\*i. +* - `q0` the corresponding exponent of `q[0]`. Note that the exponent for `q[i]` would be `q0-24*i`. * -* PIo2\[\] double precision array, obtained by cutting pi/2 into 24 bits chunks. +* - `PIo2[]` double precision array, obtained by cutting `π/2` into 24 bits chunks. * -* f\[\] ipso2\[\] in floating point +* - `f[]` `ipso2[]` in floating point * -* iq\[\] integer array by breaking up q\[\] in 24-bits chunk. +* - `iq[]` integer array by breaking up `q[]` in 24-bits chunk. * -* fq\[\] final product of x\*(2/pi) in fq\[0\],..,fq\[jk\] +* - `fq[]` final product of `x*(2/π)` in `fq[0],..,fq[jk]` * -* ih integer. If >0 it indicates q\[\] is >= 0.5, hence it also indicates the \*sign\* of the result. +* - `ih` integer. If >0 it indicates `q[]` is >= 0.5, hence it also indicates the _sign_ of the result. * * - Constants: * -* The hexadecimal values are the intended ones for the following constants. The decimal values may be used, provided that the compiler will convert from decimal to binary accurately enough to produce the hexadecimal values shown. +* - The hexadecimal values are the intended ones for the following constants. The decimal values may be used, provided that the compiler will convert from decimal to binary accurately enough to produce the hexadecimal values shown. * * @private * @param x input value