diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/README.md b/lib/node_modules/@stdlib/lapack/base/dlasq2/README.md new file mode 100644 index 000000000000..8b3121a75639 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/README.md @@ -0,0 +1,263 @@ + + +# dlasq2 + +> Compute all eigenvalues of the symmetric positive definite tridiagonal matrix associated with the QD Array `Z` to high relative accuracy. + +
+ +The `dlasq2` routine computes all the eigenvalues of the symmetric positive definite tri-diagonal matrix associated with the QD (quotient-difference) array `Z` to high relative accuracy, in the absence of denormalization, underflow and overflow. + +To see the relation of Z to the tridiagonal matrix, let `L` be a unit lower bidiagonal matrix with subdiagonals `Z(2,4,6,,..)` and let `U` be an upper bidiagonal matrix with 1's above and diagonal `Z(1,3,5,,..)`. The tridiagonal is `L*U` or, if you prefer, the symmetric tridiagonal to which it is similar. + +The QD array `Z` of length `4*N` contains pairs `(q_i, qq_i, e_i, ee_i)` for `i = 1, 2, ..., N`, where `q_i` and `e_i` are derived from the diagonal and off-diagonal elements of `T` and `qq_i` and `ee_i` are working storage. + +Upon successful completion, the first `N` indexed elements contain the eigenvalues in descending order, and `Z( 2*N+1 )` holds the trace, and `Z( 2*N+2 )` holds the sum of the eigenvalues. + +If `N > 2`, then `Z( 2*N+3 )` holds the iteration count, `Z( 2*N+4 )` holds `NDIVS/NIN^2`, and `Z( 2*N+5 )` holds the percentage of shifts that failed. + +
+ + + +
+ +## Usage + +```javascript +var dlasq2 = require( '@stdlib/lapack/base/dlasq2' ); +``` + +#### dlasq2( N, Z ) + +Computes all eigenvalues of the symmetric positive definite tridiagonal matrix associated with the QD Array `Z` to high relative accuracy. + + + +```javascript +var Float64Array = require( '@stdlib/array/float64' ); + +var Z = new Float64Array( [ 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); + +var out = dlasq2( 4, Z ); +// Z => [ ~115.713, ~83.153, ~62.17, ~42.464, ~83.153, ~20.987, ~0.0, ~0.0, 303.5, ~303.5, 15, 4.0625, 0.0, ~0.0, ~0.172, ~0.096 ] +// out => 0 +``` + +The function has the following parameters: + +- **N**: number of rows/columns in the tridiagonal matrix. +- **Z**: the QD array as a [`Float64Array`][@stdlib/array/float64]. Must have at least `4*N` indexed elements. `Z` is overwritten during computation. Upon successful exit, the first `N` indexed elements of `Z` contain the eigenvalues in decreasing order. + +Note that indexing is relative to the first index. To introduce an offset, use [`typed array`][mdn-typed-array] views. + + + +```javascript +var Float64Array = require( '@stdlib/array/float64' ); + +// Initial arrays... +var Z0 = new Float64Array( [ 0.0, 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); + +// Create offset views... +var Z1 = new Float64Array( Z0.buffer, Z0.BYTES_PER_ELEMENT*1 ); // start at 2nd element + +dlasq2( 4, Z1 ); +// Z0 => [ 0.0, ~115.713, ~83.153, ~62.17, ~42.464, ~83.153, ~20.987, ~0.0, ~0.0, 303.5, ~303.5, 15, 4.0625, 0.0, ~0.0, ~0.172, ~0.096 ] +``` + + + +#### dlasq2.ndarray( N, Z, strideZ, offsetZ ) + +Computes all eigenvalues of the symmetric positive definite tridiagonal matrix associated with the QD Array `Z` to high relative accuracy using alternative indexing semantics. + +```javascript +var Float64Array = require( '@stdlib/array/float64' ); + +var Z = new Float64Array( [ 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); + +dlasq2.ndarray( 4, Z, 1, 0 ); +// Z => [ ~115.713, ~83.153, ~62.17, ~42.464, ~83.153, ~20.987, ~0.0, ~0.0, 303.5, ~303.5, 15, 4.0625, 0.0, ~0.0, ~0.172, ~0.096 ] +``` + +The function has the following additional parameters: + +- **strideZ**: stride length for `Z`. +- **offsetZ**: starting index for `Z`. + +While [`typed array`][mdn-typed-array] views mandate a view offset based on the underlying buffer, the offset parameters support indexing semantics based on starting indices. For example, + + + +```javascript +var Float64Array = require( '@stdlib/array/float64' ); + +var Z = new Float64Array( [ 0.0, 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); + +dlasq2.ndarray( 4, Z, 1, 1 ); +// Z => [ 0.0, ~115.713, ~83.153, ~62.17, ~42.464, ~83.153, ~20.987, ~0.0, ~0.0, 303.5, ~303.5, 15, 4.0625, 0.0, ~0.0, ~0.172, ~0.096 ] +``` + +
+ + + +
+ +## Notes + +- Both functions mutate the input array `Z`. + +- Both functions return a status code indicating success or failure. A status code indicates the following conditions: + + - If `=0`, then the algorithm completed successfully. + - If `<0`: + + - If the i-th argument is a scalar and had an illegal value, then `INFO = -i`. + - If the i-th argument is an array and the j-entry had an illegal value, then `INFO = -(i*100+j)`. + + - If `>0`, the algorithm failed: + + - `=1`, a split was marked by a positive value in `E`. + - `=2`, current block of `Z` not diagonalized after `100*N` iterations (in inner while loop). On exit `Z` holds a QD array with the same eigenvalues as the given `Z`. + - `=3`, termination criterion of outer while loop not met (program created more than N unreduced blocks). + +- `dlasq2()` corresponds to the [LAPACK][LAPACK] routine [`dlasq2`][lapack-dlasq2]. + +
+ + + +
+ +## Examples + + + +```javascript +var Float64Array = require( '@stdlib/array/float64' ); +var dlasq2 = require( '@stdlib/lapack/base/dlasq2' ); + +var Z = new Float64Array( [ 144, 6, 121, 5, 100, 4, 81, 3, 64, 2, 49, 1.5, 36, 1, 25, 0.5 ] ); + +var info = dlasq2( 4, Z ); +console.log( Z ); +console.log( info ); +``` + +
+ + + + + +* * * + +
+ +## C APIs + + + +
+ +
+ + + + + +
+ +### Usage + +```c +TODO +``` + +#### TODO + +TODO. + +```c +TODO +``` + +TODO + +```c +TODO +``` + +
+ + + + + +
+ +
+ + + + + +
+ +### Examples + +```c +TODO +``` + +
+ + + +
+ + + + + + + + + + + + + + diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/benchmark/benchmark.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/benchmark/benchmark.js new file mode 100644 index 000000000000..8b481e34d83c --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/benchmark/benchmark.js @@ -0,0 +1,102 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 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 bench = require( '@stdlib/bench' ); +var uniform = require( '@stdlib/random/array/uniform' ); +var isnan = require( '@stdlib/math/base/assert/is-nan' ); +var pow = require( '@stdlib/math/base/special/pow' ); +var format = require( '@stdlib/string/format' ); +var pkg = require( './../package.json' ).name; +var dlasq2 = require( './../lib/dlasq2.js' ); + + +// VARIABLES // + +var options = { + 'dtype': 'float64' +}; + + +// FUNCTIONS // + +/** +* Creates a benchmark function. +* +* @private +* @param {PositiveInteger} len - array length +* @returns {Function} benchmark function +*/ +function createBenchmark( len ) { + var Z = uniform( len*4, 0.0, 100.0, options ); + return benchmark; + + /** + * Benchmark function. + * + * @private + * @param {Benchmark} b - benchmark instance + */ + function benchmark( b ) { + var i; + + b.tic(); + for ( i = 0; i < b.iterations; i++ ) { + dlasq2( len, Z ); + if ( isnan( Z[ i%Z.length ] ) ) { + b.fail( 'unexpected NaN value' ); + } + } + b.toc(); + if ( isnan( Z[ i%Z.length ] ) ) { + b.fail( 'unexpected NaN value' ); + } + b.pass( 'benchmark finished' ); + b.end(); + } +} + + +// MAIN // + +/** +* Main execution sequence. +* +* @private +*/ +function main() { + var len; + var min; + var max; + var f; + var i; + + min = 1; // 10^min + max = 6; // 10^max + + for ( i = min; i <= max; i++ ) { + len = pow( 10, i ); + f = createBenchmark( len ); + bench( format( '%s:len=%d', pkg, len ), f ); + } +} + +main(); diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/benchmark/benchmark.ndarray.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/benchmark/benchmark.ndarray.js new file mode 100644 index 000000000000..51500e22b5c4 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/benchmark/benchmark.ndarray.js @@ -0,0 +1,102 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 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 bench = require( '@stdlib/bench' ); +var uniform = require( '@stdlib/random/array/uniform' ); +var isnan = require( '@stdlib/math/base/assert/is-nan' ); +var pow = require( '@stdlib/math/base/special/pow' ); +var format = require( '@stdlib/string/format' ); +var pkg = require( './../package.json' ).name; +var dlasq2 = require( './../lib/ndarray.js' ); + + +// VARIABLES // + +var options = { + 'dtype': 'float64' +}; + + +// FUNCTIONS // + +/** +* Creates a benchmark function. +* +* @private +* @param {PositiveInteger} len - array length +* @returns {Function} benchmark function +*/ +function createBenchmark( len ) { + var Z = uniform( len*4, 0.0, 100.0, options ); + return benchmark; + + /** + * Benchmark function. + * + * @private + * @param {Benchmark} b - benchmark instance + */ + function benchmark( b ) { + var i; + + b.tic(); + for ( i = 0; i < b.iterations; i++ ) { + dlasq2( len, Z, 1, 0 ); + if ( isnan( Z[ i%Z.length ] ) ) { + b.fail( 'unexpected NaN value' ); + } + } + b.toc(); + if ( isnan( Z[ i%Z.length ] ) ) { + b.fail( 'unexpected NaN value' ); + } + b.pass( 'benchmark finished' ); + b.end(); + } +} + + +// MAIN // + +/** +* Main execution sequence. +* +* @private +*/ +function main() { + var len; + var min; + var max; + var f; + var i; + + min = 1; // 10^min + max = 6; // 10^max + + for ( i = min; i <= max; i++ ) { + len = pow( 10, i ); + f = createBenchmark( len ); + bench( format( '%s:ndarray:len=%d', pkg, len ), f ); + } +} + +main(); diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/docs/repl.txt b/lib/node_modules/@stdlib/lapack/base/dlasq2/docs/repl.txt new file mode 100644 index 000000000000..4159a6a8fd0d --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/docs/repl.txt @@ -0,0 +1,113 @@ + +{{alias}}( N, Z ) + Computes all the eigenvalues of the symmetric positive definite + tri-diagonal matrix associated with the QD Array `Z` to high relative + accuracy. + + Indexing is relative to the first index. To introduce an offset, use typed + array views. + + The function mutates `Z`. + + Parameters + ---------- + N: integer + Number of rows/columns in the tri-diagonal matrix. + + Z: Float64Array + QD array. Must have at least `4 * N` elements. + Upon return, the first `N` indexed elements contain the eigenvalues + in descending order, and `Z( 2*N+1 )` holds the trace, and `Z( 2*N+2 )` + holds the sum of the eigenvalues. + If `N > 2`, then `Z( 2*N+3 )` holds the iteration count, `Z( 2*N+4 )` + holds `NDIVS/NIN^2`, and `Z( 2*N+5 )` holds the percentage of shifts + that failed. + + Returns + ------- + info: integer + Status code. The status code indicates the following conditions: + + - if equal to zero, then the algorithm completed successfully. + - if smaller than zero: + - if the i-th argument is a scalar and had an illegal value, then + INFO = -i. + - if the i-th argument is an array and the j-entry had an + illegal value, then INFO = -(i*100+j). + - if greater than zero, the algorithm failed: + - equal to 1, a split was marked by a positive value in `E`. + - equal to 2, current block of `Z` not diagonalized after `100*N` + iterations (in inner while loop). On exit `Z` holds a QD array with + the same eigenvalues as the given `Z`. + - equal to 3, termination criterion of outer while loop not met + (program created more than N unreduced blocks). + + Examples + -------- + > var Z = new {{alias:@stdlib/array/float64}}( [ 100,4,81,3,64,2.5,49,2 ] ); + > {{alias}}( 2, Z ) + 0 + > Z + [ ~113.86, ~71.14, ~71.14, 3, 185, 185, 49, 2 ] + + +{{alias}}.ndarray( N, Z, strideZ, offsetZ ) + Computes all the eigenvalues of the symmetric positive definite + tri-diagonal matrix associated with the QD Array `Z` to high relative + accuracy using alternative indexing semantics. + + While typed array views mandate a view offset based on the underlying + buffer, the offset parameters support indexing semantics based on starting + indices. + + The function mutates `Z`. + + Parameters + ---------- + N: integer + Number of rows/columns in the tri-diagonal matrix. + + Z: Float64Array + QD array. Must have at least `4 * N` elements. + Upon return, the first `N` indexed elements contain the eigenvalues + in descending order, and `Z( 2*N+1 )` holds the trace, and `Z( 2*N+2 )` + holds the sum of the eigenvalues. + If `N > 2`, then `Z( 2*N+3 )` holds the iteration count, `Z( 2*N+4 )` + holds `NDIVS/NIN^2`, and `Z( 2*N+5 )` holds the percentage of shifts + that failed. + + strideZ: integer + Stride length for `Z`. + + offsetZ: integer + Starting index for `Z`. + + Returns + ------- + info: integer + Status code. The status code indicates the following conditions: + + - if equal to zero, then the algorithm completed successfully. + - if smaller than zero: + - if the i-th argument is a scalar and had an illegal value, then + INFO = -i. + - if the i-th argument is an array and the j-entry had an + illegal value, then INFO = -(i*100+j). + - if greater than zero, the algorithm failed: + - equal to 1, a split was marked by a positive value in `E`. + - equal to 2, current block of `Z` not diagonalized after `100*N` + iterations (in inner while loop). On exit `Z` holds a QD array with + the same eigenvalues as the given `Z`. + - equal to 3, termination criterion of outer while loop not met + (program created more than N unreduced blocks). + + Examples + -------- + > var Z = new {{alias:@stdlib/array/float64}}( [ 100,4,81,3,64,2.5,49,2 ] ); + > {{alias}}.ndarray( 2, Z, 1, 0 ) + 0 + > Z + [ ~113.86, ~71.14, ~71.14, 3, 185, 185, 49, 2 ] + + See Also + -------- diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/docs/types/index.d.ts b/lib/node_modules/@stdlib/lapack/base/dlasq2/docs/types/index.d.ts new file mode 100644 index 000000000000..db795383d128 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/docs/types/index.d.ts @@ -0,0 +1,111 @@ +/* +* @license Apache-2.0 +* +* Copyright (c) 2026 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. +*/ + +// TypeScript Version: 4.1 + +/** +* Status code. +* +* ## Notes +* +* The status code indicates the following conditions: +* +* - if equal to zero, then the algorithm completed successfully. +* - if smaller than zero: +* +* - if the i-th argument is a scalar and had an illegal value, then INFO = -i. +* - if the i-th argument is an array and the j-entry had an illegal value, then INFO = -(i*100+j). +* +* - if greater than zero, the algorithm failed: +* +* - equal to 1, a split was marked by a positive value in `E`. +* - equal to 2, current block of `Z` not diagonalized after `100*N` iterations (in inner while loop). On exit `Z` holds a QD array with the same eigenvalues as the given `Z`. +* - equal to 3, termination criterion of outer while loo +*/ +type StatusCode = number; + +/** +* Interface describing `dlasq2`. +*/ +interface Routine { + /** + * Computes all the eigenvalues of the symmetric positive definite tri-diagonal matrix associated with the QD Array `Z` to high relative accuracy. + * + * @param N - number of rows/columns in the tri-diagonal matrix + * @param Z - QD array + * @returns status code + * + * @example + * var Float64Array = require( '@stdlib/array/float64' ); + * + * var Z = new Float64Array( [ 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); + * + * dlasq2( 4, Z ); + * // Z => [ ~115.713, ~83.153, ~62.17, ~42.464, ~83.153, ~20.987, ~0.0, ~0.0, 303.5, ~303.5, 15, 4.0625, 0.0, ~0.0, ~0.172, ~0.096 ] + */ + ( N: number, Z: Float64Array ): StatusCode; + + /** + * Computes all the eigenvalues of the symmetric positive definite tri-diagonal matrix associated with the QD Array `Z` to high relative accuracy using alternative indexing semantics. + * + * @param N - number of rows/columns in the tri-diagonal matrix + * @param Z - QD array + * @param strideZ - stride length for `Z` + * @param offsetZ - starting index of `Z` + * @returns status code + * + * @example + * var Float64Array = require( '@stdlib/array/float64' ); + * + * var Z = new Float64Array( [ 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); + * + * dlasq2.ndarray( 4, Z, 1, 0 ); + * // Z => [ ~115.713, ~83.153, ~62.17, ~42.464, ~83.153, ~20.987, ~0.0, ~0.0, 303.5, ~303.5, 15, 4.0625, 0.0, ~0.0, ~0.172, ~0.096 ] + */ + ndarray( N: number, Z: Float64Array, strideZ: number, offsetZ: number ): StatusCode; +} + +/** +* LAPACK routine to compute all eigenvalues of the symmetric positive definite tri-diagonal matrix associated with the QD Array `Z` to high relative accuracy. +* +* @param N - number of rows/columns in the tri-diagonal matrix +* @param Z - QD array +* @returns status code +* +* @example +* var Float64Array = require( '@stdlib/array/float64' ); +* +* var Z = new Float64Array( [ 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); +* +* dlasq2( 4, Z ); +* // Z => [ ~115.713, ~83.153, ~62.17, ~42.464, ~83.153, ~20.987, ~0.0, ~0.0, 303.5, ~303.5, 15, 4.0625, 0.0, ~0.0, ~0.172, ~0.096 ] +* +* @example +* var Float64Array = require( '@stdlib/array/float64' ); +* +* var Z = new Float64Array( [ 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); +* +* dlasq2.ndarray( 4, Z, 1, 0 ); +* // Z => [ ~115.713, ~83.153, ~62.17, ~42.464, ~83.153, ~20.987, ~0.0, ~0.0, 303.5, ~303.5, 15, 4.0625, 0.0, ~0.0, ~0.172, ~0.096 ] +*/ +declare var dlasq2: Routine; + + +// EXPORTS // + +export = dlasq2; diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/docs/types/test.ts b/lib/node_modules/@stdlib/lapack/base/dlasq2/docs/types/test.ts new file mode 100644 index 000000000000..82a9a872935d --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/docs/types/test.ts @@ -0,0 +1,138 @@ +/* +* @license Apache-2.0 +* +* Copyright (c) 2026 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. +*/ + +import dlasq2 = require( './index' ); + + +// TESTS // + +// The function returns a number... +{ + const Z = new Float64Array( 16 ); + + dlasq2( 4, Z ); // $ExpectType number +} + +// The compiler throws an error if the function is provided a first argument which is not a number... +{ + const Z = new Float64Array( 16 ); + + dlasq2( '5', Z ); // $ExpectError + dlasq2( true, Z ); // $ExpectError + dlasq2( false, Z ); // $ExpectError + dlasq2( null, Z ); // $ExpectError + dlasq2( void 0, Z ); // $ExpectError + dlasq2( [], Z ); // $ExpectError + dlasq2( {}, Z ); // $ExpectError + dlasq2( ( x: number ): number => x, Z ); // $ExpectError +} + +// The compiler throws an error if the function is provided a second argument which is not a Float64Array... +{ + dlasq2( 4, '5' ); // $ExpectError + dlasq2( 4, 5 ); // $ExpectError + dlasq2( 4, true ); // $ExpectError + dlasq2( 4, false ); // $ExpectError + dlasq2( 4, null ); // $ExpectError + dlasq2( 4, void 0 ); // $ExpectError + dlasq2( 4, [] ); // $ExpectError + dlasq2( 4, {} ); // $ExpectError + dlasq2( 4, ( x: number ): number => x ); // $ExpectError +} + +// The compiler throws an error if the function is provided an unsupported number of arguments... +{ + const Z = new Float64Array( 16 ); + + dlasq2(); // $ExpectError + dlasq2( 4 ); // $ExpectError + dlasq2( 4, Z, 1 ); // $ExpectError +} + +// Attached to main export is an `ndarray` method which returns a number... +{ + const Z = new Float64Array( 16 ); + + dlasq2.ndarray( 4, Z, 1, 0 ); // $ExpectType number +} + +// The compiler throws an error if the `ndarray` method is provided a first argument which is not a number... +{ + const Z = new Float64Array( 16 ); + + dlasq2.ndarray( '5', Z, 1, 0 ); // $ExpectError + dlasq2.ndarray( true, Z, 1, 0 ); // $ExpectError + dlasq2.ndarray( false, Z, 1, 0 ); // $ExpectError + dlasq2.ndarray( null, Z, 1, 0 ); // $ExpectError + dlasq2.ndarray( void 0, Z, 1, 0 ); // $ExpectError + dlasq2.ndarray( [], Z, 1, 0 ); // $ExpectError + dlasq2.ndarray( {}, Z, 1, 0 ); // $ExpectError + dlasq2.ndarray( ( x: number ): number => x, Z, 1, 0 ); // $ExpectError +} + +// The compiler throws an error if the `ndarray` method is provided a second argument which is not a Float64Array... +{ + dlasq2.ndarray( 4, '5', 1, 0 ); // $ExpectError + dlasq2.ndarray( 4, 5, 1, 0 ); // $ExpectError + dlasq2.ndarray( 4, true, 1, 0 ); // $ExpectError + dlasq2.ndarray( 4, false, 1, 0 ); // $ExpectError + dlasq2.ndarray( 4, null, 1, 0 ); // $ExpectError + dlasq2.ndarray( 4, void 0, 1, 0 ); // $ExpectError + dlasq2.ndarray( 4, [], 1, 0 ); // $ExpectError + dlasq2.ndarray( 4, {}, 1, 0 ); // $ExpectError + dlasq2.ndarray( 4, ( x: number ): number => x, 1, 0 ); // $ExpectError +} + +// The compiler throws an error if the `ndarray` method is provided a third argument which is not a number... +{ + const Z = new Float64Array( 16 ); + + dlasq2.ndarray( 4, Z, '5', 0 ); // $ExpectError + dlasq2.ndarray( 4, Z, true, 0 ); // $ExpectError + dlasq2.ndarray( 4, Z, false, 0 ); // $ExpectError + dlasq2.ndarray( 4, Z, null, 0 ); // $ExpectError + dlasq2.ndarray( 4, Z, void 0, 0 ); // $ExpectError + dlasq2.ndarray( 4, Z, [], 0 ); // $ExpectError + dlasq2.ndarray( 4, Z, {}, 0 ); // $ExpectError + dlasq2.ndarray( 4, Z, ( x: number ): number => x, 0 ); // $ExpectError +} + +// The compiler throws an error if the `ndarray` method is provided a fourth argument which is not a number... +{ + const Z = new Float64Array( 16 ); + + dlasq2.ndarray( 4, Z, 1, '5' ); // $ExpectError + dlasq2.ndarray( 4, Z, 1, true ); // $ExpectError + dlasq2.ndarray( 4, Z, 1, false ); // $ExpectError + dlasq2.ndarray( 4, Z, 1, null ); // $ExpectError + dlasq2.ndarray( 4, Z, 1, void 0 ); // $ExpectError + dlasq2.ndarray( 4, Z, 1, [] ); // $ExpectError + dlasq2.ndarray( 4, Z, 1, {} ); // $ExpectError + dlasq2.ndarray( 4, Z, 1, ( x: number ): number => x ); // $ExpectError +} + +// The compiler throws an error if the `ndarray` method is provided an unsupported number of arguments... +{ + const Z = new Float64Array( 16 ); + + dlasq2.ndarray(); // $ExpectError + dlasq2.ndarray( 4 ); // $ExpectError + dlasq2.ndarray( 4, Z ); // $ExpectError + dlasq2.ndarray( 4, Z, 1 ); // $ExpectError + dlasq2.ndarray( 4, Z, 1, 0, 10 ); // $ExpectError +} diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/examples/index.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/examples/index.js new file mode 100644 index 000000000000..32a9b33640bb --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/examples/index.js @@ -0,0 +1,28 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 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'; + +var Float64Array = require( '@stdlib/array/float64' ); +var dlasq2 = require( './../lib' ); + +var Z = new Float64Array( [ 144, 6, 121, 5, 100, 4, 81, 3, 64, 2, 49, 1.5, 36, 1, 25, 0.5 ] ); + +var info = dlasq2( 4, Z ); +console.log( Z ); +console.log( info ); diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/base.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/base.js new file mode 100644 index 000000000000..b0fd0890d97b --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/base.js @@ -0,0 +1,527 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 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'; + +/* eslint-disable max-len, max-statements, max-lines-per-function */ + +// MODULES // + +var Float64Array = require( '@stdlib/array/float64' ); +var dlamch = require( '@stdlib/lapack/base/dlamch' ); +var pow = require( '@stdlib/math/base/special/pow' ); +var max = require( '@stdlib/math/base/special/max' ); +var min = require( '@stdlib/math/base/special/min' ); +var sqrt = require( '@stdlib/math/base/special/sqrt' ); +var abs = require( '@stdlib/math/base/special/abs' ); +var floor = require( '@stdlib/math/base/special/floor' ); +var dlasq3 = require( './dlasq3.js' ); +var dlasrt = require( './dlasrt.js' ); + + +// VARIABLES // + +var CBIAS = 1.50; +var EPS = dlamch( 'P' ); +var TOL = EPS * 100; +var TOL2 = pow( TOL, 2 ); +var SAFMIN = dlamch( 'S' ); +var IEEE = true; + + +// MAIN // + +/** +* Computes all the eigenvalues of the symmetric positive definite tri-diagonal matrix associated with the QD Array `Z` to high relative accuracy. +* +* @private +* @param {integer} N - number of rows/columns in `Z` +* @param {Float64Array} Z - qd array +* @param {integer} strideZ - stride length for `Z` +* @param {NonNegativeInteger} offsetZ - starting index of `Z` +* @returns {integer} status code +* +* @example +* var Float64Array = require( '@stdlib/array/float64' ); +* +* var Z = new Float64Array( [ 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); +* +* var out = dlasq2( 4, Z, 1, 0 ); +* // Z => [ ~115.713, ~83.153, ~62.17, ~42.464, ~83.153, ~20.987, ~0.0, ~0.0, 303.5, ~303.5, 15, 4.0625, 0.0, ~0.0, ~0.172, ~0.096 ] +* // out => 0 +*/ +function dlasq2( N, Z, strideZ, offsetZ ) { // eslint-disable-line stdlib/jsdoc-doctest-decimal-point + var deemin; + var oldemn; + var iwhila; + var iwhilb; + var desig; + var dmin1; + var dmin2; + var sigma; + var trace; + var tempq; + var tempe; + var ttype; + var nfail; + var emax; + var emin; + var iter; + var nbig; + var ndiv; + var ipn4; + var kmin; + var dmin; + var qmax; + var qmin; + var temp; + var zmax; + var splt; + var dee; + var tau; + var dn1; + var dn2; + var out; + var pp; + var i0; + var i1; + var i4; + var n0; + var n1; + var dn; + var d; + var e; + var g; + var k; + var s; + var t; + + out = new Float64Array( 18 ); + + if ( N === 0 ) { + return 0; + } + + if ( N === 1 ) { + // 1-by-1 case + if ( Z[ offsetZ ] < 0 ) { + return -201; + } + return 0; + } + + if ( N === 2 ) { + // 2-by-2 case + if ( Z[ offsetZ ] < 0 ) { + return -201; + } + + if ( Z[ offsetZ + strideZ ] < 0 ) { + return -202; + } + + if ( Z[ offsetZ + ( 2 * strideZ ) ] < 0 ) { + return -203; + } + + if ( Z[ offsetZ + ( 2 * strideZ ) ] > Z[ offsetZ ] ) { + d = Z[ offsetZ + ( 2 * strideZ ) ]; + Z[ offsetZ + ( 2 * strideZ ) ] = Z[ offsetZ ]; + Z[ offsetZ ] = d; + } + + Z[ offsetZ + ( 4 * strideZ ) ] = Z[ offsetZ ] + Z[ offsetZ + strideZ ] + Z[ offsetZ + ( 2 * strideZ ) ]; + + if ( Z[ offsetZ + strideZ ] > Z[ offsetZ + ( 2 * strideZ ) ] * TOL2 ) { + t = 0.5 * ( ( Z[ offsetZ ] - Z[ offsetZ + ( 2 * strideZ ) ] ) + Z[ offsetZ + strideZ ] ); + s = Z[ offsetZ + ( 2 * strideZ ) ] * ( Z[ offsetZ + strideZ ] / t ); + if ( s <= t ) { + s = Z[ offsetZ + ( 2 * strideZ ) ] * ( Z[ offsetZ + strideZ ] / ( t * ( 1 + sqrt( 1 + ( s / t ) ) ) ) ); + } else { + s = Z[ offsetZ + ( 2 * strideZ ) ] * ( Z[ offsetZ + strideZ ] / ( t + ( sqrt( t ) * sqrt( t + s ) ) ) ); + } + t = Z[ offsetZ ] + ( s + Z[ offsetZ + strideZ ] ); + Z[ offsetZ + ( 2 * strideZ ) ] = Z[ offsetZ + ( 2 * strideZ ) ] * ( Z[ offsetZ ] / t ); + Z[ offsetZ ] = t; + } + Z[ offsetZ + strideZ ] = Z[ offsetZ + ( 2 * strideZ ) ]; + Z[ offsetZ + ( 5 * strideZ ) ] = Z[ offsetZ + strideZ ] + Z[ offsetZ ]; + return 0; + } + + // Check for negative data and compute sums of q's and e's + Z[ offsetZ + ( ( ( 2 * N ) - 1 ) * strideZ ) ] = 0; + emin = Z[ offsetZ + strideZ ]; + qmax = 0; + zmax = 0; + d = 0; + e = 0; + + for ( k = 0; k < 2 * ( N - 1 ); k += 2 ) { + if ( Z[ offsetZ + ( k * strideZ ) ] < 0 ) { + return -( 200 + k ); + } + if ( Z[ offsetZ + ( k * strideZ ) + strideZ ] < 0 ) { + return -( 200 + k + 1 ); + } + d += Z[ offsetZ + ( k * strideZ ) ]; + e += Z[ offsetZ + ( k * strideZ ) + strideZ ]; + qmax = max( qmax, Z[ offsetZ + ( k * strideZ ) ] ); + emin = min( emin, Z[ offsetZ + ( k * strideZ ) + strideZ ] ); + zmax = max( qmax, max( zmax, Z[ offsetZ + ( k * strideZ ) + strideZ ] ) ); + } + if ( Z[ offsetZ + ( ( ( 2 * N ) - 2 ) * strideZ ) ] < 0 ) { + return -( 200 + ( 2 * N ) - 1 ); + } + d += Z[ offsetZ + ( ( ( 2 * N ) - 2 ) * strideZ ) ]; + qmax = max( qmax, Z[ offsetZ + ( ( ( 2 * N ) - 2 ) * strideZ ) ] ); + zmax = max( qmax, zmax ); + + // Check for diagonality + if ( e === 0 ) { + for ( k = 1; k < N; k++ ) { + Z[ offsetZ + ( k * strideZ ) ] = Z[ offsetZ + ( ( ( 2 * k ) + 1 ) * strideZ ) ]; + } + dlasrt( 'decreasing', N, Z, strideZ, offsetZ ); + Z[ offsetZ + ( ( ( 2 * N ) - 2 ) * strideZ ) ] = d; + return 0; + } + + trace = d + e; + + // Check for zero data + if ( trace === 0 ) { + Z[ offsetZ + ( ( ( 2 * N ) - 2 ) * strideZ ) ] = 0; + return 0; + } + + // Check whether the machine is IEEE conformable (In JS, always true) + IEEE = true; + + // Rearrange data for locality: Z=(q1,qq1,e1,ee1,q2,qq2,e2,ee2,...) + for ( k = ( 2 * N ) - 1; k >= 1; k -= 2 ) { + Z[ offsetZ + ( ( ( 2 * k ) + 1 ) * strideZ ) ] = 0; + Z[ offsetZ + ( ( 2 * k ) * strideZ ) ] = Z[ offsetZ + ( k * strideZ ) ]; + Z[ offsetZ + ( ( ( 2 * k ) - 1 ) * strideZ ) ] = 0; + Z[ offsetZ + ( ( ( 2 * k ) - 2 ) * strideZ ) ] = Z[ offsetZ + ( k * strideZ ) - strideZ ]; + } + + i0 = 0; + n0 = N-1; + + // Reverse the qd-array, if warranted + if ( CBIAS * Z[ offsetZ + ( ( 4 * i0 ) * strideZ ) ] < Z[ offsetZ + ( ( 4 * n0 ) * strideZ ) ] ) { + ipn4 = ( 4 * ( i0 + n0 ) ) + 6; + for ( i4 = ( 4 * i0 ) + 3; i4 <= ( 2 * ( i0 + n0 ) ) + 1; i4 += 4 ) { + temp = Z[ offsetZ + ( ( i4 - 3 ) * strideZ ) ]; + Z[ offsetZ + ( ( i4 - 3 ) * strideZ ) ] = Z[ offsetZ + ( ( ipn4 - i4 - 3 ) * strideZ ) ]; + Z[ offsetZ + ( ( ipn4 - i4 - 3 ) * strideZ ) ] = temp; + temp = Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ]; + Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ] = Z[ offsetZ + ( ( ipn4 - i4 - 5 ) * strideZ ) ]; + Z[ offsetZ + ( ( ipn4 - i4 - 5 ) * strideZ ) ] = temp; + } + } + + // Initial split checking via DQD and Li's test + pp = 0; + + for ( k = 0; k < 2; k++ ) { + d = Z[ offsetZ + ( ( ( 4 * n0 ) + pp ) * strideZ ) ]; + for ( i4 = ( 4 * n0 ) - 1 + pp; i4 >= ( 4 * i0 ) + pp + 3; i4 -= 4 ) { + if ( Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ] <= TOL2 * d ) { + Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ] = -0; + d = Z[ offsetZ + ( ( i4 - 3 ) * strideZ ) ]; + } else { + d = Z[ offsetZ + ( ( i4 - 3 ) * strideZ ) ] * ( d / ( d + Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ] ) ); + } + } + + // DQD maps Z to ZZ plus Li's test + emin = Z[ offsetZ + ( ( ( 4 * i0 ) + pp + 4 ) * strideZ ) ]; + d = Z[ offsetZ + ( ( ( 4 * i0 ) + pp ) * strideZ ) ]; + for ( i4 = ( 4 * i0 ) + pp + 3; i4 <= ( 4 * n0 ) + pp - 1; i4 += 4 ) { + Z[ offsetZ + ( ( i4 - ( 2 * pp ) - 2 ) * strideZ ) ] = d + Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ]; + if ( Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ] <= TOL2 * d ) { + Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ] = -0; + Z[ offsetZ + ( ( i4 - ( 2 * pp ) - 2 ) * strideZ ) ] = d; + Z[ offsetZ + ( ( i4 - ( 2 * pp ) ) * strideZ ) ] = 0; + d = Z[ offsetZ + ( ( i4 + 1 ) * strideZ ) ]; + } else if ( ( SAFMIN * Z[ offsetZ + ( ( i4 + 1 ) * strideZ ) ] < Z[ offsetZ + ( ( i4 - ( 2 * pp ) - 2 ) * strideZ ) ] ) && ( SAFMIN * Z[ offsetZ + ( ( i4 - ( 2 * pp ) - 2 ) * strideZ ) ] < Z[ offsetZ + ( ( i4 + 1 ) * strideZ ) ] ) ) { + temp = Z[ offsetZ + ( ( i4 + 1 ) * strideZ ) ] / Z[ offsetZ + ( ( i4 - ( 2 * pp ) - 2 ) * strideZ ) ]; + Z[ offsetZ + ( ( i4 - ( 2 * pp ) ) * strideZ ) ] = Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ] * temp; + d *= temp; + } else { + Z[ offsetZ + ( ( i4 - ( 2 * pp ) ) * strideZ ) ] = Z[ offsetZ + ( ( i4 + 1 ) * strideZ ) ] * ( Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ] / Z[ offsetZ + ( ( i4 - ( 2 * pp ) - 2 ) * strideZ ) ] ); + d = Z[ offsetZ + ( ( i4 + 1 ) * strideZ ) ] * ( d / Z[ offsetZ + ( ( i4 - ( 2 * pp ) - 2 ) * strideZ ) ] ); + } + emin = min( emin, Z[ offsetZ + ( ( i4 - ( 2 * pp ) ) * strideZ ) ] ); + } + Z[ offsetZ + ( ( ( 4 * n0 ) - pp + 1 ) * strideZ ) ] = d; + + // Now find qmax + qmax = Z[ offsetZ + ( ( ( 4 * i0 ) - pp + 1 ) * strideZ ) ]; + for ( i4 = ( 4 * i0 ) - pp + 5; i4 <= ( 4 * n0 ) - pp + 1; i4 += 4 ) { + qmax = max( qmax, Z[ offsetZ + ( i4 * strideZ ) ] ); + } + + // Prepare for the next iteration on K + pp = 1 - pp; + } + + // Initialize variables to pass to DLASQ3 + ttype = 0; + dmin1 = 0; + dmin2 = 0; + dn = 0; + dn1 = 0; + dn2 = 0; + g = 0; + tau = 0; + + iter = 2; + nfail = 0; + ndiv = 2 * ( n0 - i0 ); + + for ( iwhila = 0; iwhila <= N; iwhila++ ) { + if ( n0 < 0 ) { + // Move q's to the front + for ( k = 1; k < N; k++ ) { + Z[ offsetZ + ( k * strideZ ) ] = Z[ offsetZ + ( ( 4 * k ) * strideZ ) ]; + } + + // Sort and compute sum of eigenvalues + dlasrt( 'decreasing', N, Z, strideZ, offsetZ ); + + e = 0; + for ( k = N - 1; k >= 0; k-- ) { + e += Z[ offsetZ + ( k * strideZ ) ]; + } + + // Store trace, sum(eigenvalues) and information on performance + Z[ offsetZ + ( ( 2 * N ) * strideZ ) ] = trace; + Z[ offsetZ + ( ( ( 2 * N ) + 1 ) * strideZ ) ] = e; + Z[ offsetZ + ( ( ( 2 * N ) + 2 ) * strideZ ) ] = iter; + Z[ offsetZ + ( ( ( 2 * N ) + 3 ) * strideZ ) ] = ndiv / ( N * N ); + Z[ offsetZ + ( ( ( 2 * N ) + 4 ) * strideZ ) ] = ( 100 * nfail ) / iter; + return 0; + } + + // E(N0) holds the value of SIGMA when submatrix in I0:N0 splits from the rest of the array, but is negated. + desig = 0; + if ( n0 === N - 1 ) { + sigma = 0; + } else { + sigma = -Z[ offsetZ + ( ( ( 4 * n0 ) + 2 ) * strideZ ) ]; + } + if ( sigma < 0 ) { + return 1; + } + + // Find last unreduced submatrix's top index I0, find QMAX and EMIN. Find Gershgorin-type bound if Q's much greater than E's. + emax = 0; + if ( n0 > i0 ) { + emin = abs( Z[ offsetZ + ( ( ( 4 * n0 ) - 2 ) * strideZ ) ] ); + } else { + emin = 0; + } + qmin = Z[ offsetZ + ( ( 4 * n0 ) * strideZ ) ]; + qmax = qmin; + + for ( i4 = ( 4 * n0 ) + 3; i4 >= 7; i4 -= 4 ) { + if ( Z[ offsetZ + ( ( i4 - 5 ) * strideZ ) ] <= 0 ) { + break; + } + if ( qmin >= 4 * emax ) { + qmin = min( qmin, Z[ offsetZ + ( ( i4 - 3 ) * strideZ ) ] ); + emax = max( emax, Z[ offsetZ + ( ( i4 - 5 ) * strideZ ) ] ); + } + qmax = max( qmax, Z[ offsetZ + ( ( i4 - 7 ) * strideZ ) ] + Z[ offsetZ + ( ( i4 - 5 ) * strideZ ) ] ); + emin = min( emin, Z[ offsetZ + ( ( i4 - 5 ) * strideZ ) ] ); + } + + // If the loop completed without break, set i4 = 4 + if ( i4 < 7 ) { + i4 = 3; + } + + i0 = floor( ( ( i4 + 1 ) / 4 ) - 1 ); + pp = 0; + + if ( n0 - i0 > 1 ) { + dee = Z[ offsetZ + ( ( 4 * i0 ) * strideZ ) ]; + deemin = dee; + kmin = i0; + for ( i4 = ( 4 * i0 ) + 4; i4 <= 4 * n0; i4 += 4 ) { + dee = Z[ offsetZ + ( i4 * strideZ ) ] * ( dee / ( dee + Z[ offsetZ + ( ( i4 - 2 ) * strideZ ) ] ) ); + if ( dee <= deemin ) { + deemin = dee; + kmin = floor( i4 / 4 ); + } + } + if ( ( ( kmin - i0 ) * 2 < n0 - kmin ) && ( deemin <= 0.5 * Z[ offsetZ + ( ( 4 * n0 ) * strideZ ) ] ) ) { + ipn4 = ( 4 * ( i0 + n0 ) ) + 7; + pp = 2; + for ( i4 = ( 4 * i0 ) + 3; i4 <= ( 2 * ( i0 + n0 ) ) + 1; i4 += 4 ) { + temp = Z[ offsetZ + ( ( i4 - 3 ) * strideZ ) ]; + Z[ offsetZ + ( ( i4 - 3 ) * strideZ ) ] = Z[ offsetZ + ( ( ipn4 - i4 - 3 ) * strideZ ) ]; + Z[ offsetZ + ( ( ipn4 - i4 - 3 ) * strideZ ) ] = temp; + temp = Z[ offsetZ + ( ( i4 - 2 ) * strideZ ) ]; + Z[ offsetZ + ( ( i4 - 2 ) * strideZ ) ] = Z[ offsetZ + ( ( ipn4 - i4 - 2 ) * strideZ ) ]; + Z[ offsetZ + ( ( ipn4 - i4 - 2 ) * strideZ ) ] = temp; + temp = Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ]; + Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ] = Z[ offsetZ + ( ( ipn4 - i4 - 5 ) * strideZ ) ]; + Z[ offsetZ + ( ( ipn4 - i4 - 5 ) * strideZ ) ] = temp; + temp = Z[ offsetZ + ( i4 * strideZ ) ]; + Z[ offsetZ + ( i4 * strideZ ) ] = Z[ offsetZ + ( ( ipn4 - i4 - 4 ) * strideZ ) ]; + Z[ offsetZ + ( ( ipn4 - i4 - 4 ) * strideZ ) ] = temp; + } + } + } + + // Put -(initial shift) into DMIN + dmin = -max( 0, qmin - ( 2 * sqrt( qmin ) * sqrt( emax ) ) ); + + /* + Now I0:N0 is unreduced. + PP = 0 for ping, + PP = 1 for pong, + PP = 2 indicates that flipping was applied to the Z array and that the tests for deflation upon entry in DLASQ3 should not be performed. + */ + nbig = 100 * ( n0 - i0 + 1 ); + for ( iwhilb = 0; iwhilb < nbig; iwhilb++ ) { + if ( i0 > n0 ) { + break; + } + + out[ 0 ] = n0; + out[ 1 ] = pp; + out[ 2 ] = dmin; + out[ 3 ] = sigma; + out[ 4 ] = desig; + out[ 5 ] = qmax; + out[ 6 ] = nfail; + out[ 7 ] = iter; + out[ 8 ] = ndiv; + out[ 9 ] = ttype; + out[ 10 ] = dmin1; + out[ 11 ] = dmin2; + out[ 12 ] = dn; + out[ 13 ] = dn1; + out[ 14 ] = dn2; + out[ 15 ] = g; + out[ 16 ] = tau; + + // Call dlasq3 + dlasq3( i0, Z, strideZ, offsetZ, IEEE, out, 1, 0 ); + + n0 = out[ 0 ]; + pp = out[ 1 ]; + dmin = out[ 2 ]; + sigma = out[ 3 ]; + desig = out[ 4 ]; + qmax = out[ 5 ]; + nfail = out[ 6 ]; + iter = out[ 7 ]; + ndiv = out[ 8 ]; + ttype = out[ 9 ]; + dmin1 = out[ 10 ]; + dmin2 = out[ 11 ]; + dn = out[ 12 ]; + dn1 = out[ 13 ]; + dn2 = out[ 14 ]; + g = out[ 15 ]; + tau = out[ 16 ]; + + pp = 1 - pp; + + // When EMIN is very small check for splits + if ( pp === 0 && n0 - i0 >= 3 ) { + if ( ( Z[ offsetZ + ( ( ( 4 * n0 ) + 3 ) * strideZ ) ] <= TOL2 * qmax ) || ( Z[ offsetZ + ( ( ( 4 * n0 ) + 2 ) * strideZ ) ] <= TOL2 * sigma ) ) { + splt = i0 - 1; + qmax = Z[ offsetZ + ( ( 4 * i0 ) * strideZ ) ]; + emin = Z[ offsetZ + ( ( ( 4 * i0 ) + 2 ) * strideZ ) ]; + oldemn = Z[ offsetZ + ( ( ( 4 * i0 ) + 3 ) * strideZ ) ]; + for ( i4 = ( 4 * i0 ) + 3; i4 <= ( 4 * n0 ) - 9; i4 += 4 ) { + if ( ( Z[ offsetZ + ( i4 * strideZ ) ] <= TOL2 * Z[ offsetZ + ( ( i4 - 3 ) * strideZ ) ] ) || ( Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ] <= TOL2 * sigma) ) { + Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ] = -sigma; + splt = floor( ( i4 + 1 ) / 4 ) - 1; + qmax = 0; + emin = Z[ offsetZ + ( ( i4 + 3 ) * strideZ ) ]; + oldemn = Z[ offsetZ + ( ( i4 + 4 ) * strideZ ) ]; + } else { + qmax = max( qmax, Z[ offsetZ + ( ( i4 + 1 ) * strideZ ) ] ); + emin = min( emin, Z[ offsetZ + ( ( i4 - 1 ) * strideZ ) ] ); + oldemn = min( oldemn, Z[ offsetZ + ( i4 * strideZ ) ] ); + } + } + Z[ offsetZ + ( ( ( 4 * n0 ) + 2 ) * strideZ ) ] = emin; + Z[ offsetZ + ( ( ( 4 * n0 ) + 3 ) * strideZ ) ] = oldemn; + i0 = splt + 1; + } + } + } + + // If inner loop exhausted without breaking, we have INFO = 2 + if ( iwhilb > nbig ) { + // Maximum number of iterations exceeded, restore the shift SIGMA and place the new d's and e's in a qd array. This might need to be done for several blocks + + i1 = i0; + n1 = n0; + + // Label 145 loop + while ( true ) { + tempq = Z[ offsetZ + ( ( 4 * i0 ) * strideZ ) ]; + Z[ offsetZ + ( ( 4 * i0 ) * strideZ ) ] += sigma; + for ( k = i0 + 1; k <= n0; k++ ) { + tempe = Z[ offsetZ + ( ( ( 4 * k ) - 2 ) * strideZ ) ]; + Z[ offsetZ + ( ( ( 4 * k ) - 2 ) * strideZ ) ] *= tempq / Z[ offsetZ + ( ( ( 4 * k ) - 4 ) * strideZ ) ]; + tempq = Z[ offsetZ + ( ( 4 * k ) * strideZ ) ]; + Z[ offsetZ + ( ( 4 * k ) * strideZ ) ] += sigma + tempe - Z[ offsetZ + ( ( ( 4 * k ) - 2 ) * strideZ ) ]; + } + + // Prepare to do this on the previous block if there is one + if ( i1 > 0 ) { + n1 = i1 - 1; + while ( i1 >= 1 && Z[ offsetZ + ( ( ( 4 * i1 ) - 2 ) * strideZ ) ] >= 0 ) { + i1 -= 1; + } + sigma = -Z[ offsetZ + ( ( ( 4 * n1 ) + 2 ) * strideZ ) ]; + } else { + break; + } + } + + for ( k = 0; k < N; k++ ) { + Z[ offsetZ + ( ( 2 * k ) * strideZ ) ] = Z[ offsetZ + ( ( 4 * k ) * strideZ ) ]; + if ( k < n0 ) { + Z[ offsetZ + ( ( ( 2 * k ) + 1 ) * strideZ ) ] = Z[ offsetZ + ( ( ( 4 * k ) + 2 ) * strideZ ) ]; + } else { + Z[ offsetZ + ( ( ( 2 * k ) + 1 ) * strideZ ) ] = 0; + } + } + return 2; + } + } + + return 3; +} + + +// EXPORTS // + +module.exports = dlasq2; diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq2.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq2.js new file mode 100644 index 000000000000..ad9b7055a707 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq2.js @@ -0,0 +1,52 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 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 base = require( './base.js' ); + + +// MAIN // + +/** +* Computes all the eigenvalues of the symmetric positive definite tri-diagonal matrix associated with the QD Array `Z` to high relative accuracy. +* +* @private +* @param {integer} N - number of rows/columns in `Z` +* @param {Float64Array} Z - qd array +* @returns {integer} status code +* +* @example +* var Float64Array = require( '@stdlib/array/float64' ); +* +* var Z = new Float64Array( [ 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); +* +* var out = dlasq2( 4, Z ); +* // Z => [ ~115.713, ~83.153, ~62.17, ~42.464, ~83.153, ~20.987, ~0.0, ~0.0, 303.5, ~303.5, 15, 4.0625, 0.0, ~0.0, ~0.172, ~0.096 ] +* // out => 0 +*/ +function dlasq2( N, Z ) { // eslint-disable-line stdlib/jsdoc-doctest-decimal-point + return base( N, Z, 1, 0 ); +} + + +// EXPORTS // + +module.exports = dlasq2; diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq3.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq3.js new file mode 100644 index 000000000000..5b4afe58b8e1 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq3.js @@ -0,0 +1,379 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 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'; + +/* eslint-disable max-len, max-statements, max-lines-per-function */ + +// MODULES // + +var Float64Array = require( '@stdlib/array/float64' ); +var dlamch = require( '@stdlib/lapack/base/dlamch' ); +var isnan = require( '@stdlib/assert/is-nan' ); +var pow = require( '@stdlib/math/base/special/pow' ); +var max = require( '@stdlib/math/base/special/max' ); +var min = require( '@stdlib/math/base/special/min' ); +var sqrt = require( '@stdlib/math/base/special/sqrt' ); +var abs = require( '@stdlib/math/base/special/abs' ); +var dlasq4 = require( './dlasq4.js' ); +var dlasq5 = require( './dlasq5.js' ); +var dlasq6 = require( './dlasq6.js' ); + + +// VARIABLES // + +var CBIAS = 1.50; +var EPS = dlamch( 'P' ); +var TOL = EPS * 100; +var TOL2 = pow( TOL, 2 ); + + +// MAIN // + +/** +* Checks for deflation, computes a shift (`TAU`) and calls DQDS. In case of failure it changes shifts, and tries again until output is positive. +* +* @private +* @param {integer} I0 - first index +* @param {Float64Array} Z - qd array +* @param {integer} strideZ - stride length for `Z` +* @param {NonNegativeInteger} offsetZ - starting index of `Z` +* @param {boolean} IEEE - IEEE arithmetic flag +* @param {Float64Array} out - output array containing `N0`, `PP`, `DMIN`, `SIGMA`, `DESIG`, `QMAX`, `NFAIL`, `ITER`, `NDIV`, `TTYPE`, `DMIN1`, `DMIN2`, `DN`, `DN1`, `DN2`, `G` and `TAU` +* @param {integer} strideOut - stride length for `out` +* @param {NonNegativeInteger} offsetOut - starting index of `out` +* @returns {Float64Array} output array +* +* @example +* var Float64Array = require( '@stdlib/array/float64' ); +* +* var out = new Float64Array( 17 ); +* out[ 0 ] = 3; +* var Z = new Float64Array( [ 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); +* +* dlasq3( 0, Z, 1, 0, true, out, 1, 0 ); +* // out => [ 3.0, 0.0, ~6.022, 0.0, 0.0, 0.0, 0.0, 1.0, 5.0, -1.0, ~15.089, ~35.359, ~6.022, ~15.089, ~35.359, 0.0, 0.0 ] +* // Z => [ 100.0, 181.0, 81.0, ~28.641, 64.0, ~84.359, 49.0, ~20.911, ~36.0, ~40.089, 25.0, ~9.978, 16.0, ~6.022, 9.0, ~28.641 ] +*/ +function dlasq3( I0, Z, strideZ, offsetZ, IEEE, out, strideOut, offsetOut ) { + var goto80; + var desig; + var dmin1; + var dmin2; + var nfail; + var sigma; + var ttype; + var dmin; + var idx1; + var idx2; + var idx3; + var idx4; + var ipn4; + var iter; + var n0in; + var ndiv; + var out1; + var temp; + var qmax; + var dn1; + var dn2; + var idx; + var tau; + var dn; + var j4; + var N0; + var nn; + var PP; + var g; + var s; + var t; + + // Read input values from output array + idx = offsetOut; + N0 = out[ idx ]; + idx += strideOut; + PP = out[ idx ]; + idx += strideOut; + dmin = out[ idx ]; + idx += strideOut; + sigma = out[ idx ]; + idx += strideOut; + desig = out[ idx ]; + idx += strideOut; + qmax = out[ idx ]; + idx += strideOut; + nfail = out[ idx ]; + idx += strideOut; + iter = out[ idx ]; + idx += strideOut; + ndiv = out[ idx ]; + idx += strideOut; + ttype = out[ idx ]; + idx += strideOut; + dmin1 = out[ idx ]; + idx += strideOut; + dmin2 = out[ idx ]; + idx += strideOut; + dn = out[ idx ]; + idx += strideOut; + dn1 = out[ idx ]; + idx += strideOut; + dn2 = out[ idx ]; + idx += strideOut; + g = out[ idx ]; + idx += strideOut; + tau = out[ idx ]; + + n0in = N0; + idx3 = offsetZ + ( strideZ * ( 4 * N0 ) ); + idx4 = offsetZ + ( strideZ * ( 4 * I0 ) ); + + // Check for deflation + while ( true ) { + idx3 = offsetZ + ( strideZ * ( 4 * N0 ) ); + if ( N0 < I0 ) { + out[ offsetOut ] = N0; + out[ offsetOut + strideOut ] = PP; + return out; + } + if ( N0 === I0 ) { + // Deflate 1 eigenvalue + Z[ idx3 ] = Z[ idx3 + ( strideZ * PP ) ] + sigma; + N0 -= 1; + continue; + } + nn = offsetZ + ( strideZ * ( ( 4 * N0 ) + PP + 3 ) ); + if ( N0 !== ( I0 + 1 ) ) { + if ( Z[ nn - ( 5 * strideZ ) ] > TOL2 * ( sigma + Z[ nn - ( 3 * strideZ ) ] ) && Z[ nn - ( strideZ * ( ( 2 * PP ) + 4 ) ) ] > TOL2 * Z[ nn - ( 7 * strideZ ) ] ) { + // Check whether E(N0-1) is negligible, 1 eigenvalue + if ( Z[ nn - ( 9 * strideZ ) ] > TOL2 * sigma && + Z[ nn - ( strideZ * ( ( 2 * PP ) + 8 ) ) ] > TOL2 * Z[ nn - ( 11 * strideZ ) ] ) { + break; // No deflation possible + } + } else { + Z[ idx3 ] = Z[ idx3 + ( strideZ * PP ) ] + sigma; + N0 -= 1; + continue; + } + } + + // Check whether E(N0-2) is negligible, 2 eigenvalues + if ( Z[ nn - ( 3 * strideZ ) ] > Z[ nn - ( 7 * strideZ ) ] ) { + s = Z[ nn - ( 3 * strideZ ) ]; + Z[ nn - ( 3 * strideZ ) ] = Z[ nn - ( 7 * strideZ ) ]; + Z[ nn - ( 7 * strideZ ) ] = s; + } + t = 0.5 * ( ( Z[ nn - ( 7 * strideZ ) ] - Z[ nn - ( 3 * strideZ ) ] ) + Z[ nn - ( 5 * strideZ ) ] ); + if ( Z[ nn - ( 5 * strideZ ) ] > Z[ nn - ( 3 * strideZ ) ] * TOL2 && t !== 0.0 ) { + s = Z[ nn - ( 3 * strideZ ) ] * ( Z[ nn - ( 5 * strideZ ) ] / t ); + if ( s <= t ) { + s = Z[ nn - ( 3 * strideZ ) ] * ( Z[ nn - ( 5 * strideZ ) ] / ( t * ( 1.0 + sqrt( 1.0 + ( s / t ) ) ) ) ); + } else { + s = Z[ nn - ( 3 * strideZ ) ] * ( Z[ nn - ( 5 * strideZ ) ] / ( t + ( sqrt( t ) * sqrt( t + s ) ) ) ); + } + t = Z[ nn - ( 7 * strideZ ) ] + ( s + Z[ nn - ( 5 * strideZ ) ] ); + Z[ nn - ( 3 * strideZ ) ] = Z[ nn - ( 3 * strideZ ) ] * ( Z[ nn - ( 7 * strideZ ) ] / t ); + Z[ nn - ( 7 * strideZ ) ] = t; + } + Z[ idx3 - ( strideZ * 4 ) ] = Z[ nn - ( 7 * strideZ ) ] + sigma; + Z[ idx3 ] = Z[ nn - ( 3 * strideZ ) ] + sigma; + N0 -= 2; + } + + // 50: No deflation, need shift + if ( PP === 2 ) { + PP = 0; + } + + // Reverse the qd-array, if warranted + if ( dmin <= 0.0 || N0 < n0in ) { + if ( CBIAS * Z[ idx4 + ( strideZ * PP ) ] < Z[ idx3 + ( strideZ * PP ) ] ) { + ipn4 = ( 4 * ( I0 + N0 ) ) + 7; + for ( j4 = ( 4 * I0 ) + 3; j4 <= ( 2 * ( I0 + N0 ) ) + 1; j4 += 4 ) { + idx1 = offsetZ + ( strideZ * ( j4 - 3 ) ); + idx2 = offsetZ + ( strideZ * ( ipn4 - j4 - 4 ) ); + + temp = Z[ idx1 ]; + Z[ idx1 ] = Z[ idx2 ]; + Z[ idx2 ] = temp; + + idx1 += strideZ; + idx2 += strideZ; + + temp = Z[ idx1 ]; + Z[ idx1 ] = Z[ idx2 ]; + Z[ idx2 ] = temp; + + idx1 += strideZ; + idx2 -= 3*strideZ; + + temp = Z[ idx1 ]; + Z[ idx1 ] = Z[ idx2 ]; + Z[ idx2 ] = temp; + + idx1 += strideZ; + idx2 += strideZ; + + temp = Z[ idx1 ]; + Z[ idx1 ] = Z[ idx2 ]; + Z[ idx2 ] = temp; + } + if ( ( N0 - I0 ) <= 4 ) { + Z[ idx3 + ( strideZ * ( PP + 2 ) ) ] = Z[ idx4 + ( strideZ * ( PP + 2 ) ) ]; + Z[ idx3 + ( strideZ * ( 3 - PP ) ) ] = Z[ idx4 + ( strideZ * ( 3 - PP ) ) ]; + } + dmin2 = min( dmin2, Z[ idx3 + ( strideZ * ( PP + 2 ) ) ] ); + Z[ idx3 + ( strideZ * ( PP + 2 ) ) ] = min( Z[ idx3 + ( strideZ * ( PP + 2 ) ) ], min( Z[ idx4 + ( strideZ * ( PP + 2 ) ) ], Z[ idx4 + ( strideZ * ( PP + 6 ) ) ] ) ); + Z[ idx3 + ( strideZ * ( 3 - PP ) ) ] = min( Z[ idx3 + ( strideZ * ( 3 - PP ) ) ], min( Z[ idx4 + ( strideZ * ( 3 - PP ) ) ], Z[ idx4 + ( strideZ * ( 7 - PP ) ) ] ) ); + qmax = max( qmax, Z[ idx4 + ( strideZ * ( PP ) ) ], Z[ idx4 + ( strideZ * ( PP + 4 ) ) ] ); + dmin = -0.0; + } + } + + // Choose a shift + out1 = new Float64Array( 3 ); + out1[ 1 ] = ttype; + out1[ 2 ] = g; + dlasq4( I0, N0, Z, strideZ, offsetZ, PP, n0in, dmin, dmin1, dmin2, dn, dn1, dn2, out1, 1, 0 ); + tau = out1[ 0 ]; + ttype = out1[ 1 ]; + g = out1[ 2 ]; + + // Call DQDS until DMIN > 0 + out1 = new Float64Array( 6 ); + while ( true ) { + dlasq5( I0, N0, Z, strideZ, offsetZ, PP, tau, sigma, IEEE, EPS, out1, 1, 0 ); + dmin = out1[ 0 ]; + dmin1 = out1[ 1 ]; + dmin2 = out1[ 2 ]; + dn = out1[ 3 ]; + dn1 = out1[ 4 ]; + dn2 = out1[ 5 ]; + + ndiv += ( N0 - I0 + 2 ); + iter += 1; + + // Check status + if ( dmin >= 0.0 && dmin1 >= 0.0 ) { + // Success + goto80 = false; + break; + } else if ( dmin < 0.0 && dmin1 > 0.0 && Z[ idx3 - ( strideZ * ( PP + 1 ) ) ] < TOL * ( sigma + dn1 ) && abs( dn ) < TOL * sigma ) { + // Convergence hidden by negative DN + Z[ idx3 + ( strideZ * ( 1 - PP ) ) ] = 0.0; + dmin = 0.0; + goto80 = false; + break; + } else if ( dmin < 0.0 ) { + // TAU too big. Select new TAU and try again. + nfail += 1; + if ( ttype < -22 ) { + // Failed twice. Play it safe. + tau = 0.0; + } else if ( dmin1 > 0.0 ) { + // Late failure. Gives excellent shift. + tau = ( tau + dmin ) * ( 1.0 - ( 2.0 * EPS ) ); + ttype -= 11; + } else { + // Early failure. Divide by 4. + tau *= 0.25; + ttype -= 12; + } + continue; + } else if ( isnan( dmin ) ) { + // NaN + if ( tau === 0.0 ) { + goto80 = true; + break; + } + tau = 0.0; + continue; + } + // Possible underflow. Play it safe. + goto80 = true; + break; + } + + // Risk of underflow + if ( goto80 ) { + out1 = new Float64Array( 6 ); + dlasq6( I0, N0, Z, strideZ, offsetZ, PP, out1, 1, 0 ); + dmin = out1[ 0 ]; + dmin1 = out1[ 1 ]; + dmin2 = out1[ 2 ]; + dn = out1[ 3 ]; + dn1 = out1[ 4 ]; + dn2 = out1[ 5 ]; + ndiv += ( N0 - I0 + 2 ); + iter += 1; + tau = 0.0; + } + + if ( tau < sigma ) { + desig += tau; + t = sigma + desig; + desig -= ( t - sigma ); + } else { + t = sigma + tau; + desig = sigma - ( t - tau ) + desig; + } + sigma = t; + + // Store output values + idx = offsetOut; + out[ idx ] = N0; + idx += strideOut; + out[ idx ] = PP; + idx += strideOut; + out[ idx ] = dmin; + idx += strideOut; + out[ idx ] = sigma; + idx += strideOut; + out[ idx ] = desig; + idx += strideOut; + out[ idx ] = qmax; + idx += strideOut; + out[ idx ] = nfail; + idx += strideOut; + out[ idx ] = iter; + idx += strideOut; + out[ idx ] = ndiv; + idx += strideOut; + out[ idx ] = ttype; + idx += strideOut; + out[ idx ] = dmin1; + idx += strideOut; + out[ idx ] = dmin2; + idx += strideOut; + out[ idx ] = dn; + idx += strideOut; + out[ idx ] = dn1; + idx += strideOut; + out[ idx ] = dn2; + idx += strideOut; + out[ idx ] = g; + idx += strideOut; + out[ idx ] = tau; + + return out; +} + + +// EXPORTS // + +module.exports = dlasq3; diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq4.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq4.js new file mode 100644 index 000000000000..53117a0a6373 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq4.js @@ -0,0 +1,401 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 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 max = require( '@stdlib/math/base/special/max' ); +var min = require( '@stdlib/math/base/special/min' ); +var sqrt = require( '@stdlib/math/base/special/sqrt' ); + + +// MAIN // + +/** +* Computes an approximation to the smallest eigenvalue using values of d from the previous transform. +* +* ## Notes +* +* - `Z` is a 1-D array of length >= `4*N0` storing interleaved q/e values. +* - `PP` is `0` for ping, `1` for pong. +* - `TAU` is approximation to the smallest eigenvalue and is used as a shift to accelerate convergence of the dqds iteration. +* - `TTYPE` is an integer flag describing how TAU was computed. +* - `G` is a state variable that is preserved across successive calls and is used to regulate the magnitude of fallback shifts. +* +* @private +* @param {integer} I0 - first index +* @param {integer} N0 - last index +* @param {Float64Array} Z - qd array +* @param {integer} strideZ - stride length for `Z` +* @param {NonNegativeInteger} offsetZ - starting index for `Z` +* @param {integer} PP - ping-pong flag (0 or 1) +* @param {integer} N0IN - value of `N0` at the start of `EIGTEST` +* @param {number} DMIN - minimum value of `d` +* @param {number} DMIN1 - minimum value of `d`, excluding `D(N0)` +* @param {number} DMIN2 - minimum value of `d`, excluding `D(N0)` and `D(N0-1)` +* @param {number} DN - `d(N)` +* @param {number} DN1 - `d(N-1)` +* @param {number} DN2 - `d(N-2)` +* @param {Float64Array} out - output array containing `tau`, `ttype`, and `G` +* @param {integer} strideOut - stride length for `out` +* @param {NonNegativeInteger} offsetOut - starting index for `out` +* @returns {Float64Array} output array +* +* @example +* var Float64Array = require( '@stdlib/array/float64' ); +* +* var out = new Float64Array( [ 0, 0, 0.25 ] ); +* var Z = new Float64Array( [ 5, 0, 7, 0, 9, 0, 3, 0, 10, 0, 4.5, 0, 18, 0, 0, 0 ] ); +* +* dlasq4( 0, 3, Z, 1, 0, 0, 3, 0.2, 0.15, 0.1, 0.8, 0.7, 0.6, out, 1, 0 ); +* // out => [ ~0.05, -6, 0.25 ] +*/ +function dlasq4( I0, N0, Z, strideZ, offsetZ, PP, N0IN, DMIN, DMIN1, DMIN2, DN, DN1, DN2, out, strideOut, offsetOut ) { // eslint-disable-line max-len, max-params + var CNST1; + var CNST2; + var CNST3; + var ttype; + var idx2; + var gap1; + var gap2; + var gam; + var idx; + var tau; + var a2; + var b1; + var b2; + var i4; + var nn; + var np; + var G; + var s; + + idx2 = offsetOut; + tau = out[ idx2 ]; + idx2 += strideOut; + ttype = out[ idx2 ]; + idx2 += strideOut; + G = out[ idx2 ]; + + CNST1 = 0.563; + CNST2 = 1.010; + CNST3 = 1.050; + + // A negative DMIN forces the shift to take that absolute value TTYPE records the type of shift. + if ( DMIN <= 0 ) { + tau = -DMIN; + ttype = -1; + + idx2 = offsetOut; + out[ idx2 ] = tau; + idx2 += strideOut; + out[ idx2 ] = ttype; + idx2 += strideOut; + out[ idx2 ] = G; + return out; + } + + nn = offsetZ + ( strideZ*( ( 4*N0 ) + PP + 3 ) ); + + if ( N0IN === N0 ) { + // No eigenvalues deflated. + if ( DMIN === DN || DMIN === DN1 ) { + b1 = sqrt( Z[ nn - ( 3*strideZ ) ] )*sqrt( Z[ nn - ( 5*strideZ ) ] ); + b2 = sqrt( Z[ nn - ( 7*strideZ ) ] )*sqrt( Z[ nn - ( 9*strideZ ) ] ); + a2 = Z[ nn - ( 7*strideZ ) ] + Z[ nn - ( 5*strideZ ) ]; + + // Cases 2 and 3. + if ( DMIN === DN && DMIN1 === DN1 ) { + gap2 = DMIN2 - a2 - ( DMIN2*0.25 ); + + if ( gap2 > 0 && gap2 > b2 ) { + gap1 = a2 - DN - ( ( b2 / gap2 )*b2 ); + } else { + gap1 = a2 - DN - ( b1 + b2 ); + } + if ( gap1 > 0 && gap1 > b1 ) { + s = max( DN - ( ( b1 / gap1 )*b1 ), 0.5*DMIN ); + ttype = -2; + } else { + s = 0; + if ( DN > b1 ) { + s = DN - b1; + } + if ( a2 > ( b1 + b2 ) ) { + s = min( s, a2 - ( b1 + b2 ) ); + } + s = max( s, 0.333*DMIN ); + ttype = -3; + } + } else { + // Case 4. + ttype = -4; + s = 0.25*DMIN; + + if ( DMIN === DN ) { + gam = DN; + a2 = 0; + if ( Z[ nn - ( 5*strideZ ) ] > Z[ nn - ( 7*strideZ ) ] ) { + idx2 = offsetOut + strideOut; + out[ idx2 ] = ttype; + idx2 += strideOut; + out[ idx2 ] = G; + return out; + } + b2 = Z[ nn - ( 5*strideZ ) ] / Z[ nn - ( 7*strideZ ) ]; + np = nn - ( 9*strideZ ); + } else { + np = nn - ( 2*PP*strideZ ); + gam = DN1; + if ( Z[ np - ( 4*strideZ ) ] > Z[ np - ( 2*strideZ ) ] ) { + idx2 = offsetOut + strideOut; + out[ idx2 ] = ttype; + idx2 += strideOut; + out[ idx2 ] = G; + return out; + } + a2 = Z[ np - ( 4*strideZ ) ] / Z[ np - ( 2*strideZ ) ]; + if ( Z[ nn - ( 9*strideZ ) ] > Z[ nn - ( 11*strideZ ) ] ) { + idx2 = offsetOut + strideOut; + out[ idx2 ] = ttype; + idx2 += strideOut; + out[ idx2 ] = G; + return out; + } + b2 = Z[ nn - ( 9*strideZ ) ] / Z[ nn - ( 11*strideZ ) ]; + np = nn - ( 13*strideZ ); + } + + // Approximate contribution to norm squared from I < NN-1. + a2 += b2; + i4 = np; + for ( idx = ( np - offsetZ ) / strideZ; idx >= ( 4*I0 ) + PP + 2; idx -= 4 ) { + if ( b2 === 0 ) { + break; + } + b1 = b2; + if ( Z[ i4 ] > Z[ i4 - ( 2*strideZ ) ] ) { + idx2 = offsetOut + strideOut; + out[ idx2 ] = ttype; + idx2 += strideOut; + out[ idx2 ] = G; + return out; + } + b2 *= Z[ i4 ] / Z[ i4 - ( 2*strideZ ) ]; + a2 += b2; + if ( ( 100*max( b2, b1 ) ) < a2 || CNST1 < a2 ) { + break; + } + i4 -= 4*strideZ; + } + a2 *= CNST3; + + // Rayleigh quotient residual bound. + if ( a2 < CNST1 ) { + s = gam*( 1 - sqrt( a2 ) ) / ( 1 + a2 ); + } + } + } else if ( DMIN === DN2 ) { + // Case 5. + ttype = -5; + s = 0.25*DMIN; + + // Compute contribution to norm squared from I > NN-2. + np = nn - ( 2*PP*strideZ ); + b1 = Z[ np - ( 2*strideZ ) ]; + b2 = Z[ np - ( 6*strideZ ) ]; + gam = DN2; + if ( Z[ np - ( 8*strideZ ) ] > b2 || Z[ np - ( 4*strideZ ) ] > b1 ) { + idx2 = offsetOut + strideOut; + out[ idx2 ] = ttype; + idx2 += strideOut; + out[ idx2 ] = G; + return out; + } + a2 = ( Z[ np - ( 8*strideZ ) ] / b2 )*( 1 + ( Z[ np - ( 4*strideZ ) ] / b1 ) ); + + // Approximate contribution to norm squared from I < NN-2. + if ( ( N0 - I0 ) > 2 ) { + b2 = Z[ nn - ( 13*strideZ ) ] / Z[ nn - ( 15*strideZ ) ]; + a2 += b2; + + i4 = nn - ( strideZ*17 ); + for ( idx = ( ( nn - offsetZ ) / strideZ ) - 17; idx >= ( 4*I0 ) + PP + 2; idx -= 4 ) { + if ( b2 === 0 ) { + break; + } + b1 = b2; + if ( Z[ i4 ] > Z[ i4 - ( 2*strideZ ) ] ) { + idx2 = offsetOut + strideOut; + out[ idx2 ] = ttype; + idx2 += strideOut; + out[ idx2 ] = G; + return out; + } + b2 *= ( Z[ i4 ] / Z[ i4 - ( 2*strideZ ) ] ); + a2 += b2; + if ( ( 100*max( b2, b1 ) ) < a2 || CNST1 < a2 ) { + break; + } + i4 -= 4*strideZ; + } + a2 *= CNST3; + } + + if ( a2 < CNST1 ) { + s = gam * ( 1 - sqrt( a2 ) ) / ( 1 + a2 ); + } + } else { + // Case 6, no information to guide us. + if ( ttype === -6 ) { // Case when `ttype` is previously set by another routine. + G += 0.333*( 1 - G ); + } else if ( ttype === -18 ) { // Case when `ttype` is previously set by another routine. + G += 0.25*0.333; + } else { + G = 0.25; + } + s = G*DMIN; + ttype = -6; + } + } else if ( N0IN === ( N0 + 1 ) ) { + // 1 eigenvalue just deflated. Use DMIN1, DN1 for DMIN and DN. + if ( DMIN1 === DN1 && DMIN2 === DN2 ) { + // Cases 7 and 8. + ttype = -7; + s = 0.3330*DMIN1; + + if ( Z[ nn - ( 5*strideZ ) ] > Z[ nn - ( 7*strideZ ) ] ) { + idx2 = offsetOut + strideOut; + out[ idx2 ] = ttype; + idx2 += strideOut; + out[ idx2 ] = G; + return out; + } + + b1 = Z[ nn - ( 5*strideZ ) ] / Z[ nn - ( 7*strideZ ) ]; + b2 = b1; + + if ( b2 !== 0 ) { + i4 = offsetZ + ( strideZ*( ( 4*N0 ) - 6 + PP ) ); + for ( idx = 2; idx >= 0; idx-- ) { + a2 = b1; + if ( Z[ i4 ] > Z[ i4 - ( 2*strideZ ) ] ) { + idx2 = offsetOut + strideOut; + out[ idx2 ] = ttype; + idx2 += strideOut; + out[ idx2 ] = G; + return out; + } + b1 *= Z[ i4 ] / Z[ i4 - ( 2*strideZ ) ]; + b2 += b1; + + if ( ( 100*max( b1, a2 ) ) < b2 ) { + break; + } + i4 -= 4*strideZ; + } + } + + b2 = sqrt( CNST3*b2 ); + a2 = DMIN1 / ( 1 + ( b2*b2 ) ); + gap2 = ( 0.5*DMIN2 ) - a2; + + if ( gap2 > 0 && gap2 > ( b2*a2 ) ) { + s = max( s, a2*( 1 - ( CNST2*a2*( b2 / gap2 )*b2 ) ) ); + } else { + s = max( s, a2*( 1 - ( CNST2*b2 ) ) ); + ttype = -8; + } + } else { + // Case 9. + s = 0.25*DMIN1; + if ( DMIN1 === DN1 ) { + s = 0.5*DMIN1; + } + ttype = -9; + } + } else if ( N0IN === ( N0 + 2 ) ) { + // 2 eigenvalues deflated. Use DMIN2, DN2 for DMIN and DN. + + // Cases 10 and 11. + if ( DMIN2 === DN2 && ( 2*Z[ nn - ( 5*strideZ ) ] ) < Z[ nn - ( 7*strideZ ) ] ) { + ttype = -10; + s = 0.3330*DMIN2; + if ( Z[ nn - ( 5*strideZ ) ] > Z[ nn - ( 7*strideZ ) ] ) { + idx2 = offsetOut + strideOut; + out[ idx2 ] = ttype; + idx2 += strideOut; + out[ idx2 ] = G; + return out; + } + b1 = Z[ nn - ( 5*strideZ ) ] / Z[ nn - ( 7*strideZ ) ]; + b2 = b1; + if ( b2 !== 0 ) { + i4 = offsetZ + ( strideZ*( ( 4*N0 ) - 6 + PP ) ); + for ( idx = 1; idx >= 0; idx-- ) { + if ( Z[ i4 ] > Z[ i4 - ( 2*strideZ ) ] ) { + idx2 = offsetOut + strideOut; + out[ idx2 ] = ttype; + idx2 += strideOut; + out[ idx2 ] = G; + return out; + } + b1 *= ( Z[ i4 ] / Z[ i4 - ( 2*strideZ ) ] ); + b2 += b1; + if ( ( 100*b1 ) < b2 ) { + break; + } + i4 -= 4*strideZ; + } + } + + b2 = sqrt( CNST3*b2 ); + a2 = DMIN2 / ( 1 + ( b2*b2 ) ); + gap2 = Z[ nn - ( 7*strideZ ) ] + Z[ nn - ( 9*strideZ ) ] - ( sqrt( Z[ nn - ( 11*strideZ ) ] )*sqrt( Z[ nn - ( 9*strideZ ) ] ) ) - a2; + if ( gap2 > 0 && gap2 > ( b2*a2 ) ) { + s = max( s, a2*( 1 - ( CNST2*a2*( b2 / gap2 )*b2 ) ) ); + } else { + s = max( s, a2*( 1 - ( CNST2*b2 ) ) ); + } + } else { + // Case 11. + s = 0.25*DMIN2; + ttype = -11; + } + } else if ( N0IN > ( N0 + 2 ) ) { + // Case 12, more than 2 eigenvalues deflated. No information. + s = 0; + ttype = -12; + } + + tau = s; + + idx2 = offsetOut; + out[ idx2 ] = tau; + idx2 += strideOut; + out[ idx2 ] = ttype; + idx2 += strideOut; + out[ idx2 ] = G; + return out; +} + + +// EXPORTS // + +module.exports = dlasq4; diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq5.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq5.js new file mode 100644 index 000000000000..8c938fe35377 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq5.js @@ -0,0 +1,331 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 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 min = require( '@stdlib/math/base/special/min' ); + + +// MAIN // + +/** +* Computes one dqds transform in ping-pong form. +* +* ## Notes +* +* - `Z` is a 1-D array of length >= `4*N0` storing interleaved q/e values. +* - `PP` is `0` for ping, `1` for pong. +* +* @private +* @param {integer} I0 - the first index +* @param {integer} N0 - the last index +* @param {Float64Array} Z - the QD array +* @param {integer} strideZ - stride length for `Z` +* @param {NonNegativeInteger} offsetZ - starting index for `Z` +* @param {integer} PP - ping-pong flag (0 or 1) +* @param {number} TAU - the shift +* @param {number} SIGMA - the accumulated shift +* @param {boolean} IEEE - IEEE arithmetic flag +* @param {number} EPS - epsilon used by the routine +* @param {Float64Array} out - output array containing `DMIN`, `DMIN1`, `DMIN2`, `DN`, `DNM1`, and `DNM2` respectively +* @param {integer} strideOut - stride length for `out` +* @param {NonNegativeInteger} offsetOut - starting index of `out` +* @returns {Float64Array} output array +* +* @example +* var Float64Array = require( '@stdlib/array/float64' ); +* +* var out = new Float64Array( 6 ); +* var Z = new Float64Array( [ 5, 0, 7, 0, 9, 0, 3, 0, 11, 0, 4, 0, 20, 0, 0, 0 ] ); +* +* dlasq5( 0, 3, Z, 1, 0, 0, 0.1, 0.0, true, 2.220446049250313e-16, out, 1, 0 ); +* // out => [ ~3.606, ~3.606, ~3.606, ~11.823, ~5.904, ~3.606 ] +*/ +function dlasq5( I0, N0, Z, strideZ, offsetZ, PP, TAU, SIGMA, IEEE, EPS, out, strideOut, offsetOut ) { // eslint-disable-line max-len, max-params + var dthresh; + var dmin1; + var dmin2; + var dnm1; + var dnm2; + var emin; + var dmin; + var temp; + var j4p2; + var idx; + var j4; + var dn; + var d; + + // Quick return... + if ( ( N0 - I0 - 1 ) <= 0 ) { + return out; + } + + dthresh = EPS * ( SIGMA + TAU ); + if ( TAU < ( dthresh * 0.5 ) ) { + TAU = 0; + } + + if ( TAU > 0 || TAU < 0 ) { + j4 = offsetZ + ( strideZ*( ( 4*I0 ) + PP ) ); + emin = Z[ j4 + ( 4*strideZ ) ]; + d = Z[ j4 ] - TAU; + dmin = d; + dmin1 = -Z[ j4 ]; + + if ( IEEE ) { + if ( PP === 0 ) { + j4 = offsetZ + ( strideZ * ( ( 4*I0 ) + 3 ) ); + for ( idx = ( 4*I0 ) + 3; idx <= ( 4*N0 ) - 9; idx += 4 ) { + Z[ j4 - ( 2*strideZ ) ] = d + Z[ j4 - strideZ ]; + temp = Z[ j4 + strideZ ] / Z[ j4 - ( 2*strideZ ) ]; + d = ( d*temp ) - TAU; + dmin = min( dmin, d ); + Z[ j4 ] = Z[ j4 - strideZ ]*temp; + emin = min( Z[ j4 ], emin ); + j4 += 4*strideZ; + } + } else { + j4 = offsetZ + ( strideZ * ( ( 4*I0 ) + 3 ) ); + for ( idx = ( 4*I0 ) + 3; idx <= ( 4*N0 ) - 9; idx += 4 ) { + Z[ j4 - ( 3*strideZ ) ] = d + Z[ j4 ]; + temp = Z[ j4 + ( 2*strideZ ) ] / Z[ j4 - ( 3*strideZ ) ]; + d = ( d*temp ) - TAU; + dmin = min( dmin, d ); + Z[ j4 - strideZ ] = Z[ j4 ] * temp; + emin = min( Z[ j4 - strideZ ], emin ); + j4 += 4*strideZ; + } + } + + // Unroll last two steps. + dnm2 = d; + dmin2 = dmin; + + j4 = offsetZ + ( strideZ * ( ( 4*N0 ) - 5 - PP ) ); + j4p2 = j4 + ( strideZ * ( ( 2*PP ) - 1 ) ); + Z[ j4 - ( 2*strideZ ) ] = dnm2 + Z[ j4p2 ]; + Z[ j4 ] = Z[ j4p2 + ( 2*strideZ ) ] * ( Z[ j4p2 ] / Z[ j4 - ( 2*strideZ ) ] ); + dnm1 = ( Z[ j4p2 + ( 2*strideZ ) ]*( dnm2 / Z[ j4 - ( 2*strideZ ) ] ) ) - TAU; + dmin = min( dmin, dnm1 ); + + dmin1 = dmin; + j4 += 4*strideZ; + j4p2 = j4 + ( strideZ * ( ( 2*PP ) - 1 ) ); + Z[ j4 - ( 2*strideZ ) ] = dnm1 + Z[ j4p2 ]; + Z[ j4 ] = Z[ j4p2 + ( 2*strideZ ) ] * ( Z[ j4p2 ] / Z[ j4 - ( 2*strideZ ) ] ); + dn = ( Z[ j4p2 + ( 2*strideZ ) ]*( dnm1 / Z[ j4 - ( 2*strideZ ) ] ) ) - TAU; + dmin = min( dmin, dn ); + } else { + if ( PP === 0 ) { + j4 = offsetZ + ( strideZ * ( ( 4*I0 ) + 3 ) ); + for ( idx = ( 4*I0 ) + 3; idx <= ( 4*N0 ) - 9; idx += 4 ) { + Z[ j4 - ( 2*strideZ ) ] = d + Z[ j4 - strideZ ]; + if ( d < 0 ) { + return out; + } + Z[ j4 ] = Z[ j4 + strideZ ] * ( Z[ j4 - strideZ ] / Z[ j4 - ( 2*strideZ ) ] ); + d = ( Z[ j4 + strideZ ]*( d / Z[ j4 - ( 2*strideZ ) ] ) ) - TAU; + dmin = min( dmin, d ); + emin = min( emin, Z[ j4 ] ); + j4 += 4*strideZ; + } + } else { + j4 = offsetZ + ( strideZ * ( ( 4*I0 ) + 3 ) ); + for ( idx = ( 4*I0 ) + 3; idx <= ( 4*N0 ) - 9; idx += 4 ) { + Z[ j4 - ( 3*strideZ ) ] = d + Z[ j4 ]; + if ( d < 0 ) { + return out; + } + Z[ j4 - strideZ ] = Z[ j4 + ( 2*strideZ ) ] * ( Z[ j4 ] / Z[ j4 - ( 3*strideZ ) ] ); + d = ( Z[ j4 + ( 2*strideZ ) ]*( d / Z[ j4 - ( 3*strideZ ) ] ) ) - TAU; + dmin = min( dmin, d ); + emin = min( emin, Z[ j4 - strideZ ] ); + j4 += 4*strideZ; + } + } + + // Unroll last two steps. + dnm2 = d; + dmin2 = dmin; + + j4 = offsetZ + ( strideZ * ( ( 4*N0 ) - 5 - PP ) ); + j4p2 = j4 + ( strideZ * ( ( 2*PP ) - 1 ) ); + Z[ j4 - ( 2*strideZ ) ] = dnm2 + Z[ j4p2 ]; + if ( dnm2 < 0 ) { + return out; + } + Z[ j4 ] = Z[ j4p2 + ( 2*strideZ ) ] * ( Z[ j4p2 ] / Z[ j4 - ( 2*strideZ ) ] ); + dnm1 = ( Z[ j4p2 + ( 2*strideZ ) ]*( dnm2 / Z[ j4 - ( 2*strideZ ) ] ) ) - TAU; + dmin = min( dmin, dnm1 ); + + dmin1 = dmin; + j4 += 4*strideZ; + j4p2 = j4 + ( strideZ * ( ( 2*PP ) - 1 ) ); + Z[ j4 - ( 2*strideZ ) ] = dnm1 + Z[ j4p2 ]; + if ( dnm1 < 0 ) { + return out; + } + Z[ j4 ] = Z[ j4p2 + ( 2*strideZ ) ] * ( Z[ j4p2 ] / Z[ j4 - ( 2*strideZ ) ] ); + dn = ( Z[ j4p2 + ( 2*strideZ ) ]*( dnm1 / Z[ j4 - ( 2*strideZ ) ] ) ) - TAU; + dmin = min( dmin, dn ); + } + } else { + // This is the version that sets d's to zero if they are small enough + j4 = offsetZ + ( strideZ * ( ( 4*I0 ) + PP ) ); + emin = Z[ j4 + ( 4*strideZ ) ]; + d = Z[ j4 ]; + dmin = d; + dmin1 = -Z[ j4 ]; + + if ( IEEE ) { + if ( PP === 0 ) { + j4 = offsetZ + ( strideZ * ( ( 4*I0 ) + 3 ) ); + for ( idx = ( 4*I0 ) + 3; idx <= ( 4*N0 ) - 9; idx += 4 ) { + Z[ j4 - ( 2*strideZ ) ] = d + Z[ j4 - strideZ ]; + temp = Z[ j4 + strideZ ] / Z[ j4 - ( 2*strideZ ) ]; + d = ( d*temp ) - TAU; + if ( d < dthresh ) { + d = 0; + } + dmin = min( dmin, d ); + Z[ j4 ] = Z[ j4 - strideZ ] * temp; + emin = min( Z[ j4 ], emin ); + j4 += 4*strideZ; + } + } else { + j4 = offsetZ + ( strideZ * ( ( 4*I0 ) + 3 ) ); + for ( idx = ( 4*I0 ) + 3; idx <= ( 4*N0 ) - 9; idx += 4 ) { + Z[ j4 - ( 3*strideZ ) ] = d + Z[ j4 ]; + temp = Z[ j4 + ( 2*strideZ ) ] / Z[ j4 - ( 3*strideZ ) ]; + d = ( d*temp ) - TAU; + if ( d < dthresh ) { + d = 0; + } + dmin = min( dmin, d ); + Z[ j4 - strideZ ] = Z[ j4 ] * temp; + emin = min( Z[ j4 - strideZ ], emin ); + j4 += 4*strideZ; + } + } + + // Unroll last two steps. + dnm2 = d; + dmin2 = dmin; + + j4 = offsetZ + ( strideZ * ( ( 4*N0 ) - 5 - PP ) ); + j4p2 = j4 + ( strideZ * ( ( 2*PP ) - 1 ) ); + Z[ j4 - ( 2*strideZ ) ] = dnm2 + Z[ j4p2 ]; + Z[ j4 ] = Z[ j4p2 + ( 2*strideZ ) ] * ( Z[ j4p2 ] / Z[ j4 - ( 2*strideZ ) ] ); + dnm1 = ( Z[ j4p2 + ( 2*strideZ ) ]*( dnm2 / Z[ j4 - ( 2*strideZ ) ] ) ) - TAU; + dmin = min( dmin, dnm1 ); + + dmin1 = dmin; + j4 += 4*strideZ; + j4p2 = j4 + ( strideZ * ( ( 2*PP ) - 1 ) ); + Z[ j4 - ( 2*strideZ ) ] = dnm1 + Z[ j4p2 ]; + Z[ j4 ] = Z[ j4p2 + ( 2*strideZ ) ] * ( Z[ j4p2 ] / Z[ j4 - ( 2*strideZ ) ] ); + dn = ( Z[ j4p2 + ( 2*strideZ ) ]*( dnm1 / Z[ j4 - ( 2*strideZ ) ] ) ) - TAU; + dmin = min( dmin, dn ); + } else { + if ( PP === 0 ) { + j4 = offsetZ + ( strideZ * ( ( 4*I0 ) + 3 ) ); + for ( idx = ( 4*I0 ) + 3; idx <= ( 4*N0 ) - 9; idx += 4 ) { + Z[ j4 - ( 2*strideZ ) ] = d + Z[ j4 - strideZ ]; + if ( d < 0 ) { + return out; + } + Z[ j4 ] = Z[ j4 + strideZ ] * ( Z[ j4 - strideZ ] / Z[ j4 - ( 2*strideZ ) ] ); + d = ( Z[ j4 + strideZ ]*( d / Z[ j4 - ( 2*strideZ ) ] ) ) - TAU; + if ( d < dthresh ) { + d = 0; + } + dmin = min( dmin, d ); + emin = min( emin, Z[ j4 ] ); + j4 += 4*strideZ; + } + } else { + j4 = offsetZ + ( strideZ * ( ( 4*I0 ) + 3 ) ); + for ( idx = ( 4*I0 ) + 3; idx <= ( 4*N0 ) - 9; idx += 4 ) { + Z[ j4 - ( 3*strideZ ) ] = d + Z[ j4 ]; + if ( d < 0 ) { + return out; + } + Z[ j4 - strideZ ] = Z[ j4 + ( 2*strideZ ) ] * ( Z[ j4 ] / Z[ j4 - ( 3*strideZ ) ] ); + d = ( Z[ j4 + ( 2*strideZ ) ]*( d / Z[ j4 - ( 3*strideZ ) ] ) ) - TAU; + if ( d < dthresh ) { + d = 0; + } + dmin = min( dmin, d ); + emin = min( emin, Z[ j4 - strideZ ] ); + j4 += 4*strideZ; + } + } + + // Unroll last two steps. + dnm2 = d; + dmin2 = dmin; + + j4 = offsetZ + ( strideZ * ( ( 4*N0 ) - 5 - PP ) ); + j4p2 = j4 + ( strideZ * ( ( 2*PP ) - 1 ) ); + Z[ j4 - ( 2*strideZ ) ] = dnm2 + Z[ j4p2 ]; + if ( dnm2 < 0 ) { + return out; + } + Z[ j4 ] = Z[ j4p2 + ( 2*strideZ ) ] * ( Z[ j4p2 ] / Z[ j4 - ( 2*strideZ ) ] ); + dnm1 = ( Z[ j4p2 + ( 2*strideZ ) ]*( dnm2 / Z[ j4 - ( 2*strideZ ) ] ) ) - TAU; + + dmin = min( dmin, dnm1 ); + dmin1 = dmin; + j4 += 4*strideZ; + j4p2 = j4 + ( strideZ * ( ( 2*PP ) - 1 ) ); + Z[ j4 - ( 2*strideZ ) ] = dnm1 + Z[ j4p2 ]; + if ( dnm1 < 0 ) { + return out; + } + Z[ j4 ] = Z[ j4p2 + ( 2*strideZ ) ] * ( Z[ j4p2 ] / Z[ j4 - ( 2*strideZ ) ] ); + dn = ( Z[ j4p2 + ( 2*strideZ ) ]*( dnm1 / Z[ j4 - ( 2*strideZ ) ] ) ) - TAU; + dmin = min( dmin, dn ); + } + } + + Z[ j4 + ( 2*strideZ ) ] = dn; + Z[ offsetZ + ( strideZ * ( ( 4*N0 ) + 3 - PP ) ) ] = emin; + + idx = offsetOut; + out[ idx ] = dmin; + idx += strideOut; + out[ idx ] = dmin1; + idx += strideOut; + out[ idx ] = dmin2; + idx += strideOut; + out[ idx ] = dn; + idx += strideOut; + out[ idx ] = dnm1; + idx += strideOut; + out[ idx ] = dnm2; + + return out; +} + + +// EXPORTS // + +module.exports = dlasq5; diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq6.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq6.js new file mode 100644 index 000000000000..8c925518502d --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasq6.js @@ -0,0 +1,220 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 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 dlamch = require( '@stdlib/lapack/base/dlamch' ); +var min = require( '@stdlib/math/base/special/min' ); + + +// VARIABLES // + +var SFMIN = dlamch( 'S' ); + + +// MAIN // + +/** +* Computes one dqd transform in ping-pong form. +* +* ## Notes +* +* - Z is a 1-D array of length >= 4*N0 storing interleaved q/e values. +* - PP is 0 for ping, 1 for pong. +* +* @private +* @param {integer} I0 - the first index +* @param {integer} N0 - the last index +* @param {Float64Array} Z - the QD array +* @param {integer} strideZ - stride length for `z` +* @param {NonNegativeInteger} offsetZ - starting index for `z` +* @param {boolean} PP - ping-pong flag (0 or 1) +* @param {Float64Array} out - output array containing `DMIN`, `DMIN1`, `DMIN2`, `DN`, `DNM1`, and `DNM2` respectively +* @param {integer} strideOut - stride length for `out` +* @param {NonNegativeInteger} offsetOut - starting index of `out` +* @returns {Float64Array} output array +* +* @example +* var Float64Array = require( '@stdlib/array/float64' ); +* +* var out = new Float64Array( 6 ); +* var Z = new Float64Array( [ 5, 0, 7, 0, 9, 0, 3, 0, 11, 0, 4, 0, 20, 0, 0, 0 ] ); +* +* dlasq6( 0, 3, Z, 1, 0, 0, out, 1, 0 ); +* // out => [ 3.75, 3.75, 3.75, ~12.088, ~6.111, 3.75 ] +*/ +function dlasq6( I0, N0, Z, strideZ, offsetZ, PP, out, strideOut, offsetOut ) { + var dmin1; + var dmin2; + var dnm1; + var dnm2; + var emin; + var dmin; + var temp; + var j4p2; + var idx; + var iz1; + var iz2; + var iz3; + var iz4; + var j4; + var dn; + var d; + + // Quick return... + if ( ( N0 - I0 - 1 ) <= 0 ) { + return out; + } + + idx = offsetZ; + + j4 = offsetZ + ( strideZ*( ( 4*I0 ) + PP ) ); + emin = Z[ j4 + ( 4*strideZ ) ]; + d = Z[ j4 ]; + dmin = d; + + j4 = offsetZ + ( strideZ*( ( 4*I0 ) + 3 ) ); + + if ( PP === 0 ) { + for ( idx = ( 4*I0 ) + 3; idx <= ( 4*N0 ) - 9; idx += 4 ) { + iz1 = j4 + strideZ; + iz2 = j4 - strideZ; + iz3 = j4 - ( 2*strideZ ); + + Z[ iz3 ] = d + Z[ iz2 ]; + if ( Z[ iz3 ] === 0 ) { + Z[ j4 ] = 0; + d = Z[ iz1 ]; + dmin = d; + emin = 0; + } else if ( ( SFMIN*Z[ iz1 ] < Z[ iz3 ] ) && + ( SFMIN*Z[ iz3 ] < Z[ iz1 ] ) ) { + temp = Z[ iz1 ] / Z[ iz3 ]; + Z[ j4 ] = Z[ iz2 ]*temp; + d *= temp; + } else { + Z[ j4 ] = Z[ iz1 ]*( Z[ iz2 ] / Z[ iz3 ] ); + d = Z[ iz1 ]*( d / Z[ iz3 ] ); + } + dmin = min( dmin, d ); + emin = min( emin, Z[ j4 ] ); + j4 += 4*strideZ; + } + } else { + for ( idx = ( 4*I0 ) + 3; idx <= ( 4*N0 ) - 9; idx += 4 ) { + iz1 = j4 + ( 2*strideZ ); + iz2 = j4 - strideZ; + iz3 = j4 - ( 2*strideZ ); + iz4 = j4 - ( 3*strideZ ); + + Z[ iz4 ] = d + Z[ j4 ]; + if ( Z[ iz4 ] === 0 ) { + Z[ iz2 ] = 0; + d = Z[ iz1 ]; + dmin = d; + emin = 0; + } else if ( ( SFMIN*Z[ iz1 ] < Z[ iz4 ] ) && + ( SFMIN*Z[ iz4 ] < Z[ iz1 ] ) ) { + temp = Z[ iz1 ] / Z[ iz4 ]; + Z[ j4 - strideZ ] = Z[ j4 ]*temp; + d *= temp; + } else { + Z[ iz2 ] = Z[ iz1 ]*( Z[ j4 ] / Z[ iz4 ] ); + d = Z[ iz1 ]*( d / Z[ iz4 ] ); + } + dmin = min( dmin, d ); + emin = min( emin, Z[ iz2 ] ); + j4 += 4*strideZ; + } + } + + // Unroll last two steps. + dnm2 = d; + dmin2 = dmin; + j4 = offsetZ + ( strideZ*( ( 4*N0 ) - 5 - PP ) ); + j4p2 = j4 + ( strideZ*( ( 2*PP ) - 1 ) ); + + iz1 = j4 - ( 2*strideZ ); + iz2 = j4p2 + ( 2*strideZ ); + + Z[ j4 - ( 2*strideZ ) ] = dnm2 + Z[ j4p2 ]; + if ( Z[ iz1 ] === 0 ) { + Z[ j4 ] = 0; + dnm1 = Z[ iz2 ]; + dmin = dnm1; + emin = 0; + } else if ( ( SFMIN*Z[ iz2 ] < Z[ iz1 ] ) && + ( SFMIN*Z[ iz1 ] < Z[ iz2 ] ) ) { + temp = Z[ iz2 ] / Z[ iz1 ]; + Z[ j4 ] = Z[ j4p2 ]*temp; + dnm1 = dnm2*temp; + } else { + Z[ j4 ] = Z[ iz2 ]*( Z[ j4p2 ] / Z[ iz1 ] ); + dnm1 = Z[ iz2 ]*( dnm2 / Z[ iz1 ] ); + } + dmin = min( dmin, dnm1 ); + + dmin1 = dmin; + j4 += 4*strideZ; + j4p2 = j4 + ( strideZ*( ( 2*PP ) - 1 ) ); + + iz1 = j4 - ( 2*strideZ ); + iz2 = j4p2 + ( 2*strideZ ); + + Z[ iz1 ] = dnm1 + Z[ j4p2 ]; + if ( Z[ iz1 ] === 0 ) { + Z[ j4 ] = 0; + dn = Z[ iz2 ]; + dmin = dn; + emin = 0; + } else if ( ( SFMIN*Z[ iz2 ] < Z[ iz1 ] ) && + ( SFMIN*Z[ iz1 ] < Z[ iz2 ] ) ) { + temp = Z[ iz2 ] / Z[ iz1 ]; + Z[ j4 ] = Z[ j4p2 ]*temp; + dn = dnm1*temp; + } else { + Z[ j4 ] = Z[ iz2 ]*( Z[ j4p2 ] / Z[ iz1 ] ); + dn = Z[ iz2 ]*( dnm1 / Z[ iz1 ] ); + } + dmin = min( dmin, dn ); + + Z[ j4 + ( 2*strideZ ) ] = dn; + Z[ offsetZ + ( strideZ*( ( 4*( N0 ) ) + 3 - PP ) ) ] = emin; + + idx = offsetOut; + out[ idx ] = dmin; + idx += strideOut; + out[ idx ] = dmin1; + idx += strideOut; + out[ idx ] = dmin2; + idx += strideOut; + out[ idx ] = dn; + idx += strideOut; + out[ idx ] = dnm1; + idx += strideOut; + out[ idx ] = dnm2; + + return out; +} + + +// EXPORTS // + +module.exports = dlasq6; diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasrt.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasrt.js new file mode 100644 index 000000000000..01ad9afb80a9 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/dlasrt.js @@ -0,0 +1,250 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 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 lowercase = require( '@stdlib/string/base/lowercase' ); +var floor = require( '@stdlib/math/base/special/floor' ); +var Float64Array = require( '@stdlib/array/float64' ); + + +// VARIABLES // + +var SELECT = 20; + + +// MAIN // + +/** +* Sort an array of doubles in increasing or decreasing order using quicksort, with insertion sort for small partitions (size <= 20). +* +* @private +* @param {string} ID - sort direction: 'I' or 'D' +* @param {NonNegativeInteger} N - number of elements to sort +* @param {Float64Array} D - array to sort in-place +* @param {integer} strideD - stride length for `d` +* @param {NonNegativeInteger} offsetD - starting index for `d` +* @returns {integer} status code (0 = success, -i = the i-th argument had an illegal value) +* +* @example +* var Float64Array = require( '@stdlib/array/float64' ); +* +* var D = new Float64Array( [ 5, 7, 9, 3, 11 ] ); +* +* dlasrt( 'I', 5, D, 1, 0 ); +* // D => [ 3, 5, 7, 9, 11 ] +*/ +function dlasrt( ID, N, D, strideD, offsetD ) { + var stkpnt; + var stack; + var start; + var dmnmx; + var endd; + var idx1; + var idx2; + var idx3; + var dir; + var tmp; + var d1; + var d2; + var d3; + var i; + var j; + + // Determine sort direction + dir = -1; + if ( lowercase( ID ) === 'd' ) { + dir = 0; + } else if ( lowercase( ID ) === 'i' ) { + dir = 1; + } + + if ( dir === -1 ) { + return -1; + } + if ( N < 0 ) { + return -2; + } + + // Quick return + if ( N <= 1 ) { + return 0; + } + + // Initialize the stack with the full range (using 0-based indices)... + stkpnt = 0; + stack = new Float64Array( 64 ); + stack[ 0 ] = 0; + stack[ 1 ] = N - 1; + + while ( stkpnt >= 0 ) { + start = stack[ 2*stkpnt ]; + endd = stack[ ( 2*stkpnt ) + 1 ]; + stkpnt -= 1; + + if ( endd - start <= SELECT && endd - start > 0 ) { + // Insertion sort for small partitions + if ( dir === 0 ) { + // Sort in decreasing order + idx1 = offsetD + ( ( start + 1 )*strideD ); + for ( i = start + 1; i <= endd; i++ ) { + idx2 = idx1; // index of i + for ( j = i; j >= start + 1; j-- ) { + idx3 = idx2 - strideD; // index of j - 1 + if ( D[ idx2 ] > D[ idx3 ] ) { + dmnmx = D[ idx2 ]; + D[ idx2 ] = D[ idx3 ]; + D[ idx3 ] = dmnmx; + } else { + break; + } + idx2 -= strideD; + } + idx1 += strideD; + } + } else { + // Sort in increasing order + idx1 = offsetD + ( ( start + 1 )*strideD ); + for ( i = start + 1; i <= endd; i++ ) { + idx2 = idx1; // index of i + for ( j = i; j >= start + 1; j-- ) { + idx3 = idx2 - strideD; // index of j - 1 + if ( D[ idx2 ] < D[ idx3 ] ) { + dmnmx = D[ idx2 ]; + D[ idx2 ] = D[ idx3 ]; + D[ idx3 ] = dmnmx; + } else { + break; + } + idx2 -= strideD; + } + idx1 += strideD; + } + } + } else if ( endd - start > SELECT ) { + // Quicksort partition using median-of-three pivot + d1 = D[ offsetD + ( start*strideD ) ]; + d2 = D[ offsetD + ( endd*strideD ) ]; + i = floor( ( start + endd ) / 2 ); + d3 = D[ offsetD + ( i*strideD ) ]; + + // Find median of d1, d2, d3 + if ( d1 < d2 ) { + if ( d3 < d1 ) { + dmnmx = d1; + } else if ( d3 < d2 ) { + dmnmx = d3; + } else { + dmnmx = d2; + } + } else if ( d3 < d2 ) { + dmnmx = d2; + } else if ( d3 < d1 ) { + dmnmx = d3; + } else { + dmnmx = d1; + } + + if ( dir === 0 ) { + // Partition for decreasing order + i = start; + j = endd; + while ( true ) { + idx1 = offsetD + ( j*strideD ); + while ( D[ idx1 ] < dmnmx ) { + idx1 -= strideD; + j -= 1; + } + + idx2 = offsetD + ( i*strideD ); + while ( D[ idx2 ] > dmnmx ) { + idx2 += strideD; + i += 1; + } + if ( i < j ) { + tmp = D[ idx2 ]; + D[ idx2 ] = D[ idx1 ]; + D[ idx1 ] = tmp; + i += 1; + j -= 1; + } else { + break; + } + } + } else { + // Partition for increasing order + i = start; + j = endd; + while ( true ) { + idx1 = offsetD + ( j*strideD ); + while ( D[ idx1 ] > dmnmx ) { + idx1 -= strideD; + j -= 1; + } + + idx2 = offsetD + ( i*strideD ); + while ( D[ idx2 ] < dmnmx ) { + idx2 += strideD; + i += 1; + } + if ( i < j ) { + tmp = D[ idx2 ]; + D[ idx2 ] = D[ idx1 ]; + D[ idx1 ] = tmp; + i += 1; + j -= 1; + } else { + break; + } + } + } + + // Push sub-partitions onto stack (larger first for bounded stack depth) + idx1 = 2*stkpnt; + if ( j - start > endd - j - 1 ) { + stkpnt += 1; + idx1 += 2; + stack[ idx1 ] = start; + stack[ idx1 + 1 ] = j; + + stkpnt += 1; + idx1 += 2; + stack[ idx1 ] = j + 1; + stack[ idx1 + 1 ] = endd; + } else { + stkpnt += 1; + idx1 += 2; + stack[ idx1 ] = j + 1; + stack[ idx1 + 1 ] = endd; + + stkpnt += 1; + idx1 += 2; + stack[ idx1 ] = start; + stack[ idx1 + 1 ] = j; + } + } + } + return 0; +} + + +// EXPORTS // + +module.exports = dlasrt; diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/index.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/index.js new file mode 100644 index 000000000000..ca906e2dce89 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/index.js @@ -0,0 +1,72 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 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'; + +/** +* Compute all the eigenvalues of the symmetric positive definite tri-diagonal matrix associated with the QD Array `Z` to high relative accuracy. +* +* @module @stdlib/lapack/base/dlasq2 +* +* @example +* var Float64Array = require( '@stdlib/array/float64' ); +* var dlasq2 = require( '@stdlib/lapack/base/dlasq2' ); +* +* var Z = new Float64Array( [ 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); +* +* var info = dlasq2( 4, Z ); +* // info => 0 +* // Z => [ ~115.713, ~83.153, ~62.17, ~42.464, ~83.153, ~20.99, ~0.0, ~0.0, 303.5, ~303.5, 14, 3.8125, 0.0, ~0.0, ~0.172, ~0.096 ] +* +* @example +* var Float64Array = require( '@stdlib/array/float64' ); +* var dlasq2 = require( '@stdlib/lapack/base/dlasq2' ); +* +* var Z = new Float64Array( [ 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); +* +* var info = dlasq2.ndarray( 4, Z, 1, 0 ); +* // info => 0 +* // Z => [ ~115.713, ~83.153, ~62.17, ~42.464, ~83.153, ~20.99, ~0.0, ~0.0, 303.5, ~303.5, 14, 3.8125, 0.0, ~0.0, ~0.172, ~0.096 ] +* +*/ + + +// MODULES // + +var join = require( 'path' ).join; +var tryRequire = require( '@stdlib/utils/try-require' ); +var isError = require( '@stdlib/assert/is-error' ); +var main = require( './main.js' ); + + +// MAIN // + +var dlasq2; +var tmp = tryRequire( join( __dirname, './native.js' ) ); +if ( isError( tmp ) ) { + dlasq2 = main; +} else { + dlasq2 = tmp; +} + + +// EXPORTS // + +module.exports = dlasq2; + +// exports: { "ndarray": "dlasq2.ndarray" } diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/main.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/main.js new file mode 100644 index 000000000000..17b79e773a48 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/main.js @@ -0,0 +1,35 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 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 setReadOnly = require( '@stdlib/utils/define-nonenumerable-read-only-property' ); +var dlasq2 = require( './dlasq2.js' ); +var ndarray = require( './ndarray.js' ); + + +// MAIN // + +setReadOnly( dlasq2, 'ndarray', ndarray ); + + +// EXPORTS // + +module.exports = dlasq2; diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/ndarray.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/ndarray.js new file mode 100644 index 000000000000..1012c8e3d0bc --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/lib/ndarray.js @@ -0,0 +1,54 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 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 base = require( './base.js' ); + + +// MAIN // + +/** +* Computes all the eigenvalues of the symmetric positive definite tri-diagonal matrix associated with the QD Array `Z` to high relative accuracy using alternative indexing semantics. +* +* @private +* @param {integer} N - number of rows/columns in `Z` +* @param {Float64Array} Z - qd array +* @param {integer} strideZ - stride length for `Z` +* @param {NonNegativeInteger} offsetZ - starting index of `Z` +* @returns {integer} status code +* +* @example +* var Float64Array = require( '@stdlib/array/float64' ); +* +* var Z = new Float64Array( [ 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ] ); +* +* var info = dlasq2( 4, Z, 1, 0 ); +* // Z => [ ~115.713, ~83.153, ~62.17, ~42.464, ~83.153, ~20.987, ~0.0, ~0.0, 303.5, ~303.5, 15, 4.0625, 0.0, ~0.0, ~0.172, ~0.096 ] +* // info => 0 +*/ +function dlasq2( N, Z, strideZ, offsetZ ) { // eslint-disable-line stdlib/jsdoc-doctest-decimal-point + return base( N, Z, strideZ, offsetZ ); +} + + +// EXPORTS // + +module.exports = dlasq2; diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/package.json b/lib/node_modules/@stdlib/lapack/base/dlasq2/package.json new file mode 100644 index 000000000000..4ba71eeb0825 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/package.json @@ -0,0 +1,72 @@ +{ + "name": "@stdlib/lapack/base/dlasq2", + "version": "0.0.0", + "description": "LAPACK routine to compute all the eigenvalues of the symmetric positive definite tri-diagonal matrix associated with the QD Array `Z` to high relative accuracy.", + "license": "Apache-2.0", + "author": { + "name": "The Stdlib Authors", + "url": "https://github.com/stdlib-js/stdlib/graphs/contributors" + }, + "contributors": [ + { + "name": "The Stdlib Authors", + "url": "https://github.com/stdlib-js/stdlib/graphs/contributors" + } + ], + "main": "./lib", + "directories": { + "benchmark": "./benchmark", + "doc": "./docs", + "example": "./examples", + "lib": "./lib", + "test": "./test" + }, + "types": "./docs/types", + "scripts": {}, + "homepage": "https://github.com/stdlib-js/stdlib", + "repository": { + "type": "git", + "url": "git://github.com/stdlib-js/stdlib.git" + }, + "bugs": { + "url": "https://github.com/stdlib-js/stdlib/issues" + }, + "dependencies": {}, + "devDependencies": {}, + "engines": { + "node": ">=0.10.0", + "npm": ">2.7.0" + }, + "os": [ + "aix", + "darwin", + "freebsd", + "linux", + "macos", + "openbsd", + "sunos", + "win32", + "windows" + ], + "keywords": [ + "stdlib", + "stdmath", + "mathematics", + "math", + "lapack", + "svd", + "decomposition", + "dlasq2", + "exchange", + "permute", + "permutedims", + "linear", + "algebra", + "subroutines", + "array", + "ndarray", + "float64", + "double", + "float64array" + ] +} diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_2_s_gt_t.json b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_2_s_gt_t.json new file mode 100644 index 000000000000..3b15eee74f82 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_2_s_gt_t.json @@ -0,0 +1,28 @@ +{ + "N": 2, + + "Z": [ 100.0, 50.0, 90.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ], + "strideZ": 1, + "offsetZ": 0, + + "info": 0, + + "Z_out": [ + 1.93484692283495349E+002, + 4.65153077165046582E+001, + 4.65153077165046582E+001, + 0.00000000000000000E+000, + 2.40000000000000000E+002, + 2.40000000000000000E+002, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000 + ] +} diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_2_swap.json b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_2_swap.json new file mode 100644 index 000000000000..473ba984155b --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_2_swap.json @@ -0,0 +1,28 @@ +{ + "N": 2, + + "Z": [ 50.0, 1.0, 100.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ], + "strideZ": 1, + "offsetZ": 0, + + "info": 0, + + "Z_out": [ + 1.01962237244798459E+002, + 4.90377627552015198E+001, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 4.90377627552015198E+001, + 9.80392156862745026E+001, + 0.00000000000000000E+000, + -0.00000000000000000E+000, + 1.51000000000000000E+002, + 1.50999999999999972E+002, + 2.00000000000000000E+000, + 3.75000000000000000E-001, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000 + ] +} diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_3.json b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_3.json new file mode 100644 index 000000000000..d3feca019af5 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_3.json @@ -0,0 +1,32 @@ +{ + "N": 3, + + "Z": [ 5.0, 1.0, 0.0, 1.0, 3.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ], + "strideZ": 1, + "offsetZ": 0, + + "info": 0, + + "Z_out": [ + 6.00000000000000000E+000, + 4.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 4.00000000000000000E+000, + 1.00000000000000000E+000, + 0.00000000000000000E+000, + 3.00000000000000000E+000, + 1.00000000000000000E+001, + 1.00000000000000000E+001, + 2.00000000000000000E+000, + 3.75000000000000000E-001, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000 + ] +} diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_4.json b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_4.json new file mode 100644 index 000000000000..12666e1e052c --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_4.json @@ -0,0 +1,28 @@ +{ + "N": 4, + + "Z": [ 100, 4, 81, 3, 64, 2.5, 49, 2, 36, 1.5, 25, 1, 16, 0.5, 9, 0 ], + "strideZ": 1, + "offsetZ": 0, + + "info": 0, + + "Z_out": [ + 1.15713345308519422E+002, + 8.31528656136522670E+001, + 6.21697557821649909E+001, + 4.24640332956633912E+001, + 8.31528656136522670E+001, + 2.09860054467078996E+001, + 4.26761123604735024E-050, + 4.82003247816166265E-023, + 3.03500000000000000E+002, + 3.03500000000000057E+002, + 1.40000000000000000E+001, + 3.81250000000000000E+000, + 0.00000000000000000E+000, + 6.01321014246954851E-020, + 1.71992034899211055E-001, + 9.56139715683959612E-002 + ] +} diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_5.json b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_5.json new file mode 100644 index 000000000000..649b9a70df27 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_5.json @@ -0,0 +1,40 @@ +{ + "N": 5, + + "Z": [ 25.0, 4.0, 16.0, 3.0, 9.0, 2.0, 4.0, 1.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ], + "strideZ": 1, + "offsetZ": 0, + + "info": 0, + + "Z_out": [ + 3.37308331471525236E+001, + 1.77381457249531636E+001, + 9.04235592736787375E+000, + 3.78556961657009072E+000, + 7.03095583956349190E-001, + 8.69578979758529158E+000, + 5.65616647756791637E-020, + 1.21793792362522094E-044, + 9.04235592736787375E+000, + 1.87246548818156069E-024, + 6.50000000000000000E+001, + 6.50000000000000000E+001, + 1.50000000000000000E+001, + 2.95999999999999996E+000, + 0.00000000000000000E+000, + 2.71122498593797051E-003, + 7.03095583956349190E-001, + 1.02042861520919738E-019, + 6.39182685795127159E-002, + 3.14824862268408834E-002, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000 + ] +} diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_6.json b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_6.json new file mode 100644 index 000000000000..7b85664cd666 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/fixtures/n_eq_6.json @@ -0,0 +1,52 @@ +{ + "N": 6, + + "Z": [ 1.0, 1.033074, 1.1, 0.962918, 1.2, 0.365277, 1.3, 0.392244, 1.4, 1.670586, 101.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ], + "strideZ": 1, + "offsetZ": 0, + + "info": 0, + + "Z_out": [ + 3.52326429814598185E+000, + 2.50257208096376882E+000, + 1.69523907948892338E+000, + 8.41422012060362934E-001, + 1.91015529340962331E-001, + 8.07333001474845546E-001, + 4.79264697420401755E-023, + 9.02086614587196387E-049, + 1.69523907948892338E+000, + 1.51966164781775714E-026, + 8.75351299999999988E+000, + 8.75351299999999988E+000, + 1.90000000000000000E+001, + 3.79999999999999982E+000, + 0.00000000000000000E+000, + 2.31392126923531767E-003, + 1.91015529340962331E-001, + 5.28783478380648471E-019, + 2.11738417771642137E-002, + 1.34821119757860090E-002, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000, + 0.00000000000000000E+000 + ] +} diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/test/test.dlasq2.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/test.dlasq2.js new file mode 100644 index 000000000000..d7796012dbf4 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/test.dlasq2.js @@ -0,0 +1,205 @@ +/* +* @license Apache-2.0 +* +* Copyright (c) 2026 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'; + +var tape = require( 'tape' ); +var Float64Array = require( '@stdlib/array/float64' ); +var isAlmostEqual = require( '@stdlib/assert/is-almost-equal-float64array' ); +var dlasq2 = require( './../lib/dlasq2.js' ); + + +// FIXTURES // + +// All tests have been verified with the LAPACK fortran code. +var N_EQ_4 = require( './fixtures/n_eq_4.json' ); +var N_EQ_5 = require( './fixtures/n_eq_5.json' ); +var N_EQ_6 = require( './fixtures/n_eq_6.json' ); +var N_EQ_3 = require( './fixtures/n_eq_3.json' ); +var N_EQ_2_S_GT_T = require( './fixtures/n_eq_2_s_gt_t.json' ); +var N_EQ_2_SWAP = require( './fixtures/n_eq_2_swap.json' ); + + +// TESTS // + +tape( 'main export is a function', function test( t ) { + t.ok( true, __filename ); + t.strictEqual( typeof dlasq2, 'function', 'main export is a function' ); + t.end(); +}); + +tape( 'the function has an arity of 2', function test( t ) { + t.strictEqual( dlasq2.length, 2, 'returns expected value' ); + t.end(); +}); + +tape( 'the function quick returns (N=0)', function test( t ) { + var Z = new Float64Array( 10 ); + t.strictEqual( dlasq2( 0, Z ), 0, 'returns expected value' ); + t.end(); +}); + +tape( 'the function quick returns (N=1)', function test( t ) { + var Z = new Float64Array( 10 ); + t.strictEqual( dlasq2( 1, Z ), 0, 'returns expected value' ); + + Z[ 0 ] = -1; + t.strictEqual( dlasq2( 1, Z ), -201, 'returns expected value' ); + t.end(); +}); + +tape( 'the function quick returns (N=2)', function test( t ) { + var Z = new Float64Array( [ -1.0, 1.0, 3.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ] ); + t.strictEqual( dlasq2( 2, Z ), -201, 'returns expected value' ); + + Z = new Float64Array( [ 1.0, -1.0, 3.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ] ); + t.strictEqual( dlasq2( 2, Z ), -202, 'returns expected value' ); + + Z = new Float64Array( [ 1.0, 1.0, -3.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ] ); + t.strictEqual( dlasq2( 2, Z ), -203, 'returns expected value' ); + t.end(); +}); + +tape( 'the function quick returns (N=3)', function test( t ) { + var Z = new Float64Array( [ -1.0, 1.0, 3.0, 1.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ] ); + t.strictEqual( dlasq2( 3, Z ), -200, 'returns expected value' ); + + Z = new Float64Array( [ 1.0, -1.0, 3.0, 1.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ] ); + t.strictEqual( dlasq2( 3, Z ), -201, 'returns expected value' ); + + Z = new Float64Array( [ 1.0, 1.0, -3.0, 1.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ] ); + t.strictEqual( dlasq2( 3, Z ), -202, 'returns expected value' ); + + Z = new Float64Array( [ 1.0, 1.0, 3.0, -1.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ] ); + t.strictEqual( dlasq2( 3, Z ), -203, 'returns expected value' ); + + Z = new Float64Array( [ 1.0, 1.0, 3.0, 1.0, -1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ] ); + t.strictEqual( dlasq2( 3, Z ), -205, 'returns expected value' ); + + Z = new Float64Array( [ 3.0, 0.0, 2.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ] ); + t.strictEqual( dlasq2( 3, Z ), 0, 'returns expected value' ); + t.end(); +}); + +tape( 'the function computes all the eigenvalues (N=4)', function test( t ) { + var expectedZ; + var data; + var info; + var Z; + + data = N_EQ_4; + + Z = new Float64Array( data.Z ); + expectedZ = new Float64Array( data.Z_out ); + + info = dlasq2( data.N, Z ); + + t.strictEqual( info, data.info, 'returns expected value' ); + t.strictEqual( isAlmostEqual( Z, expectedZ, 1 ), true, 'returns expected value' ); + t.end(); +}); + +tape( 'the function computes all the eigenvalues (N=2, s > t)', function test( t ) { + var expectedZ; + var data; + var info; + var Z; + + data = N_EQ_2_S_GT_T; + + Z = new Float64Array( data.Z ); + expectedZ = new Float64Array( data.Z_out ); + + info = dlasq2( data.N, Z ); + + t.strictEqual( info, data.info, 'returns expected value' ); + t.strictEqual( isAlmostEqual( Z, expectedZ, 1 ), true, 'returns expected value' ); + t.end(); +}); + +tape( 'the function computes all the eigenvalues (N=2, swap)', function test( t ) { + var expectedZ; + var data; + var info; + var Z; + + data = N_EQ_2_SWAP; + + Z = new Float64Array( data.Z ); + expectedZ = new Float64Array( data.Z_out ); + + info = dlasq2( data.N, Z ); + + t.strictEqual( info, data.info, 'returns expected value' ); + t.strictEqual( isAlmostEqual( Z, expectedZ, 1 ), true, 'returns expected value' ); + t.end(); +}); + +tape( 'the function computes all the eigenvalues (N=3, q2=0)', function test( t ) { + var expectedZ; + var data; + var info; + var Z; + + data = N_EQ_3; + + Z = new Float64Array( data.Z ); + expectedZ = new Float64Array( data.Z_out ); + + info = dlasq2( data.N, Z ); + + t.strictEqual( info, data.info, 'returns expected value' ); + t.strictEqual( isAlmostEqual( Z, expectedZ, 1 ), true, 'returns expected value' ); + t.end(); +}); + +tape( 'the function computes all the eigenvalues (N=5, decreasing)', function test( t ) { + var expectedZ; + var data; + var info; + var Z; + + data = N_EQ_5; + + Z = new Float64Array( data.Z ); + expectedZ = new Float64Array( data.Z_out ); + + info = dlasq2( data.N, Z ); + + t.strictEqual( info, data.info, 'returns expected value' ); + t.strictEqual( isAlmostEqual( Z, expectedZ, 1 ), true, 'returns expected value' ); + t.end(); +}); + +tape( 'the function computes all the eigenvalues (N=6, clustered)', function test( t ) { + var expectedZ; + var data; + var info; + var Z; + + data = N_EQ_6; + + Z = new Float64Array( data.Z ); + expectedZ = new Float64Array( data.Z_out ); + + info = dlasq2( data.N, Z ); + + t.strictEqual( info, data.info, 'returns expected value' ); + t.strictEqual( isAlmostEqual( Z, expectedZ, 1 ), true, 'returns expected value' ); + t.end(); +}); diff --git a/lib/node_modules/@stdlib/lapack/base/dlasq2/test/test.js b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/test.js new file mode 100644 index 000000000000..7a5dc50eb809 --- /dev/null +++ b/lib/node_modules/@stdlib/lapack/base/dlasq2/test/test.js @@ -0,0 +1,82 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2026 The Stdlib Authors. +* +* Licensed under the Apache License, Version 2.0 (the "License"); +* you may not use this file except in compliance with the License. +* You may obtain a copy of the License at +* +* http://www.apache.org/licenses/LICENSE-2.0 +* +* Unless required by applicable law or agreed to in writing, software +* distributed under the License is distributed on an "AS IS" BASIS, +* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +* See the License for the specific language governing permissions and +* limitations under the License. +*/ + +'use strict'; + +// MODULES // + +var tape = require( 'tape' ); +var proxyquire = require( 'proxyquire' ); +var IS_BROWSER = require( '@stdlib/assert/is-browser' ); +var dlasq2 = require( './../lib' ); + + +// VARIABLES // + +var opts = { + 'skip': IS_BROWSER +}; + + +// TESTS // + +tape( 'main export is a function', function test( t ) { + t.ok( true, __filename ); + t.strictEqual( typeof dlasq2, 'function', 'main export is a function' ); + t.end(); +}); + +tape( 'attached to the main export is a method providing an ndarray interface', function test( t ) { + t.strictEqual( typeof dlasq2.ndarray, 'function', 'method is a function' ); + t.end(); +}); + +tape( 'if a native implementation is available, the main export is the native implementation', opts, function test( t ) { + var dlasq2 = proxyquire( './../lib', { + '@stdlib/utils/try-require': tryRequire + }); + + t.strictEqual( dlasq2, mock, 'returns expected value' ); + t.end(); + + function tryRequire() { + return mock; + } + + function mock() { + // Mock... + } +}); + +tape( 'if a native implementation is not available, the main export is a JavaScript implementation', opts, function test( t ) { + var dlasq2; + var main; + + main = require( './../lib/dlasq2.js' ); + + dlasq2 = proxyquire( './../lib', { + '@stdlib/utils/try-require': tryRequire + }); + + t.strictEqual( dlasq2, main, 'returns expected value' ); + t.end(); + + function tryRequire() { + return new Error( 'Cannot find module' ); + } +});