Skip to content

Commit caff697

Browse files
aayush0325aayush0325
committed
feat: lapack/base/dorm2r
--- type: pre_commit_static_analysis_report description: Results of running static analysis checks when committing changes. report: - task: lint_filenames status: passed - task: lint_editorconfig status: passed - task: lint_markdown status: na - task: lint_package_json status: na - task: lint_repl_help status: na - task: lint_javascript_src status: passed - task: lint_javascript_cli status: na - task: lint_javascript_examples status: na - task: lint_javascript_tests status: na - task: lint_javascript_benchmarks status: na - task: lint_python status: na - task: lint_r status: na - task: lint_c_src status: na - task: lint_c_examples status: na - task: lint_c_benchmarks status: na - task: lint_c_tests_fixtures status: na - task: lint_shell status: na - task: lint_typescript_declarations status: na - task: lint_typescript_tests status: na - task: lint_license_headers status: passed ---
1 parent f79bf70 commit caff697

File tree

6 files changed

+585
-0
lines changed

6 files changed

+585
-0
lines changed

‎lib/node_modules/@stdlib/lapack/base/dorm2r/lib/base.js

Lines changed: 143 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,143 @@
1+
/**
2+
* @license Apache-2.0
3+
*
4+
* Copyright (c) 2025 The Stdlib Authors.
5+
*
6+
* Licensed under the Apache License, Version 2.0 (the "License");
7+
* you may not use this file except in compliance with the License.
8+
* You may obtain a copy of the License at
9+
*
10+
* http://www.apache.org/licenses/LICENSE-2.0
11+
*
12+
* Unless required by applicable law or agreed to in writing, software
13+
* distributed under the License is distributed on an "AS IS" BASIS,
14+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
15+
* See the License for the specific language governing permissions and
16+
* limitations under the License.
17+
*/
18+
19+
'use strict';
20+
21+
// MODULES //
22+
23+
var dlarf1f = require( './dlarf1f.js' );
24+
25+
26+
// MAIN //
27+
28+
/**
29+
* Multiplies a matrix `C` by an orthogonal matrix `Q` from a QR factorization.
30+
*
31+
* ## Notes
32+
*
33+
* `dorm2r` overwrites the general real M by N matrix C with
34+
*
35+
* - `Q * C` if side = 'left' and trans = 'no-transpose'
36+
* - `Q^T * C` if side = 'left' and trans = 'transpose'
37+
* - `C * Q` if side = 'right' and trans = 'no-transpose'
38+
* - `C * Q^T` if side = 'right' and trans = 'transpose'
39+
*
40+
* where Q is a real orthogonal matrix defined as the product of K elementary reflectors `Q = H(1) H(2) ... H(K)` as returned by `dgeqrf`. Q is of order M if side = 'left' and of order N if side = 'right'.
41+
*
42+
* @private
43+
* @param {string} side - specifies the side of multiplication with `C`
44+
* @param {string} trans - specifies the operation to be performed
45+
* @param {NonNegativeInteger} M - number of rows in matrix `C`
46+
* @param {NonNegativeInteger} N - number of columns in matrix `C`
47+
* @param {NonNegativeInteger} K - number of elementary reflectors whose product defines the matrix `Q`
48+
* @param {Float64Array} A - input matrix containing the elementary reflectors
49+
* @param {integer} strideA1 - stride length for the first dimension of `A`
50+
* @param {integer} strideA2 - stride length for the second dimension of `A`
51+
* @param {NonNegativeInteger} offsetA - starting index for `A`
52+
* @param {Float64Array} tau - array containing the scalar factors of the elementary reflectors
53+
* @param {integer} strideTau - stride length for `tau`
54+
* @param {NonNegativeInteger} offsetTau - starting index for `tau`
55+
* @param {Float64Array} C - input/output matrix to be multiplied
56+
* @param {integer} strideC1 - stride length for the first dimension of `C`
57+
* @param {integer} strideC2 - stride length for the second dimension of `C`
58+
* @param {NonNegativeInteger} offsetC - starting index for `C`
59+
* @param {Float64Array} work - workspace array for intermediate calculations
60+
* @param {integer} strideWork - stride length for `work`
61+
* @param {NonNegativeInteger} offsetWork - starting index for `work`
62+
* @returns {Float64Array} the modified matrix `C`
63+
*
64+
* @example
65+
* var Float64Array = require( '@stdlib/array/float64' );
66+
*
67+
* var A = new Float64Array( [ 1.0, 0.0, 0.0, 2.0, 4.0, 0.0, 3.0, 5.0, 6.0 ] );
68+
* var tau = new Float64Array( [ 7.0, 8.0, 9.0 ] );
69+
* var C = new Float64Array( [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0 ] );
70+
* var work = new Float64Array( 3 );
71+
*
72+
* dorm2r( 'left', 'no-transpose', 3, 3, 3, A, 3, 1, 0, tau, 1, 0, C, 3, 1, 0, work, 1, 0 );
73+
* // C => <Float64Array>[ -261638.0, -298618.0, -335598.0, -521066.0, -594715.0, -668364.0, -773933.0, -883324.0, -992715.0 ]
74+
*/
75+
function dorm2r( side, trans, M, N, K, A, strideA1, strideA2, offsetA, tau, strideTau, offsetTau, C, strideC1, strideC2, offsetC, work, strideWork, offsetWork ) { // eslint-disable-line max-params, max-len
76+
var aOffset;
77+
var cOffset;
78+
var tauVal;
79+
var notran;
80+
var left;
81+
var i1;
82+
var i2;
83+
var i3;
84+
var mi;
85+
var ni;
86+
var ic;
87+
var jc;
88+
var i;
89+
90+
// Quick return if possible
91+
if ( M === 0 || N === 0 || K === 0 ) {
92+
return C;
93+
}
94+
95+
left = ( side === 'left' );
96+
notran = ( trans === 'no-transpose' );
97+
98+
// Determine loop indices based on side and trans
99+
if ( ( left && !notran ) || ( !left && notran ) ) {
100+
i1 = 1;
101+
i2 = K;
102+
i3 = 1;
103+
} else {
104+
i1 = K;
105+
i2 = 1;
106+
i3 = -1;
107+
}
108+
109+
// Initialize dimensions
110+
if ( left ) {
111+
ni = N;
112+
jc = 0; // 0-based indexing
113+
} else {
114+
mi = M;
115+
ic = 0; // 0-based indexing
116+
}
117+
118+
// Main loop
119+
for ( i = i1; ( i3 > 0 ) ? ( i <= i2 ) : ( i >= i2 ); i += i3 ) {
120+
if ( left ) {
121+
// H(i) is applied to C(i:m,1:n)
122+
mi = M - i + 1;
123+
ic = i - 1; // Convert to 0-based indexing
124+
} else {
125+
// H(i) is applied to C(1:m,i:n)
126+
ni = N - i + 1;
127+
jc = i - 1; // Convert to 0-based indexing
128+
}
129+
130+
// Apply H(i) using `dlarf1f`
131+
tauVal = tau[ offsetTau + ( ( i - 1 ) * strideTau ) ];
132+
aOffset = offsetA + ( ( i - 1 ) * strideA1 ) + ( ( i - 1 ) * strideA2 );
133+
cOffset = offsetC + ( ic * strideC1 ) + ( jc * strideC2 );
134+
dlarf1f( side, mi, ni, A, strideA1, aOffset, tauVal, C, strideC1, strideC2, cOffset, work, strideWork, offsetWork ); // eslint-disable-line max-len
135+
}
136+
137+
return C;
138+
}
139+
140+
141+
// EXPORTS //
142+
143+
module.exports = dorm2r;

