Skip to content

refactor: improve C implementation of math/base/special/hyp2f1 #11353

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
kgryte merged 2 commits into from
Apr 10, 2026
+93 −2
Merged

refactor: improve C implementation of math/base/special/hyp2f1 #11353

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
430e38a
refactor: improve C implementation of math/base/special/hyp2f1
officiallyanee Apr 10, 2026
c7e1747
style: follow conventions in for loop
officiallyanee Apr 10, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 3 additions & 0 deletions lib/node_modules/@stdlib/math/base/special/hyp2f1/manifest.json
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -47,6 +47,7 @@
"@stdlib/math/base/special/ln",
"@stdlib/math/base/special/abs",
"@stdlib/math/base/special/pow",
"@stdlib/math/base/special/max",
"@stdlib/constants/float64/pinf",
"@stdlib/constants/float64/nan"
]
Expand All @@ -72,6 +73,7 @@
"@stdlib/math/base/special/ln",
"@stdlib/math/base/special/abs",
"@stdlib/math/base/special/pow",
"@stdlib/math/base/special/max",
"@stdlib/constants/float64/pinf",
"@stdlib/constants/float64/nan"
]
Expand All @@ -97,6 +99,7 @@
"@stdlib/math/base/special/ln",
"@stdlib/math/base/special/abs",
"@stdlib/math/base/special/pow",
"@stdlib/math/base/special/max",
"@stdlib/constants/float64/pinf",
"@stdlib/constants/float64/nan"
]
Expand Down
46 changes: 44 additions & 2 deletions lib/node_modules/@stdlib/math/base/special/hyp2f1/src/main.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -39,6 +39,7 @@
#include "stdlib/math/base/special/ln.h"
#include "stdlib/math/base/special/abs.h"
#include "stdlib/math/base/special/pow.h"
#include "stdlib/math/base/special/max.h"
#include "stdlib/math/base/assert/is_nan.h"
#include "stdlib/constants/float64/pinf.h"
#include "stdlib/constants/float64/nan.h"
Expand Down Expand Up @@ -235,6 +236,38 @@ static double hys2f1( double a, double b, const double c, const double x, double
return s;
}

/**
* Evaluates 2F1(a, b; b; x) when `b = c` is a negative integer using AMS55 #15.4.2.
*
* @param a first parameter
* @param b second parameter (equals c, a non-positive integer)
* @param x argument
* @return function value, or NaN if precision is insufficient
*/
static double hyp2f1NegCEqualBC( const double a, const double b, const double x ) {
double collectorMax;
double collector;
double sum;
double k;

collectorMax = 1.0;
collector = 1.0;
sum = 1.0;

if ( stdlib_base_abs( b ) >= 1.0e5 ) {
return STDLIB_CONSTANT_FLOAT64_NAN;
}
for ( k = 1.0; k <= -b; k ++ ) {
collector *= ( ( a + k - 1.0 ) * x ) / k;
collectorMax = stdlib_base_max( stdlib_base_abs( collector ), collectorMax );
sum += collector;
}
if ( 1.0e-16 * ( 1.0 + ( collectorMax / stdlib_base_abs( sum ) ) ) > 1.0e-7 ) {
return STDLIB_CONSTANT_FLOAT64_NAN;
}
return sum;
}

/**
* Applies transformations for `|x|` near unity before performing a power series expansion.
*
Expand Down Expand Up @@ -289,7 +322,7 @@ static double hyt2f1( const double a, const double b, const double c, const doub
d = c - a - b;
id = stdlib_base_round( d );

if ( x > 0.9 && !negIntA && !negIntB ) {
if ( x > 0.85 && !negIntA && !negIntB ) {
if ( isInteger( d ) == false ) {
// Try the power series first:
y = hys2f1( a, b, c, x, &err );
Expand Down Expand Up @@ -386,7 +419,11 @@ static double hyt2f1( const double a, const double b, const double c, const doub
y = -y;
}
q = stdlib_base_pow( s, id );
y = ( id > 0.0 ) ? y*q : y1*q;
if ( id > 0.0 ) {
y *= q;
} else {
y1 *= q;
}
y += y1;
}
*loss = err;
Expand Down Expand Up @@ -477,6 +514,11 @@ double stdlib_base_hyp2f1( const double a, const double b, const double c, const
if ( ax < 1.0 || x == -1.0 ) {
if ( b == c ) {
// 2F1(a,b;b;x) = (1-x)**(-a):
if ( negIntB ) {
// For negative integer b=c use the finite polynomial (AMS55 #15.4.2):
y = hyp2f1NegCEqualBC( a, b, x );
return y;
}
y = stdlib_base_pow( s, -a );
return y;
}
Expand Down
46 changes: 46 additions & 0 deletions lib/node_modules/@stdlib/math/base/special/hyp2f1/test/test.native.js
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -43,6 +43,8 @@ var edgeCases1 = require( './fixtures/python/edge_cases1.json' );
var edgeCases2 = require( './fixtures/python/edge_cases2.json' );
var edgeCases3 = require( './fixtures/python/edge_cases3.json' );
var edgeCases4 = require( './fixtures/python/edge_cases4.json' );
var edgeCases5 = require( './fixtures/python/edge_cases5.json' );
var edgeCases6 = require( './fixtures/python/edge_cases6.json' );
var outliers = require( './fixtures/python/outliers.json' );


Expand Down Expand Up @@ -254,6 +256,50 @@ tape( 'the function correctly evaluates the hypergeometric function', opts, func
t.end();
});

tape( 'the function correctly evaluates the hypergeometric function', function test( t ) {
var expected;
var a;
var b;
var c;
var x;
var v;
var i;

a = edgeCases5.a;
b = edgeCases5.b;
c = edgeCases5.c;
x = edgeCases5.x;
expected = edgeCases5.expected;

for ( i = 0; i < x.length; i++ ) {
v = hyp2f1( a[ i ], b[ i ], c[ i ], x[ i ] );
t.strictEqual( isAlmostSameValue( v, expected[ i ], 12 ), true, 'returns expected value.' );
Copy link
Copy Markdown
Contributor Author

@officiallyanee officiallyanee Apr 10, 2026

Choose a reason for hiding this comment

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

A note on both tol I did the following to replicate compiler optimization:

CFLAGS="-ffp-contract=on" make install-node-addons NODE_ADDONS_PATTERN="math/base/special/hyp2f1"
make test TESTS_FILTER=".*/math/base/special/hyp2f1/.*"

and they were still passing.

}
t.end();
});

tape( 'the function correctly evaluates the hypergeometric function (edge case 6)', function test( t ) {
var expected;
var a;
var b;
var c;
var x;
var v;
var i;

a = edgeCases6.a;
b = edgeCases6.b;
c = edgeCases6.c;
x = edgeCases6.x;
expected = edgeCases6.expected;

for ( i = 0; i < x.length; i++ ) {
v = hyp2f1( a[ i ], b[ i ], c[ i ], x[ i ] );
t.strictEqual( isAlmostSameValue( v, expected[ i ], 1 ), true, 'returns expected value.' );
}
t.end();
});

tape( 'the function correctly evaluates the hypergeometric function', opts, function test( t ) {
var expected;
var a;
Expand Down