Skip to content

Commit

Permalink
Fix inverse_discrete_quantile for large guess
Browse files Browse the repository at this point in the history
If `guess` passed to `inverse_discrete_quantile` cannot be represented
as floating point number, it is possible that `guess + 1` or `guess - 1`
does not change the value at all and we are stuck in infinite loop
inside `round_to_floor` or `round_to_ceil`.  Fix this by
increase/decrease more than 1 in these cases.

Example code to reproduce this:
```c++
boost::math::binomial_distribution<> dist(9079765771874083840, 0.561815);
boost::math::quantile(dist, 0.0365346);
```
  • Loading branch information
Yuhta committed Jul 31, 2023
1 parent 281f491 commit 320db24
Showing 1 changed file with 14 additions and 10 deletions.
24 changes: 14 additions & 10 deletions include/boost/math/distributions/detail/inv_discrete_quantile.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -302,15 +302,17 @@ inline typename Dist::value_type round_to_floor(const Dist& d, typename Dist::va
//
while(result != 0)
{
cc = result - 1;
cc = result;
for (typename Dist::value_type d = 1; cc == result; d += 1)
{
cc = result - d;
}
if(cc < support(d).first)
break;
pp = c ? cdf(complement(d, cc)) : cdf(d, cc);
if(pp == p)
result = cc;
else if(c ? pp > p : pp < p)
if(c ? pp > p : pp < p)
break;
result -= 1;
result = cc;
}

return result;
Expand All @@ -336,15 +338,17 @@ inline typename Dist::value_type round_to_ceil(const Dist& d, typename Dist::val
//
while(true)
{
cc = result + 1;
cc = result;
for (typename Dist::value_type d = 1; cc == result; d += 1)
{
cc = result + d;
}
if(cc > support(d).second)
break;
pp = c ? cdf(complement(d, cc)) : cdf(d, cc);
if(pp == p)
result = cc;
else if(c ? pp < p : pp > p)
if(c ? pp < p : pp > p)
break;
result += 1;
result = cc;
}

return result;
Expand Down

0 comments on commit 320db24

Please sign in to comment.