‎lib/node_modules/@stdlib/lapack/base/dorm2r/lib/dlarf1f.js

Lines changed: 128 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,128 @@
1+
/**
2+
* @license Apache-2.0
3+
*
4+
* Copyright (c) 2025 The Stdlib Authors.
5+
*
6+
* Licensed under the Apache License, Version 2.0 (the "License");
7+
* you may not use this file except in compliance with the License.
8+
* You may obtain a copy of the License at
9+
*
10+
* http://www.apache.org/licenses/LICENSE-2.0
11+
*
12+
* Unless required by applicable law or agreed to in writing, software
13+
* distributed under the License is distributed on an "AS IS" BASIS,
14+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
15+
* See the License for the specific language governing permissions and
16+
* limitations under the License.
17+
*/
18+
19+
'use strict';
20+
21+
/* eslint-disable max-len */
22+
23+
// MODULES //
24+
25+
var iladlc = require( '@stdlib/lapack/base/iladlc' ).ndarray;
26+
var iladlr = require( '@stdlib/lapack/base/iladlr' ).ndarray;
27+
var dgemv = require( '@stdlib/blas/base/dgemv' ).ndarray;
28+
var dger = require( '@stdlib/blas/base/dger' ).ndarray;
29+
var daxpy = require( '@stdlib/blas/base/daxpy' ).ndarray;
30+
var dscal = require( '@stdlib/blas/base/dscal' ).ndarray;
31+
32+
33+
// MAIN //
34+
35+
/**
36+
* Applies a real elementary reflector `H = I - tau * v * v ^ T` to a real M by N matrix `C`.
37+
*
38+
* ## Notes
39+
*
40+
* - `work` should have `N` indexed elements if side = `left` and `M` indexed elements if side = `right`.
41+
* - `V` should have `1 + (M-1) * abs(strideV)` indexed elements if side = `left` and `1 + (N-1) * abs(strideV)` indexed elements if side = `right`.
42+
* - `C` is overwritten by `H * C` if side = `left` and `C * H` if side = `right`.
43+
*
44+
* @private
45+
* @param {string} side - specifies the side of multiplication with `C`. Use `left` to form `H * C` and `right` to from `C * H`.
46+
* @param {NonNegativeInteger} M - number of rows in `C`
47+
* @param {NonNegativeInteger} N - number of columns in `C`
48+
* @param {Float64Array} V - the vector `v` in the representation of `H`
49+
* @param {integer} strideV - stride length for `V`
50+
* @param {NonNegativeInteger} offsetV - starting index for `V`
51+
* @param {number} tau - the value of `tau` in representation of `H`
52+
* @param {Float64Array} C - input matrix
53+
* @param {integer} strideC1 - stride of the first dimension of `C`
54+
* @param {integer} strideC2 - stride of the second dimension of `C`
55+
* @param {NonNegativeInteger} offsetC - starting index for `C`
56+
* @param {Float64Array} work - workspace array
57+
* @param {integer} strideWork - stride length for `work`
58+
* @param {NonNegativeInteger} offsetWork - starting index for `work`
59+
* @returns {Float64Array} `C * H` or `H * C`
60+
*
61+
* @example
62+
* var Float64Array = require( '@stdlib/array/float64' );
63+
*
64+
* var C = new Float64Array( [ 1.0, 5.0, 9.0, 2.0, 6.0, 10.0, 3.0, 7.0, 11.0, 4.0, 8.0, 12.0 ] );
65+
* var V = new Float64Array( [ 0.5, 0.5, 0.5, 0.5 ] );
66+
* var work = new Float64Array( 3 );
67+
*
68+
* var out = dlarf1f( 'left', 4, 3, V, 1, 0, 1.0, C, 3, 1, 0, work, 1, 0 );
69+
* // returns <Float64Array>[ -4.5, -10.5, -16.5, -0.75, -1.75, -2.75, 0.25, -0.75, -1.75, 1.25, 0.25, -0.75 ]
70+
*/
71+
function dlarf1f( side, M, N, V, strideV, offsetV, tau, C, strideC1, strideC2, offsetC, work, strideWork, offsetWork ) { // eslint-disable-line max-params
72+
var lastv;
73+
var lastc;
74+
var i;
75+
76+
lastv = 1;
77+
lastc = 0;
78+
if ( tau !== 0.0 ) {
79+
if ( side === 'left' ) {
80+
lastv = M;
81+
} else {
82+
lastv = N;
83+
}
84+
85+
// i points to the last element in V
86+
i = offsetV + ( ( lastv - 1 ) * strideV );
87+
88+
// Move i to the last non-zero element in V
89+
while ( lastv > 0 && V[ i ] === 0.0 ) {
90+
lastv -= 1;
91+
i -= strideV;
92+
}
93+
if ( side === 'left' ) {
94+
lastc = iladlc( lastv + 1, N, C, strideC1, strideC2, offsetC ) + 1; // to account for the difference between zero-based and one-based indexing
95+
} else {
96+
lastc = iladlr( M, lastv + 1, C, strideC1, strideC2, offsetC ) + 1; // to account for the difference between zero-based and one-based indexing
97+
}
98+
// lastc is zero if a matrix is filled with zeros
99+
}
100+
if ( lastc === 0 ) {
101+
// Returns C unchanged if tau is zero or all elements in C are zero
102+
return C;
103+
}
104+
if ( side === 'left' ) {
105+
if ( lastv === 0 ) {
106+
dscal( lastc, 1.0 - tau, C, strideC2, offsetC ); // scale the first row
107+
} else {
108+
dgemv( 'transpose', lastv-1, lastc, 1.0, C, strideC1, strideC2, offsetC + strideC1, V, strideV, offsetV + strideV, 0.0, work, strideWork, offsetWork ); // C( 1, 0 ) is accessed here
109+
daxpy( lastc, 1.0, C, strideC2, offsetC, work, strideWork, offsetWork ); // operates on the first row of C
110+
daxpy( lastc, -tau, work, strideWork, offsetWork, C, strideC2, offsetC ); // operates on the first row of C
111+
dger( lastv-1, lastc, -tau, V, strideV, offsetV + strideV, work, strideWork, offsetWork, C, strideC1, strideC2, offsetC + strideC1 ); // C( 1, 0 ) is accessed here
112+
}
113+
} else if ( lastv === 0 ) {
114+
dscal( lastc, 1.0 - tau, C, strideC1, offsetC ); // scale the first column
115+
} else {
116+
dgemv( 'no-transpose', lastc, lastv-1, 1.0, C, strideC1, strideC2, offsetC + strideC2, V, strideV, offsetV + strideV, 0.0, work, strideWork, offsetWork ); // C( 0, 1 ) is accessed here
117+
daxpy( lastc, 1.0, C, strideC1, offsetC, work, strideWork, offsetWork ); // operates on the first column of C
118+
daxpy( lastc, -tau, work, strideWork, offsetWork, C, strideC1, offsetC ); // operates on the first column of C
119+
dger( lastc, lastv-1, -tau, work, strideWork, offsetWork, V, strideV, offsetV + strideV, C, strideC1, strideC2, offsetC + strideC2 ); // C( 0, 1 ) is accessed here
120+
}
121+
122+
return C;
123+
}
124+
125+
126+
// EXPORTS //
127+
128+
module.exports = dlarf1f;

