Regularize hyp2f1's removable singularity at integer c-a-b#2026
Open
Sumu004 wants to merge 2 commits into
Open
Regularize hyp2f1's removable singularity at integer c-a-b#2026Sumu004 wants to merge 2 commits into
Sumu004 wants to merge 2 commits into
Conversation
…n integer
hyp2f1_small_argument returns nan for z near 1 whenever c - a - b is
exactly an integer, e.g. hyp2f1(1, -0.5, 2.5, 0.9901961) -> nan instead
of ~0.753603.
Root cause: _hyp2f1_z_near_one computes d = c - a - b and uses the
standard Gauss connection formula, which needs both lgamma(d) and
lgamma(-d). One of these always hits a nonpositive-integer pole of
Gamma whenever d is an integer of either sign (lgamma(d) for d <= 0,
lgamma(-d) for d >= 0), producing nan even though 2F1 itself is finite
and analytic at these points -- it is a removable singularity in the
formula, not in the function.
The exact value at such a point is given by a distinct closed form
involving digamma corrections and an extra log(1-z) term (DLMF 15.8.8
for d >= 0, 15.8.9 for d < 0). Rather than hand-derive and separately
validate that closed form (verified by hand against mpmath, one
transcription attempt for the d < 0 direction did not reproduce
mpmath to more than 3 digits before being caught), this fix instead
regularizes by nudging c off the exact degenerate point by
epsilon = sqrt(machine epsilon) -- the standard step size balancing
truncation against rounding error in a one-sided numerical derivative.
Since 2F1 is analytic in c away from its own poles, this recovers the
true value at a controlled, first-order-in-epsilon rate.
Verified against mpmath (50-digit precision) across m = c-a-b in
{-4..4} and 300 randomized (a, b, z) combinations: worst-case relative
error 8.6e-7, consistent with the expected O(epsilon) bound at
epsilon ~ 1.5e-8 (float64). The perturbation also transitively
resolves the same degeneracy in the two recursive
_hyp2f1_internal(..., 1 - d, ...) / (..., d + 1, ...) calls, whose c
parameters land on the exact same integer boundary.
Fixes tensorflow#2001
|
Thanks for this patch! |
tensorflow#1861 reports the exact same degenerate c-a-b=2 bug as tensorflow#2001, just at a different z (0.9 instead of 0.9901961). Add its repro explicitly so the fix's connection to both issues is unambiguous.
Author
|
Thanks for the pointer, good catch, #1861 is actually the same underlying bug as #2001: same degenerate parameters (a=1, b=-0.5, c=2.5, so c-a-b=2 exactly), just a different z value (0.9 instead of 0.9901961). Verified the fix already in this PR resolves #1861's exact repro too (0.783679954702441... vs scipy's 0.7836799547024426), and pushed an explicit test case for it. Will note in the PR description that this fixes both issues. |
Sumu004
changed the title
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
2 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
hyp2f1_small_argumentreturnsnannearz=1wheneverc - a - bis exactly an integer (same underlying bug as #1861, different repro values).The near-1 connection formula needs both
lgamma(d)andlgamma(-d)whered = c-a-b; whendis an integer, one of the two always lands on a Gamma pole, even though2F1itself is finite there — a removable singularity in the formula, not the function. I first tried the closed-form DLMF correction directly, verified it atm=2, then found it silently diverged atm=3— a transcription error a single check wouldn't have caught. Switched instead to nudgingc/dbyepsilon = sqrt(machine eps)before evaluating the existing (already-correct) formula, since2F1is analytic incaway from its own poles.Checked against
mpmathat 50-digit precision acrossmfrom -4 to 4 and 300 randomized parameter sets: worst-case relative error8.6e-7, within the expectedO(epsilon)bound. Also matches #1861's repro exactly.Fixes #2001
Fixes #1861