[Math] Fill in some gaps in custom derivative support for Clad - #23078
Merged
Conversation
Clad computes Hessians in reverse-over-forward mode, so the forward pass looks up `<name>_pushforward` custom derivatives. For ROOT::Math::inc_gamma, inc_gamma_c, landau_pdf, and landau_cdf no pushforward existed, and Clad silently fell back to numerical differentiation, producing badly wrong second derivatives (e.g. the Poisson NLL Hessian was off by more than 10 %). Implement the pushforwards with derivatives that stay exact under Clad's reverse pass over their bodies: - dP/dx of inc_gamma is the gamma density exp((a-1) ln x - x - lgamma(a)); its hand-written pullback gives the closed-form d2P/dx2 and d2P/dxda (via the digamma function). - dP/da is evaluated exactly by seeding the existing algorithm-differentiated inc_gamma_pullback. Only d2P/da2 has no closed form; it uses a central difference of the exact first derivative (accuracy ~1e-10, far better than differentiating the original function numerically). It is not exercised by RooFit, where the first argument is data. - landau_cdf derivatives are fully closed-form via the pdf (dF/dx = p(v)/xi with v = (x - x0)/xi). - landau_pdf needs p'(v), evaluated exactly by seeding the existing landau_pdf_pullback, and p''(v), again a central difference of the exact p'. Helpers live in clad::custom_derivatives::ROOT::Math so that Clad resolves their custom pullbacks in the same namespace. Full 2x2 Hessians of inc_gamma / inc_gamma_c were validated against finite differences at points covering both Cephes branches (<= 4e-5 relative), and the Landau NLL Hessians against finite differences of the exact Clad gradient (<= 1e-4 relative). 🤖 Done with the help of AI
The central difference for the second derivative of the standardized Landau density in landau_pdf_dv_pullback() straddled the branch boundaries of the piecewise-rational DENLAN approximation. The exact first derivative has small jumps at those seams, and dividing a jump by the ~1e-5 step size ruined the difference quotient: up to ~50 % error right at v = 1 (~3 % at v = -1, ~0.3 % at v = 5), far outside the validated 1e-4 accuracy. Shift the stencil sideways when it would straddle a seam so both points stay on the branch that contains v; the error at v = 1 drops to ~0.3 %, which is the intrinsic inter-branch consistency of DENLAN itself. Away from the seams the stencil is unchanged. Also document two approximations that the code presented as exact: - landau_cdf_pushforward() takes dF/dx = landau_pdf, which is exact only for the mathematical Landau distribution; ROOT's cdf (DISLAN) and density (DENLAN) are independent rational approximations consistent to ~1e-7, so the forward derivative differs at that level from the exact derivative of the implemented cdf and from landau_cdf_pullback(), which differentiates the DISLAN algorithm. - In inc_gamma_da_pullback(), for a <= h the lower stencil point leaves the domain where inc_gamma_da() silently returns zero, making d2P/da2 unreliable there (not exercised by RooFit, where a is data). 🤖 Done with the help of AI
The pullbacks for landau_pdf and landau_cdf were pasted machine-generated reverse-mode code: eight branches of ~30 lines each whose only real output is the derivative of the standardized function with respect to v = (x - x0) / xi, buried in adjoint bookkeeping that could not be compared against the primals by eye. Replace them with scalar functions landau_pdf_dv() and landau_cdf_dv() that differentiate each branch of the CERNLIB DENLAN/DISLAN rational approximations directly, built on small horner()/horner_deriv()/ rational_deriv() helpers. The branch structure and coefficient tables now mirror the primals in PdfFuncMathCore.cxx and ProbFuncMathCore.cxx line by line, and the pullbacks collapse to the three-line chain rule in x, x0 and xi. The existing landau_pdf_dv() used by the Hessian helpers becomes this direct implementation instead of seeding the pullback and discarding two of its three outputs. Also drop provably dead code from the inc_gamma pullbacks: constants the generated code never uses, and the adjoint machinery for `t` in inc_gamma_c_pullback(), which is identically zero because t only controls the loop exit (its only updates are self-cancelling). Validated against the previous code over ~1000 samples covering all branches, seams and domain guards: the inc_gamma outputs are bitwise identical, the Landau outputs agree to 3e-14 relative (evaluation-order rounding only), the new derivative functions match finite differences of the primals at the FD noise floor, and the Hessian branch-seam behavior is unchanged. 🤖 Done with the help of AI
Test Results 23 files 23 suites 3d 15h 6m 8s ⏱️ Results for commit e59a1e3. |
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
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.
This is improving the code for existing custom pullbacks, and also adds missing pullbacks and pushforwards for the planned RooFit tests that cover Hessians generated with Clad.