‎lib/node_modules/@stdlib/lapack/base/dorm2r/lib/dorm2r.js

Lines changed: 116 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,116 @@
1+
/**
2+
* @license Apache-2.0
3+
*
4+
* Copyright (c) 2025 The Stdlib Authors.
5+
*
6+
* Licensed under the Apache License, Version 2.0 (the "License");
7+
* you may not use this file except in compliance with the License.
8+
* You may obtain a copy of the License at
9+
*
10+
* http://www.apache.org/licenses/LICENSE-2.0
11+
*
12+
* Unless required by applicable law or agreed to in writing, software
13+
* distributed under the License is distributed on an "AS IS" BASIS,
14+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
15+
* See the License for the specific language governing permissions and
16+
* limitations under the License.
17+
*/
18+
19+
'use strict';
20+
21+
// MODULES //
22+
23+
var isLayout = require( '@stdlib/blas/base/assert/is-layout' );
24+
var isOperationSide = require( '@stdlib/blas/base/assert/is-operation-side' );
25+
var isTransposeOperation = require( '@stdlib/blas/base/assert/is-transpose-operation' );
26+
var isRowMajor = require( '@stdlib/ndarray/base/assert/is-row-major-string' );
27+
var isColumnMajor = require( '@stdlib/ndarray/base/assert/is-column-major-string' );
28+
var max = require( '@stdlib/math/base/special/max' );
29+
var format = require( '@stdlib/string/format' );
30+
var base = require( './base.js' );
31+
32+
33+
// MAIN //
34+
35+
/**
36+
* Multiplies a matrix `C` by an orthogonal matrix `Q` from a QR factorization.
37+
*
38+
* ## Notes
39+
*
40+
* `dorm2r` overwrites the general real M by N matrix C with
41+
*
42+
* - `Q * C` if side = 'left' and trans = 'no-transpose'
43+
* - `Q^T * C` if side = 'left' and trans = 'transpose'
44+
* - `C * Q` if side = 'right' and trans = 'no-transpose'
45+
* - `C * Q^T` if side = 'right' and trans = 'transpose'
46+
*
47+
* where Q is a real orthogonal matrix defined as the product of K elementary reflectors `Q = H(1) H(2) ... H(K)` as returned by `dgeqrf`. Q is of order M if side = 'left' and of order N if side = 'right'.
48+
*
49+
* @param {string} order - storage layout
50+
* @param {string} side - specifies the side of multiplication with `C`
51+
* @param {string} trans - specifies the operation to be performed
52+
* @param {NonNegativeInteger} M - number of rows in matrix `C`
53+
* @param {NonNegativeInteger} N - number of columns in matrix `C`
54+
* @param {NonNegativeInteger} K - number of elementary reflectors whose product defines the matrix `Q`
55+
* @param {Float64Array} A - input matrix containing the elementary reflectors
56+
* @param {PositiveInteger} LDA - stride of the first dimension of `A` (a.k.a., leading dimension of the matrix `A`)
57+
* @param {Float64Array} tau - array containing the scalar factors of the elementary reflectors
58+
* @param {Float64Array} C - input/output matrix to be multiplied
59+
* @param {PositiveInteger} LDC - stride of the first dimension of `C` (a.k.a., leading dimension of the matrix `C`)
60+
* @param {Float64Array} work - workspace array for intermediate calculations
61+
* @throws {TypeError} first argument must be a valid order
62+
* @throws {TypeError} second argument must be a valid operation side
63+
* @throws {TypeError} third argument must be a valid transpose operation
64+
* @throws {RangeError} eighth argument must be greater than or equal to max(1,K) when side is 'left' or max(1,M) when side is 'right'
65+
* @throws {RangeError} twelfth argument must be greater than or equal to max(1,M)
66+
* @returns {Float64Array} the modified matrix `C`
67+
*
68+
* @example
69+
* var Float64Array = require( '@stdlib/array/float64' );
70+
*
71+
* var A = new Float64Array( [ 1.0, 0.0, 0.0, 2.0, 4.0, 0.0, 3.0, 5.0, 6.0 ] );
72+
* var tau = new Float64Array( [ 7.0, 8.0, 9.0 ] );
73+
* var C = new Float64Array( [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0 ] );
74+
* var work = new Float64Array( 3 );
75+
*
76+
* dorm2r( 'row-major', 'left', 'no-transpose', 3, 3, 3, A, 3, tau, C, 3, work );
77+
* // C => <Float64Array>[ -261638.0, -298618.0, -335598.0, -521066.0, -594715.0, -668364.0, -773933.0, -883324.0, -992715.0 ]
78+
*/
79+
function dorm2r( order, side, trans, M, N, K, A, LDA, tau, C, LDC, work ) { // eslint-disable-line max-params
80+
var sa1;
81+
var sa2;
82+
var sc1;
83+
var sc2;
84+
if ( !isLayout( order ) ) {
85+
throw new TypeError( format( 'invalid argument. First argument must be a valid order. Value: `%s`.', order ) );
86+
}
87+
if ( !isOperationSide( side ) ) {
88+
throw new TypeError( format( 'invalid argument. Second argument must be a valid operation side. Value: `%s`.', side ) );
89+
}
90+
if ( !isTransposeOperation( trans ) ) {
91+
throw new TypeError( format( 'invalid argument. Third argument must be a valid transpose operation. Value: `%s`.', trans ) );
92+
}
93+
if ( isRowMajor( order ) && LDC < max( 1, N ) ) {
94+
throw new RangeError( format( 'invalid argument. Twelfth argument must be greater than or equal to max(1,%d). Value: `%d`.', N, LDC ) );
95+
}
96+
if ( isColumnMajor( order ) && LDC < max( 1, M ) ) {
97+
throw new RangeError( format( 'invalid argument. Twelfth argument must be greater than or equal to max(1,%d). Value: `%d`.', M, LDC ) );
98+
}
99+
if ( isColumnMajor( order ) ) {
100+
sa1 = 1;
101+
sa2 = LDA;
102+
sc1 = 1;
103+
sc2 = LDC;
104+
} else {
105+
sa1 = LDA;
106+
sa2 = 1;
107+
sc1 = LDC;
108+
sc2 = 1;
109+
}
110+
return base( side, trans, M, N, K, A, sa1, sa2, 0, tau, 1, 0, C, sc1, sc2, 0, work, 1, 0 ); // eslint-disable-line max-len
111+
}
112+
113+
114+
// EXPORTS //
115+
116+
module.exports = dorm2r;

0 commit comments

Comments
 (0)