Skip to content

Commit 842b5ef

Browse files
authored
feat: add /blas/ext/base/dnancusumkbn2
Signed-off-by: Shabareesh Shetty <[email protected]>
1 parent 87abb74 commit 842b5ef

File tree

1 file changed

+108
-0
lines changed
  • lib/node_modules/@stdlib/blas/ext/base/dnancusumkbn2/lib

1 file changed

+108
-0
lines changed
Lines changed: 108 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,108 @@
1+
/**
2+
* @license Apache-2.0
3+
*
4+
* Copyright (c) 2020 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 isnan = require( '@stdlib/math/base/assert/is-nan' );
24+
var abs = require( '@stdlib/math/base/special/abs' );
25+
26+
27+
// MAIN //
28+
29+
/**
30+
* Computes the cumulative sum of double-precision floating-point strided array elements, ignoring `NaN` values and using a second-order iterative Kahan–Babuška algorithm.
31+
*
32+
* ## Method
33+
*
34+
* - This implementation uses a second-order iterative Kahan–Babuška algorithm, as described by Klein (2005).
35+
*
36+
* ## References
37+
*
38+
* - Klein, Andreas. 2005. "A Generalized Kahan-Babuška-Summation-Algorithm." _Computing_ 76 (3): 279–93. doi:[10.1007/s00607-005-0139-x](https://doi.org/10.1007/s00607-005-0139-x).
39+
*
40+
* @param {PositiveInteger} N - number of indexed elements
41+
* @param {number} sum - initial sum
42+
* @param {Float64Array} x - input array
43+
* @param {integer} strideX - stride length for `x`
44+
* @param {NonNegativeInteger} offsetX - starting index for `x`
45+
* @param {Float64Array} y - output array
46+
* @param {integer} strideY - stride length for `y`
47+
* @param {NonNegativeInteger} offsetY - starting index for `y`
48+
* @returns {Float64Array} output array
49+
*
50+
* @example
51+
* var Float64Array = require( '@stdlib/array/float64' );
52+
*
53+
* var x = new Float64Array( [ 2.0, 1.0, 2.0, -2.0, -2.0, 2.0, 3.0, NaN ] );
54+
* var y = new Float64Array( x.length );
55+
*
56+
* var v = dnancusumkbn2( 4, 0.0, x, 2, 1, y, 1, 0 );
57+
* // returns <Float64Array>[ 1.0, -1.0, 1.0, 1.0, 0.0, 0.0, 0.0, 0.0 ]
58+
*/
59+
function dcusumkbn2( N, sum, x, strideX, offsetX, y, strideY, offsetY ) {
60+
var ccs;
61+
var ix;
62+
var iy;
63+
var cs;
64+
var cc;
65+
var v;
66+
var t;
67+
var c;
68+
var i;
69+
70+
if ( N <= 0 ) {
71+
return y;
72+
}
73+
ix = offsetX;
74+
iy = offsetY;
75+
76+
ccs = 0.0; // second order correction term for lost low order bits
77+
cs = 0.0; // first order correction term for lost low order bits
78+
for ( i = 0; i < N; i++ ) {
79+
v = x[ ix ];
80+
if ( isnan( v ) === false ) {
81+
t = sum + v;
82+
if ( abs( sum ) >= abs( v ) ) {
83+
c = (sum-t) + v;
84+
} else {
85+
c = (v-t) + sum;
86+
}
87+
sum = t;
88+
t = cs + c;
89+
if ( abs( cs ) >= abs( c ) ) {
90+
cc = (cs-t) + c;
91+
} else {
92+
cc = (c-t) + cs;
93+
}
94+
cs = t;
95+
ccs += cc;
96+
97+
y[ iy ] = sum + cs + ccs;
98+
}
99+
ix += strideX;
100+
iy += strideY;
101+
}
102+
return y;
103+
}
104+
105+
106+
// EXPORTS //
107+
108+
module.exports = dnancusumkbn2;

0 commit comments

Comments
 (0)