Skip to content

Commit

Permalink
feat: add C implementation for blas/base/sspr
Browse files Browse the repository at this point in the history
PR-URL: #4491
Ref: #2039
Co-authored-by: Athan Reines <[email protected]>
Reviewed-by: Athan Reines <[email protected]> 
Co-authored-by: stdlib-bot <[email protected]>
  • Loading branch information
3 people authored Feb 16, 2025
1 parent 7cfb8f7 commit 8d34c0c
Show file tree
Hide file tree
Showing 26 changed files with 2,824 additions and 14 deletions.
85 changes: 77 additions & 8 deletions lib/node_modules/@stdlib/blas/base/sspr/README.md
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@
@license Apache-2.0
Copyright (c) 2024 The Stdlib Authors.
Copyright (c) 2025 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.
Expand Down Expand Up @@ -183,21 +183,63 @@ console.log( AP );
### Usage

```c
TODO
#include "stdlib/blas/base/sspr.h"
```

#### TODO
#### c_sspr( order, uplo, N, alpha, \*X, strideX, \*AP )

TODO.
Performs the symmetric rank 1 operation `A = α*x*x^T + A` where `α` is a scalar, `x` is an `N` element vector, and `A` is an `N` by `N` symmetric matrix supplied in packed form.

```c
#include "stdlib/blas/base/shared.h"

float AP[] = { 1.0f, 2.0f, 3.0f, 1.0f, 2.0f, 1.0f };
const float x[] = { 1.0f, 2.0f, 3.0f };

c_sspr( CblasColMajor, CblasUpper, 3, 1.0f, x, 1, AP );
```
The function accepts the following arguments:
- **order**: `[in] CBLAS_LAYOUT` storage layout.
- **uplo**: `[in] CBLAS_UPLO` specifies whether the upper or lower triangular part of the symmetric matrix `A` should be referenced.
- **N**: `[in] CBLAS_INT` number of elements along each dimension of `A`.
- **alpha**: `[in] float` scalar.
- **X**: `[in] float*` input vector.
- **strideX**: `[in] CBLAS_INT` stride length for `X`.
- **AP**: `[inout] float*` packed form of a symmetric matrix `A`.
```c
TODO
void c_sspr( const CBLAS_LAYOUT order, const CBLAS_UPLO uplo, const CBLAS_INT N, const float alpha, const float *X, const CBLAS_INT strideX, float *AP )
```

TODO
#### c_sspr_ndarray( order, uplo, N, alpha, \*X, strideX, \*AP, strideAP, offsetAP )

Performs the symmetric rank 1 operation `A = α*x*x^T + A` where `α` is a scalar, `x` is an `N` element vector, and `A` is an `N` by `N` symmetric matrix supplied in packed form using alternative indexing semantics.

```c
#include "stdlib/blas/base/shared.h"

float AP[] = { 1.0f, 2.0f, 3.0f, 1.0f, 2.0f, 1.0f };
const float x[] = { 1.0f, 2.0f, 3.0f };

c_sspr_ndarray( CblasColMajor, CblasUpper, 3, 1.0f, x, 1, AP, 1, 0 );
```
The function accepts the following arguments:
- **order**: `[in] CBLAS_LAYOUT` storage layout.
- **uplo**: `[in] CBLAS_UPLO` specifies whether the upper or lower triangular part of the symmetric matrix `A` should be referenced.
- **N**: `[in] CBLAS_INT` number of elements along each dimension of `A`.
- **alpha**: `[in] float` scalar.
- **X**: `[in] float*` input vector.
- **strideX**: `[in] CBLAS_INT` stride length for `X`.
- **AP**: `[inout] float*` packed form of a symmetric matrix `A`.
- **strideAP**: `[in] CBLAS_INT` stride length for `AP`.
- **offsetAP**: `[in] CBLAS_INT` starting index for `AP`.
```c
TODO
void c_sspr_ndarray( const CBLAS_LAYOUT order, const CBLAS_UPLO uplo, const CBLAS_INT N, const float alpha, const float *X, const CBLAS_INT strideX, const CBLAS_INT offsetX, float *AP, const CBLAS_INT strideAP, const CBLAS_INT offsetAP )
```

</section>
Expand All @@ -219,7 +261,34 @@ TODO
### Examples

