diff --git a/CHANGELOG.md b/CHANGELOG.md index 1193ba921..770419064 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -2143,6 +2143,8 @@ A total of 35 people contributed to this release. Thank you to the following con
+- [`c3895df`](https://github.com/stdlib-js/stdlib/commit/c3895df672126473f5803e93b529bcdd0775c75a) - **refactor:** use utility to resolve an index offset _(by Athan Reines)_ +- [`1654659`](https://github.com/stdlib-js/stdlib/commit/1654659445a6dee281706379770c9cb0498c36c7) - **refactor:** update implementation to reduce code duplication [(#2480)](https://github.com/stdlib-js/stdlib/pull/2480) _(by Aman Bhansali, Athan Reines)_ - [`4d08374`](https://github.com/stdlib-js/stdlib/commit/4d0837401b68ebd4e5b8c38e0214158dbe410a07) - **refactor:** reduce code duplication [(#2479)](https://github.com/stdlib-js/stdlib/pull/2479) _(by Aman Bhansali, Athan Reines)_ - [`a591e05`](https://github.com/stdlib-js/stdlib/commit/a591e052cf1b1515c267781b914c6a482e150425) - **test:** fix test configuration _(by Athan Reines)_ - [`ca56638`](https://github.com/stdlib-js/stdlib/commit/ca566387ddc147c4f15fd012a09bd55713307394) - **feat:** add `blas/base/dspmv` [(#2456)](https://github.com/stdlib-js/stdlib/pull/2456) _(by Aman Bhansali, Athan Reines)_ diff --git a/base/dcopy/lib/dcopy.js b/base/dcopy/lib/dcopy.js index 0e3c85c2d..fad96d6d2 100644 --- a/base/dcopy/lib/dcopy.js +++ b/base/dcopy/lib/dcopy.js @@ -20,6 +20,7 @@ // MODULES // +var stride2offset = require( '@stdlib/strided/base/stride2offset' ); var ndarray = require( './ndarray.js' ); @@ -50,16 +51,8 @@ function dcopy( N, x, strideX, y, strideY ) { if ( N <= 0 ) { return y; } - if ( strideX < 0 ) { - ix = ( 1 - N ) * strideX; - } else { - ix = 0; - } - if ( strideY < 0 ) { - iy = ( 1 - N ) * strideY; - } else { - iy = 0; - } + ix = stride2offset( N, strideX ); + iy = stride2offset( N, strideY ); return ndarray( N, x, strideX, ix, y, strideY, iy ); } diff --git a/base/dspmv/lib/base.js b/base/dspmv/lib/base.js new file mode 100644 index 000000000..613fcd0be --- /dev/null +++ b/base/dspmv/lib/base.js @@ -0,0 +1,138 @@ +/** +* @license Apache-2.0 +* +* Copyright (c) 2024 The Stdlib Authors. +* +* Licensed under the Apache License, Version 2.0 (the "License"); +* you may not use this file except in compliance with the License. +* You may obtain a copy of the License at +* +* http://www.apache.org/licenses/LICENSE-2.0 +* +* Unless required by applicable law or agreed to in writing, software +* distributed under the License is distributed on an "AS IS" BASIS, +* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +* See the License for the specific language governing permissions and +* limitations under the License. +*/ + +'use strict'; + +// MODULES // + +var dfill = require( './../../../ext/base/dfill' ).ndarray; +var dscal = require( './../../../base/dscal' ).ndarray; + + +// MAIN // + +/** +* Performs the matrix-vector operation `y = alpha*A*x + beta*y` where `alpha` and `beta` are scalars, `x` and `y` are `N` element vectors, and `A` is an `N` by `N` symmetric matrix supplied in packed form. +* +* @private +* @param {string} order - storage layout +* @param {string} uplo - specifies whether the upper or lower triangular part of the symmetric matrix `A` is supplied +* @param {NonNegativeInteger} N - number of elements along each dimension of `A` +* @param {number} alpha - scalar constant +* @param {Float64Array} AP - packed form of a symmetric matrix `A` +* @param {Float64Array} x - first input array +* @param {integer} strideX - `x` stride length +* @param {NonNegativeInteger} offsetX - starting `x` index +* @param {number} beta - scalar constant +* @param {Float64Array} y - second input array +* @param {integer} strideY - `y` stride length +* @param {NonNegativeInteger} offsetY - starting `y` index +* @returns {Float64Array} `y` +* +* @example +* var Float64Array = require( '@stdlib/array/float64' ); +* +* var A = new Float64Array( [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0 ] ); +* var x = new Float64Array( [ 1.0, 1.0, 1.0 ] ); +* var y = new Float64Array( [ 1.0, 1.0, 1.0 ] ); +* +* dspmv( 'column-major', 'lower', 3, 1.0, A, x, 1, 0, 1.0, y, 1, 0 ); +* // y => [ ~7.0, ~12.0, ~15.0 ] +*/ +function dspmv( order, uplo, N, alpha, AP, x, strideX, offsetX, beta, y, strideY, offsetY ) { // eslint-disable-line max-params, max-len + var temp1; + var temp2; + var ix; + var iy; + var jx; + var jy; + var kk; + var kx; + var ky; + var j; + var k; + + if ( N === 0 || ( alpha === 0.0 && beta === 1.0 ) ) { + return y; + } + // Form: y = beta*y + if ( beta !== 1.0 ) { + if ( beta === 0.0 ) { + dfill( N, 0.0, y, strideY, offsetY ); + } else { + dscal( N, beta, y, strideY, offsetY ); + } + } + if ( alpha === 0.0 ) { + return y; + } + // Form: y = alpha*A*x + y + kx = offsetX; + ky = offsetY; + kk = 0; + if ( + ( order === 'row-major' && uplo === 'upper' ) || + ( order === 'column-major' && uplo === 'lower' ) + ) { + jx = kx; + jy = ky; + for ( j = 0; j < N; j++ ) { + temp1 = alpha * x[ jx ]; + temp2 = 0.0; + y[ jy ] += temp1 * AP[ kk ]; + ix = jx; + iy = jy; + for ( k = kk + 1; k < kk + N - j; k++ ) { + ix += strideX; + iy += strideY; + y[ iy ] += temp1 * AP[ k ]; + temp2 += AP[ k ] * x[ ix ]; + } + y[ jy ] += alpha * temp2; + jx += strideX; + jy += strideY; + kk += N - j; + } + return y; + } + // ( order === 'row-major' && uplo === 'lower') || ( order === 'column-major' && uplo === 'upper' ) + jx = kx; + jy = ky; + for ( j = 0; j < N; j++ ) { + temp1 = alpha * x[ jx ]; + temp2 = 0.0; + ix = kx; + iy = ky; + for ( k = kk; k < kk + j; k++ ) { + y[ iy ] += temp1 * AP[ k ]; + temp2 += AP[ k ] * x[ ix ]; + ix += strideX; + iy += strideY; + } + y[ jy ] += ( temp1 * AP[ kk + j ] ) + ( alpha * temp2 ); + jx += strideX; + jy += strideY; + kk += j + 1; + } + return y; +} + + +// EXPORTS // + +module.exports = dspmv; diff --git a/base/dspmv/lib/dspmv.js b/base/dspmv/lib/dspmv.js index 3052a4d36..4455c6b88 100644 --- a/base/dspmv/lib/dspmv.js +++ b/base/dspmv/lib/dspmv.js @@ -20,10 +20,10 @@ // MODULES // -var dfill = require( './../../../ext/base/dfill' ); -var dscal = require( './../../../base/dscal' ); var isLayout = require( './../../../base/assert/is-layout' ); var isMatrixTriangle = require( './../../../base/assert/is-matrix-triangle' ); +var format = require( '@stdlib/string/format' ); +var base = require( './base.js' ); // MAIN // @@ -59,109 +59,35 @@ var isMatrixTriangle = require( './../../../base/assert/is-matrix-triangle' ); * // y => [ ~7.0, ~12.0, ~15.0 ] */ function dspmv( order, uplo, N, alpha, AP, x, strideX, beta, y, strideY ) { - var temp1; - var temp2; - var ix; - var iy; - var jx; - var jy; - var kk; - var kx; - var ky; - var sy; - var j; - var k; + var offsetX; + var offsetY; if ( !isLayout( order ) ) { - throw new TypeError( 'invalid argument. First argument must be a valid order. Value: `%s`.', order ); + throw new TypeError( format( 'invalid argument. First argument must be a valid order. Value: `%s`.', order ) ); } if ( !isMatrixTriangle( uplo ) ) { - throw new TypeError( 'invalid argument. Second argument must specify whether the lower or upper triangular matrix is supplied. Value: `%s`.', uplo ); + throw new TypeError( format( 'invalid argument. Second argument must specify whether the lower or upper triangular matrix is supplied. Value: `%s`.', uplo ) ); } if ( N < 0 ) { - throw new RangeError( 'invalid argument. Third argument must be a nonnegative integer. Value: `%d`.', N ); + throw new RangeError( format( 'invalid argument. Third argument must be a nonnegative integer. Value: `%d`.', N ) ); } if ( strideX === 0 ) { - throw new RangeError( 'invalid argument. Seventh argument must be non-zero. Value: `%d`.', strideX ); + throw new RangeError( format( 'invalid argument. Seventh argument must be non-zero. Value: `%d`.', strideX ) ); } if ( strideY === 0 ) { - throw new RangeError( 'invalid argument. Tenth argument must be non-zero. Value: `%d`.', strideY ); - } - if ( N === 0 || ( alpha === 0.0 && beta === 1.0 ) ) { - return y; - } - // Form: y = beta*y - sy = strideY; - if ( beta !== 1.0 ) { - if ( beta === 0.0 ) { - dfill( N, 0.0, y, strideY ); - } else { - if ( strideY < 0 ) { - sy = -sy; - } - dscal( N, beta, y, sy ); - } - } - if ( alpha === 0.0 ) { - return y; + throw new RangeError( format( 'invalid argument. Tenth argument must be non-zero. Value: `%d`.', strideY ) ); } if ( strideX > 0 ) { - kx = 0; + offsetX = 0; } else { - kx = ( 1 - N ) * strideX; + offsetX = ( 1 - N ) * strideX; } if ( strideY > 0 ) { - ky = 0; + offsetY = 0; } else { - ky = ( 1 - N ) * strideY; - } - // Form: y = alpha*A*x + y - kk = 0; - if ( - ( order === 'row-major' && uplo === 'upper' ) || - ( order === 'column-major' && uplo === 'lower' ) - ) { - jx = kx; - jy = ky; - for ( j = 0; j < N; j++ ) { - temp1 = alpha * x[ jx ]; - temp2 = 0.0; - y[ jy ] += temp1 * AP[ kk ]; - ix = jx; - iy = jy; - for ( k = kk + 1; k < kk + N - j; k++ ) { - ix += strideX; - iy += strideY; - y[ iy ] += temp1 * AP[ k ]; - temp2 += AP[ k ] * x[ ix ]; - } - y[ jy ] += alpha * temp2; - jx += strideX; - jy += strideY; - kk += N - j; - } - return y; - } - // ( order === 'row-major' && uplo === 'lower') || ( order === 'column-major' && uplo === 'upper' ) - jx = kx; - jy = ky; - for ( j = 0; j < N; j++ ) { - temp1 = alpha * x[ jx ]; - temp2 = 0.0; - ix = kx; - iy = ky; - for ( k = kk; k < kk + j; k++ ) { - y[ iy ] += temp1 * AP[ k ]; - temp2 += AP[ k ] * x[ ix ]; - ix += strideX; - iy += strideY; - } - y[ jy ] += ( temp1 * AP[ kk + j ] ) + ( alpha * temp2 ); - jx += strideX; - jy += strideY; - kk += j + 1; + offsetY = ( 1 - N ) * strideY; } - return y; + return base( order, uplo, N, alpha, AP, x, strideX, offsetX, beta, y, strideY, offsetY ); // eslint-disable-line max-len } diff --git a/base/dspmv/lib/ndarray.js b/base/dspmv/lib/ndarray.js index 5f3e5a607..44b35d3c7 100644 --- a/base/dspmv/lib/ndarray.js +++ b/base/dspmv/lib/ndarray.js @@ -20,10 +20,10 @@ // MODULES // -var dfill = require( './../../../ext/base/dfill' ).ndarray; -var dscal = require( './../../../base/dscal' ).ndarray; var isLayout = require( './../../../base/assert/is-layout' ); var isMatrixTriangle = require( './../../../base/assert/is-matrix-triangle' ); +var format = require( '@stdlib/string/format' ); +var base = require( './base.js' ); // MAIN // @@ -61,96 +61,22 @@ var isMatrixTriangle = require( './../../../base/assert/is-matrix-triangle' ); * // y => [ ~7.0, ~12.0, ~15.0 ] */ function dspmv( order, uplo, N, alpha, AP, x, strideX, offsetX, beta, y, strideY, offsetY ) { // eslint-disable-line max-params, max-len - var temp1; - var temp2; - var ix; - var iy; - var jx; - var jy; - var kk; - var kx; - var ky; - var j; - var k; - if ( !isLayout( order ) ) { - throw new TypeError( 'invalid argument. First argument must be a valid order. Value: `%s`.', order ); + throw new TypeError( format( 'invalid argument. First argument must be a valid order. Value: `%s`.', order ) ); } if ( !isMatrixTriangle( uplo ) ) { - throw new TypeError( 'invalid argument. Second argument must specify whether the lower or upper triangular matrix is supplied. Value: `%s`.', uplo ); + throw new TypeError( format( 'invalid argument. Second argument must specify whether the lower or upper triangular matrix is supplied. Value: `%s`.', uplo ) ); } if ( N < 0 ) { - throw new RangeError( 'invalid argument. Third argument must be a nonnegative integer. Value: `%d`.', N ); + throw new RangeError( format( 'invalid argument. Third argument must be a nonnegative integer. Value: `%d`.', N ) ); } if ( strideX === 0 ) { - throw new RangeError( 'invalid argument. Seventh argument must be non-zero. Value: `%d`.', strideX ); + throw new RangeError( format( 'invalid argument. Seventh argument must be non-zero. Value: `%d`.', strideX ) ); } if ( strideY === 0 ) { - throw new RangeError( 'invalid argument. Eleventh argument must be non-zero. Value: `%d`.', strideY ); - } - if ( N === 0 || ( alpha === 0.0 && beta === 1.0 ) ) { - return y; - } - // Form: y = beta*y - if ( beta !== 1.0 ) { - if ( beta === 0.0 ) { - dfill( N, 0.0, y, strideY, offsetY ); - } else { - dscal( N, beta, y, strideY, offsetY ); - } - } - if ( alpha === 0.0 ) { - return y; - } - // Form: y = alpha*A*x + y - kx = offsetX; - ky = offsetY; - kk = 0; - if ( - ( order === 'row-major' && uplo === 'upper' ) || - ( order === 'column-major' && uplo === 'lower' ) - ) { - jx = kx; - jy = ky; - for ( j = 0; j < N; j++ ) { - temp1 = alpha * x[ jx ]; - temp2 = 0.0; - y[ jy ] += temp1 * AP[ kk ]; - ix = jx; - iy = jy; - for ( k = kk + 1; k < kk + N - j; k++ ) { - ix += strideX; - iy += strideY; - y[ iy ] += temp1 * AP[ k ]; - temp2 += AP[ k ] * x[ ix ]; - } - y[ jy ] += alpha * temp2; - jx += strideX; - jy += strideY; - kk += N - j; - } - return y; - } - // ( order === 'row-major' && uplo === 'lower') || ( order === 'column-major' && uplo === 'upper' ) - jx = kx; - jy = ky; - for ( j = 0; j < N; j++ ) { - temp1 = alpha * x[ jx ]; - temp2 = 0.0; - ix = kx; - iy = ky; - for ( k = kk; k < kk + j; k++ ) { - y[ iy ] += temp1 * AP[ k ]; - temp2 += AP[ k ] * x[ ix ]; - ix += strideX; - iy += strideY; - } - y[ jy ] += ( temp1 * AP[ kk + j ] ) + ( alpha * temp2 ); - jx += strideX; - jy += strideY; - kk += j + 1; + throw new RangeError( format( 'invalid argument. Eleventh argument must be non-zero. Value: `%d`.', strideY ) ); } - return y; + return base( order, uplo, N, alpha, AP, x, strideX, offsetX, beta, y, strideY, offsetY ); // eslint-disable-line max-len }