|
| 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 dlamch = require( '@stdlib/lapack/base/dlamch' ); |
| 24 | +var sqrt = require( '@stdlib/math/base/special/sqrt' ); |
| 25 | +var abs = require( '@stdlib/math/base/special/fast/abs' ); |
| 26 | +var sign = require( '@stdlib/math/base/special/copysign' ); |
| 27 | +var min = require( '@stdlib/math/base/special/fast/min' ); |
| 28 | +var max = require( '@stdlib/math/base/special/fast/max' ); |
| 29 | + |
| 30 | + |
| 31 | +// VARIABLES // |
| 32 | + |
| 33 | +var safmin = dlamch( 'safe minimum' ); |
| 34 | +var safmax = 1.0 / safmin; |
| 35 | +var rtmin = sqrt( safmin ); |
| 36 | +var rtmax = sqrt( safmax / 2.0 ); |
| 37 | + |
| 38 | + |
| 39 | +// MAIN // |
| 40 | + |
| 41 | +/** |
| 42 | +* Generates a plane rotation of a vector such that the following equation is satisfied. |
| 43 | +* |
| 44 | +* ```tex |
| 45 | +* \(\begin{bmatrix} C & S \\ -S & C \end{bmatrix} \begin{bmatrix} F \\ G \end{bmatrix} = \begin{bmatrix} R \\ 0 \end{bmatrix}\) |
| 46 | +* ``` |
| 47 | +* |
| 48 | +* @private |
| 49 | +* @param {Float64Array} F - single element array representing the first component of the vector to be rotated |
| 50 | +* @param {NonNegativeInteger} offsetF - starting index for `F` |
| 51 | +* @param {Float64Array} G - single element array representing the second component of the vector to be rotated |
| 52 | +* @param {NonNegativeInteger} offsetG - starting index for `G` |
| 53 | +* @param {Float64Array} CS - single element array overwritten by the cosine of the rotation |
| 54 | +* @param {NonNegativeInteger} offsetCS - starting index for `CS` |
| 55 | +* @param {Float64Array} SN - single element array overwritten by the sine of the rotation |
| 56 | +* @param {NonNegativeInteger} offsetSN - starting index for `SN` |
| 57 | +* @param {Float64Array} R - single element array overwritten by the non-zero component of the rotated vector |
| 58 | +* @param {NonNegativeInteger} offsetR - starting index for `R` |
| 59 | +* @returns {void} |
| 60 | +* |
| 61 | +* @example |
| 62 | +* var Float64Array = require( '@stdlib/array/float64' ); |
| 63 | +* |
| 64 | +* var F = new Float64Array( [ 3.0 ] ); |
| 65 | +* var G = new Float64Array( [ 4.0 ] ); |
| 66 | +* var CS = new Float64Array( 1 ); |
| 67 | +* var SN = new Float64Array( 1 ); |
| 68 | +* var R = new Float64Array( 1 ); |
| 69 | +* |
| 70 | +* dlartg( F, 0, G, 0, CS, 0, SN, 0, R, 0 ); |
| 71 | +* // R => <Float64Array>[ 5.0 ] |
| 72 | +* // CS => <Float64Array>[ 0.6 ] |
| 73 | +* // SN => <Float64Array>[ 0.8 ] |
| 74 | +*/ |
| 75 | +function dlartg( F, offsetF, G, offsetG, CS, offsetCS, SN, offsetSN, R, offsetR ) { // eslint-disable-line max-len |
| 76 | + var f1; |
| 77 | + var g1; |
| 78 | + var fs; |
| 79 | + var gs; |
| 80 | + var u; |
| 81 | + var d; |
| 82 | + |
| 83 | + f1 = abs( F[ offsetF ] ); |
| 84 | + g1 = abs( G[ offsetG ] ); |
| 85 | + |
| 86 | + if ( G[ offsetG ] === 0.0 ) { |
| 87 | + CS[ offsetCS ] = 1.0; |
| 88 | + SN[ offsetSN ] = 0.0; |
| 89 | + R[ offsetR ] = F[ offsetF ]; |
| 90 | + } else if ( F[ offsetF ] === 0.0 ) { |
| 91 | + CS[ offsetCS ] = 0.0; |
| 92 | + SN[ offsetSN ] = sign( 1.0, G[ offsetG ] ); |
| 93 | + R[ offsetR ] = g1; |
| 94 | + } else if ( f1 > rtmin && f1 < rtmax && g1 > rtmin && g1 < rtmax ) { |
| 95 | + d = sqrt( (F[ offsetF ]*F[ offsetF ]) + (G[ offsetG ]*G[ offsetG ]) ); |
| 96 | + CS[ offsetCS ] = f1 / d; |
| 97 | + R[ offsetR ] = sign( d, F[ offsetF ] ); |
| 98 | + SN[ offsetSN ] = G[ offsetG ] / R[ offsetR ]; |
| 99 | + } else { |
| 100 | + u = min( safmax, max( safmin, f1, g1 ) ); |
| 101 | + fs = F[ offsetF ] / u; |
| 102 | + gs = G[ offsetG ] / u; |
| 103 | + d = sqrt( (fs*fs) + (gs*gs) ); |
| 104 | + CS[ offsetCS ] = abs( fs ) / d; |
| 105 | + R[ offsetR ] = sign( d, F[ offsetF ] ); |
| 106 | + SN[ offsetSN ] = gs / R[ offsetR ]; |
| 107 | + R[ offsetR ] *= u; |
| 108 | + } |
| 109 | +} |
| 110 | + |
| 111 | + |
| 112 | +// EXPORTS // |
| 113 | + |
| 114 | +module.exports = dlartg; |
0 commit comments