```c
TODO
#include "stdlib/blas/base/sspr.h"
#include "stdlib/blas/base/shared.h"
#include <stdio.h>

int main( void ) {
// Create strided arrays:
float AP[] = { 1.0f, 2.0f, 3.0f, 1.0f, 2.0f, 1.0f };
const float x[] = { 1.0f, 2.0f, 3.0f };

// Specify the number of elements along each dimension of `A`:
const int N = 3;

// Perform the symmetric rank 1 operation `A = α*x*x^T + A`:
c_sspr( CblasRowMajor, CblasUpper, N, 1.0f, x, 1, AP );

// Print the result:
for ( int i = 0; i < N*(N+1)/2; i++ ) {
printf( "AP[ %i ] = %f\n", i, AP[ i ] );
}

// Perform the symmetric rank 1 operation `A = α*x*x^T + A` using alternative indexing semantics:
c_sspr_ndarray( CblasRowMajor, CblasUpper, N, 1.0f, x, 1, 0, AP, 1, 0 );

// Print the result:
for ( int i = 0; i < N*(N+1)/2; i++ ) {
printf( "AP[ %i ] = %f\n", i, AP[ i ] );
}
}
```
</section>
Expand Down
109 changes: 109 additions & 0 deletions lib/node_modules/@stdlib/blas/base/sspr/benchmark/benchmark.native.js
Original file line number Diff line number Diff line change
@@ -0,0 +1,109 @@
/**
* @license Apache-2.0
*
* Copyright (c) 2025 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 resolve = require( 'path' ).resolve;
var bench = require( '@stdlib/bench' );
var uniform = require( '@stdlib/random/array/uniform' );
var isnanf = require( '@stdlib/math/base/assert/is-nanf' );
var pow = require( '@stdlib/math/base/special/pow' );
var floor = require( '@stdlib/math/base/special/floor' );
var tryRequire = require( '@stdlib/utils/try-require' );
var pkg = require( './../package.json' ).name;


// VARIABLES //

var sspr = tryRequire( resolve( __dirname, './../lib/sspr.native.js' ) );
var opts = {
'skip': ( sspr instanceof Error )
};
var options = {
'dtype': 'float32'
};


// FUNCTIONS //

/**
* Creates a benchmark function.
*
* @private
* @param {PositiveInteger} N - number of elements along each dimension
* @returns {Function} benchmark function
*/
function createBenchmark( N ) {
var AP = uniform( N * ( N + 1 ) / 2, -10.0, 10.0, options );
var x = uniform( N, -10.0, 10.0, options );
return benchmark;

/**
* Benchmark function.
*
* @private
* @param {Benchmark} b - benchmark instance
*/
function benchmark( b ) {
var z;
var i;

b.tic();
for ( i = 0; i < b.iterations; i++ ) {
z = sspr( 'row-major', 'upper', N, 1.0, x, 1, AP );
if ( isnanf( z[ i%z.length ] ) ) {
b.fail( 'should not return NaN' );
}
}
b.toc();
if ( isnanf( z[ i%z.length ] ) ) {
b.fail( 'should not return NaN' );
}
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 = floor( pow( pow( 10, i ), 1.0/2.0 ) );
f = createBenchmark( len );
bench( pkg+'::native:size='+( len * ( len + 1 ) / 2 ), opts, f );
}
}

main();
Original file line number Diff line number Diff line change
@@ -0,0 +1,109 @@
/**
* @license Apache-2.0
*
* Copyright (c) 2025 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 resolve = require( 'path' ).resolve;
var bench = require( '@stdlib/bench' );
var uniform = require( '@stdlib/random/array/uniform' );
var isnanf = require( '@stdlib/math/base/assert/is-nanf' );
var pow = require( '@stdlib/math/base/special/pow' );
var floor = require( '@stdlib/math/base/special/floor' );
var tryRequire = require( '@stdlib/utils/try-require' );
var pkg = require( './../package.json' ).name;


// VARIABLES //

var sspr = tryRequire( resolve( __dirname, './../lib/ndarray.native.js' ) );
var opts = {
'skip': ( sspr instanceof Error )
};
var options = {
'dtype': 'float32'
};


// FUNCTIONS //

/**
* Creates a benchmark function.
*
* @private
* @param {PositiveInteger} N - number of elements along each dimension
* @returns {Function} benchmark function
*/
function createBenchmark( N ) {
var AP = uniform( N * ( N + 1 ) / 2, -10.0, 10.0, options );
var x = uniform( N, -10.0, 10.0, options );
return benchmark;

/**
* Benchmark function.
*
* @private
* @param {Benchmark} b - benchmark instance
*/
function benchmark( b ) {
var z;
var i;

b.tic();
for ( i = 0; i < b.iterations; i++ ) {
z = sspr( 'row-major', 'upper', N, 1.0, x, 1, 0, AP, 1, 0 );
if ( isnanf( z[ i%z.length ] ) ) {
b.fail( 'should not return NaN' );
}
}
b.toc();
if ( isnanf( z[ i%z.length ] ) ) {
b.fail( 'should not return NaN' );
}
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 = floor( pow( pow( 10, i ), 1.0/2.0 ) );
f = createBenchmark( len );
bench( pkg+'::native:ndarray:size='+( len * ( len + 1 ) / 2 ), opts, f );
}
}

main();
Loading

1 comment on commit 8d34c0c

@stdlib-bot
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Coverage Report

Package Statements Branches Functions Lines
blas/base/sspr $\color{green}533/533$
$\color{green}+100.00\%$
$\color{green}49/49$
$\color{green}+100.00\%$
$\color{green}5/5$
$\color{green}+100.00\%$
$\color{green}533/533$
$\color{green}+100.00\%$

The above coverage report was generated for the changes in this push.

Please sign in to comment.