boost/math
2
0
Fork 0
mirror of https://github.com/boostorg/math.git synced 2026-01-19 04:22:09 +00:00
Code Issues Projects Releases Wiki Activity

Fix for scipy issue 17916 (#939)

Browse Source
This commit is contained in:
Matt Borland
2023-02-03 09:10:47 -08:00
committed by GitHub
parent 01630b8b12
commit fb85bdb510
3 changed files with 36 additions and 4 deletions
Download Patch File Download Diff File Expand all files Collapse all files

8
include/boost/math/distributions/non_central_chi_squared.hpp
View File

@@ -13,7 +13,7 @@
#include <boost/math/distributions/fwd.hpp>
#include <boost/math/special_functions/gamma.hpp> // for incomplete gamma. gamma_q
#include <boost/math/special_functions/bessel.hpp> // for cyl_bessel_i
#include <boost/math/special_functions/round.hpp> // for iround
#include <boost/math/special_functions/round.hpp> // for llround
#include <boost/math/distributions/complement.hpp> // complements
#include <boost/math/distributions/chi_squared.hpp> // central distribution
#include <boost/math/distributions/detail/common_error_handling.hpp> // error checks
@@ -71,7 +71,7 @@ namespace boost
// k is chosen as the peek of the Poisson weights, which
// will occur *before* the largest term.
//
int k = iround(lambda, pol);
long long k = llround(lambda, pol);
// Forwards and backwards Poisson weights:
T poisf = boost::math::gamma_p_derivative(static_cast<T>(1 + k), lambda, pol);
T poisb = poisf * k / lambda;
@@ -215,7 +215,7 @@ namespace boost
// function, which ocurrs *after* the largest term in the
// sum.
//
int k = iround(del, pol);
long long k = llround(del, pol);
T a = n / 2 + k;
// Central chi squared term for forward iteration:
T gamkf = boost::math::gamma_p(a, x, pol);
@@ -294,7 +294,7 @@ namespace boost
T n2 = n / 2;
T l2 = lambda / 2;
T sum = 0;
int k = itrunc(l2);
long long k = lltrunc(l2);
T pois = gamma_p_derivative(static_cast<T>(k + 1), l2, pol) * gamma_p_derivative(static_cast<T>(n2 + k), x2);
if(pois == 0)
return 0;

1
test/Jamfile.v2
View File

@@ -896,6 +896,7 @@ test-suite distribution_tests :
[ run scipy_issue_14901.cpp ../../test/build//boost_unit_test_framework ]
[ run scipy_issue_17146.cpp ../../test/build//boost_unit_test_framework ]
[ run scipy_issue_17388.cpp ../../test/build//boost_unit_test_framework ]
[ run scipy_issue_17916.cpp ../../test/build//boost_unit_test_framework ]
;
test-suite mp :

31
test/scipy_issue_17916.cpp Normal file
View File

@@ -0,0 +1,31 @@
// Copyright Matt Borland, 2023
// Use, modification and distribution are subject to the
// Boost Software License, Version 1.0. (See accompanying file
// LICENSE or copy at http://www.boost.org/LICENSE_1_0.txt)
//
// See: https://github.com/scipy/scipy/issues/17916
#include <boost/math/distributions/non_central_chi_squared.hpp>
#include "math_unit_test.hpp"
int main(void)
{
auto dist = boost::math::non_central_chi_squared(2.0, 4820232647677555.0);
double test_pdf;
double test_cdf;
try
{
test_pdf = boost::math::pdf(dist, 2.0);
test_cdf = boost::math::cdf(dist, 2.0);
}
catch (...)
{
return 1;
}
CHECK_ULP_CLOSE(test_pdf, 0.0, 1);
CHECK_ULP_CLOSE(test_cdf, 0.0, 1);
return boost::math::test::report_errors();
}