Repository navigation
Read a factorial where it has no logarithm at all (#754) - #767
Conversation
The two readings of Stirling's expansion already here both need the factorial to be under a logarithm: one supplies it, by turning f^g into e^(g * ln f), and the other rewrites a logarithm that is written down. Neither reaches x! / x^x, which has none, so that shape and its relatives had no limit. The logarithm is supplied instead: for a positive expression lim H is e^(lim ln H), and ln H is where the expansion applies. The guard is the power of the factorial the expression depends on -- exactly the coefficient ln(f!) carries in ln H -- so the dropped 1/(12f) has to vanish against it, as everywhere else. (x!)^x depends on its factorial through the power x, x/(12x) does not vanish, and it is refused. **Not by substituting e^(Stirling(f)) for the factorial**, which is the obvious move and was measured to be far worse. It puts an e to a large exponent into the expression, the machinery evaluates that constant to a hundred-digit decimal, and everything downstream drags it: lim x->+oo (x!/e^x)^(1/x) went from half a second to over a minute, on an expression the library already answered. Going through the logarithm keeps the only e in the final answer. That was found by timing the calculus tests -- 46 s to 2 m 9 s -- and not by reading the code. The power is read off the shape rather than by differentiating and simplifying. The symbolic route needs a full Simplify to cancel the target back out, since InnerSimplified leaves p * (1/x^x) / (p/x^x) standing, and it is asked on every expression holding a factorial. The rewrite that reads a logarithm is broadened from ln(f!) to any logarithm holding a diverging factorial, and this is not tidiness. Without it, ln(x!/x^x) stays one opaque node, the descent assembles -oo/+oo, reads NaN off it, and three limits that had been slow non-answers become wrong answers. Measuring caught that too; the sweep now reports no NaN at all. Measured over 70 generated factorial expressions against master: 24 timeouts and 9 unevaluated became values, and nothing was lost -- 33 answered before, 66 now, with timeouts 25 -> 1. The one that remains is ln((x+1)!)/(x*ln(x)), already recorded when PR #766 landed. Values spot-checked numerically at up to x = 1e6. Calculus tests 46 s -> 55 s for fifteen more of them, casbench 113/117 with 0 wrong and 0 timeout, unit tests 5322 -> 5337 passing with 15 added, F# 130, propcheck 0 failures, rootcheck 596/596, simpsweep 0 disagreements. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
|
Force-pushed: the first version had a performance regression I nearly shipped, and the fix is a better design. Substituting — on an expression the library already answered without any of this. The calculus tests went from 46 s to 2 m 9 s. I found it by timing the suite, not by reading the code, and only after redoing a comparison I had botched (I had The fix is to never put It is better on every axis measured:
The 55 s is for fifteen more tests than master's 46 s — the regression is gone rather than reduced. Two smaller changes came with it. The guard is now read off the shape ( |
e2db5bf to
11aa29f
Compare
Completes the factorial work of #754 / #765. Follows #764 and #766, both merged.
What was wrong
The two readings of Stirling's expansion already in the tree both need the factorial to be under a logarithm: #764 supplies one, by turning
f^gintoe^(g * ln f), and #766 rewrites a logarithm that is written down. Neither reachesx! / x^x, which has none.A third guard, for a third kind of error
The factorial is replaced outright by
e^(Stirling(f)), so what may be dropped is no longer an additive error that vanishes but a relative one that tends to 1 —f!ise^(Stirling(f)) * e^(1/(12f) + ...).A factor tending to 1 leaves a limit alone only where the expression does not amplify it, and what decides that is the elasticity
d(ln H)/d(ln F)— the power of the factorial the expression effectively depends on. A constant is fine.(x!)^xhas elasticityx, so the error would arrive ase^(x/(12x)), which ise^(1/12)and not 1, and it is refused.The three guards are worth keeping apart, and that is why they are three functions:
StirlingExponent(#764)StirlingRewritten(#766)StirlingSubstituted(here)The part measuring caught
Substituting for the factorial alone leaves
ln(e^Stirling / x^x), and nothing here simplifiesln(e^A). Those expressions did not merely stay unanswered — the descent assembled-oo/+oo, readNaNoff it, and three limits that had been slow non-answers became wrong answers:The fix is to broaden #766's rewrite from
ln(f!)to any logarithm holding a diverging factorial, so those go through the logarithm's expansion instead of this one. The sweep now reports noNaNat all. This is the reason the two changes are in one PR: the substitution is not safe to land without the broadening.Measured
70 generated factorial expressions —
x!,(2*x)!,(x+1)!,x!/x^x,x!/e^x,x^x/x!,(x+1)!/x!under ten wrappers — against master:NaN20 timeouts and 9 unevaluated became values, and nothing was lost. A representative sample was checked numerically with mpmath at 50 digits up to
x = 1e7— includingln(x!/x^x)/(x*ln(x)), which is-0.062atx = 1e7and tends to0slowly enough to be easy to get wrong by eye.casbench113/117, 0 wrong, 0 error, 0 timeout — unchangedpropcheck0 failures,rootcheck596/596,simpsweep0 disagreements of 10463