Skip to content

Commit

Permalink
calculate incomplete Gamma function by an asymptotic expansion
Browse files Browse the repository at this point in the history
  • Loading branch information
@rabauke
rabauke committed Feb 21, 2024
1 parent 4cccf4a commit 1af9373
Show file tree
Hide file tree
Showing 3 changed files with 260 additions and 37 deletions.
68 changes: 48 additions & 20 deletions tests/test_special_functions.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -207,16 +207,30 @@ TEMPLATE_TEST_CASE("special_functions", "", float, double, long double) {

SECTION("GammaP") {
using tuple = arg_res_tuple<T, 2>;
auto x_yref{
GENERATE(tuple{T(2.0l), T(0.0l), T(0.0000000000000000000000000000000000000000e+00l)},
tuple{T(2.0l), T(1.0l), T(2.6424111765711535680895245967707826510838e-01l)},
tuple{T(2.0l), T(2.0l), T(5.9399415029016192431800151508254678977711e-01l)},
tuple{T(2.0l), T(3.0l), T(8.0085172652854422808263033739975289347320e-01l)},
tuple{T(2.0l), T(4.0l), T(9.0842180555632909853140989363379378894044e-01l)},
tuple{T(2.0l), T(5.0l), T(9.5957231800548719742018370946110945450690e-01l)},
tuple{T(2.0l), T(6.0l), T(9.8264873476333549103868382798428332475945e-01l)},
tuple{T(2.0l), T(7.0l), T(9.9270494427556387033597491132472573898821e-01l)},
tuple{T(2.0l), T(8.0l), T(9.9698083634887739345060749786797225082620e-01l)})};
auto x_yref{GENERATE(
tuple{T(2.0l), T(0.0l), T(0.0000000000000000000000000000000000000000e+00l)},
tuple{T(2.0l), T(1.0l), T(2.6424111765711535680895245967707826510838e-01l)},
tuple{T(2.0l), T(2.0l), T(5.9399415029016192431800151508254678977711e-01l)},
tuple{T(2.0l), T(3.0l), T(8.0085172652854422808263033739975289347320e-01l)},
tuple{T(2.0l), T(4.0l), T(9.0842180555632909853140989363379378894044e-01l)},
tuple{T(2.0l), T(5.0l), T(9.5957231800548719742018370946110945450690e-01l)},
tuple{T(2.0l), T(6.0l), T(9.8264873476333549103868382798428332475945e-01l)},
tuple{T(2.0l), T(7.0l), T(9.9270494427556387033597491132472573898821e-01l)},
tuple{T(2.0l), T(8.0l), T(9.9698083634887739345060749786797225082620e-01l)},
tuple{T(12.5l), T(9.0l), T(1.5760928441953977820468234649142238289159e-01l)},
tuple{T(12.5l), T(11.0l), T(3.6425597236137699758564203792074202059559e-01l)},
tuple{T(12.5l), T(13.0l), T(5.9240130858255329901102834319115488696740e-01l)},
tuple{T(12.5l), T(15.0l), T(7.7571099516559608601853775825027997448055e-01l)},
tuple{T(12.5l), T(17.0l), T(8.9209214616519247831579080497239130365603e-01l)},
tuple{T(40.0l), T(32.0l), T(9.5602816630230858703771575061447476503886e-02l)},
tuple{T(40.0l), T(34.0l), T(1.7173502161590030279066281224095312864764e-01l)},
tuple{T(40.0l), T(36.0l), T(2.7369639320460330667999631204057987978994e-01l)},
tuple{T(40.0l), T(38.0l), T(3.9408063211889345109193252781723965433945e-01l)},
tuple{T(40.0l), T(40.0l), T(5.2102886106105516334095689757755319315557e-01l)},
tuple{T(40.0l), T(42.0l), T(6.4193141067359395324846555212831384747500e-01l)},
tuple{T(40.0l), T(44.0l), T(7.4692196604585985651799646842674920219553e-01l)},
tuple{T(40.0l), T(46.0l), T(8.3074626260026305803119722265978017795290e-01l)},
tuple{T(40.0l), T(48.0l), T(8.9272353393334340495239865479281817396499e-01l)})};
const T s{x_yref.x[0]};
const T x{x_yref.x[1]};
const T y{trng::math::GammaP(s, x)};
Expand All @@ -225,16 +239,30 @@ TEMPLATE_TEST_CASE("special_functions", "", float, double, long double) {

SECTION("GammaQ") {
using tuple = arg_res_tuple<T, 2>;
auto x_yref{
GENERATE(tuple{T(2.0l), T(0.0l), T(1.0000000000000000000000000000000000000000e+00l)},
tuple{T(2.0l), T(1.0l), T(7.3575888234288464319104754032292173489162e-01l)},
tuple{T(2.0l), T(2.0l), T(4.0600584970983807568199848491745321022289e-01l)},
tuple{T(2.0l), T(3.0l), T(1.9914827347145577191736966260024710652680e-01l)},
tuple{T(2.0l), T(4.0l), T(9.1578194443670901468590106366206211059560e-02l)},
tuple{T(2.0l), T(5.0l), T(4.0427681994512802579816290538890545493098e-02l)},
tuple{T(2.0l), T(6.0l), T(1.7351265236664508961316172015716675240545e-02l)},
tuple{T(2.0l), T(7.0l), T(7.2950557244361296640250886752742610117898e-03l)},
tuple{T(2.0l), T(8.0l), T(3.0191636511226065493925021320277491737981e-03l)})};
auto x_yref{GENERATE(
tuple{T(2.0l), T(0.0l), T(1.0000000000000000000000000000000000000000e+00l)},
tuple{T(2.0l), T(1.0l), T(7.3575888234288464319104754032292173489162e-01l)},
tuple{T(2.0l), T(2.0l), T(4.0600584970983807568199848491745321022289e-01l)},
tuple{T(2.0l), T(3.0l), T(1.9914827347145577191736966260024710652680e-01l)},
tuple{T(2.0l), T(4.0l), T(9.1578194443670901468590106366206211059560e-02l)},
tuple{T(2.0l), T(5.0l), T(4.0427681994512802579816290538890545493098e-02l)},
tuple{T(2.0l), T(6.0l), T(1.7351265236664508961316172015716675240545e-02l)},
tuple{T(2.0l), T(7.0l), T(7.2950557244361296640250886752742610117898e-03l)},
tuple{T(2.0l), T(8.0l), T(3.0191636511226065493925021320277491737981e-03l)},
tuple{T(12.5l), T(9.0l), T(8.4239071558046022179531765350857761710841e-01l)},
tuple{T(12.5l), T(11.0l), T(6.3574402763862300241435796207925797940441e-01l)},
tuple{T(12.5l), T(13.0l), T(4.0759869141744670098897165680884511303260e-01l)},
tuple{T(12.5l), T(15.0l), T(2.2428900483440391398146224174972002551945e-01l)},
tuple{T(12.5l), T(17.0l), T(1.0790785383480752168420919502760869634397e-01l)},
tuple{T(40.0l), T(32.0l), T(9.0439718336976914129622842493855252349611e-01l)},
tuple{T(40.0l), T(34.0l), T(8.2826497838409969720933718775904687135236e-01l)},
tuple{T(40.0l), T(36.0l), T(7.2630360679539669332000368795942012021006e-01l)},
tuple{T(40.0l), T(38.0l), T(6.0591936788110654890806747218276034566055e-01l)},
tuple{T(40.0l), T(40.0l), T(4.7897113893894483665904310242244680684443e-01l)},
tuple{T(40.0l), T(42.0l), T(3.5806858932640604675153444787168615252500e-01l)},
tuple{T(40.0l), T(44.0l), T(2.5307803395414014348200353157325079780447e-01l)},
tuple{T(40.0l), T(46.0l), T(1.6925373739973694196880277734021982204710e-01l)},
tuple{T(40.0l), T(48.0l), T(1.0727646606665659504760134520718182603501e-01l)})};
const T s{x_yref.x[0]};
const T x{x_yref.x[1]};
const T y{trng::math::GammaQ(s, x)};
Expand Down
3 changes: 3 additions & 0 deletions trng/math.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -193,6 +193,9 @@ namespace trng {
return ::std::modf(x, &dummy);
}

using ::std::copysign;
using ::std::signbit;

using ::std::isfinite;
using ::std::isnan;
using ::std::isnormal;
Expand Down
Loading

0 comments on commit 1af9373

Please sign in to comment.