Skip to content

The integrator's sampled-point checks evaluate in intervals, which read no setting - #1497

Open
Rafael-SOWNet wants to merge 7 commits into
masterfrom
interval-evaluation
Open

Rafael-SOWNet wants to merge 7 commits into
masterfrom
interval-evaluation

Conversation

@Rafael-SOWNet

@Rafael-SOWNet Rafael-SOWNet commented Sep 27, 2026 •

Copy link
Copy Markdown
Member

Part of #1338, #1486 and #1019 (item 10).

The integrator checks a candidate answer numerically: it differentiates the answer, then compares the derivative with the integrand at a few sampled points. HoldsAtSampledPoints is that comparison, and about thirty integration routes call it. It used to evaluate both sides with EvalNumerical at a hundred digits, under DowncastingEnabled switched off. It now uses an internal interval evaluation that reads no global setting:

  • AngouriMath.Numerics.IntervalEvaluation: double intervals with outward rounding. An interval holds the exact value of what was evaluated, so every comparison has three outcomes: the values agree to the tolerance asked, they differ by more, or the intervals are too wide to tell.
  • AngouriMath.Numerics.PreciseEvaluation: the same evaluation in decimal intervals. The precision is a constructor argument. It is used at 40 and then 80 digits where cancellation leaves the doubles too wide.

Both are internal until v3, and EvalNumerical is unchanged beside them. They read numbers, pi, e, arithmetic, powers, exp/ln/log, the six trigonometric functions, arcsin, arccos, arctan, arccot, arcsec, arccsc (on real arguments), abs, sgn, provided and piecewise. Conditions are decided three ways: the comparisons, the connectives, and membership in RR and CC. Anything else is not read, and a point where either side isn't read is skipped, as a NaN was before.

One difference from what I described on #1486: the second tier's bounds are PeterO decimals, not the binary floats (a BigInteger mantissa and an int exponent) I proposed there. With PeterO's floor and ceiling contexts and the library's decimal functions, the second tier could be a direct port of the double one, and its verdicts could be checked against the old evaluation before anything moved. It reads no setting, and callers see only Of and Agree. The binary-float tier replaces it in the next PR, using the #1364/#1366 BigInteger series directly.

IsZeroAtPinnedSymbols, the zero test the partial-fraction code runs on symbolic coefficients, moves onto the same evaluation.

The one setting left

DerivativeHoldsAtSampledPoints pins the symbols to numbers and then simplifies the answer before differentiating it. That simplification folds number nodes. With downcasting on, a pinned 1.37 is 137/100, and simplifying the answer to u^5/((u^2 - 1)^5 (2 a u + b (u^2 + 1))) took 28 s instead of 3 ms. The downcasting scope therefore stays around that step. It belongs to the legacy number tower, not to the evaluation, and goes when number nodes stop reading the setting.

Audit

Before the switch, the pilot computed both verdicts at every check and returned the old one, so no answer moved:

checks same verdict old evaluation wrong
suite 659 659 0
Rubi 4.7.2, 4.3.10, 4.4.10 266 263 3
family 1 (sample) 651 649 2
family 4 (sample) 1203 1202 1

In every disagreement the intervals accepted an answer that the old evaluation rejected, and in every one the intervals were right. I checked five with exact rational arithmetic in Python, and the sixth, which has a square root in it, in 80-digit mpmath. The old evaluation's failures come from the limits of its decimals (#1498). Their exponents are bounded to [-100, 1000], so a decimal below about 1e-199 is zero and one above 1e1000 is infinite. In the five I bisected, a factor near 1e-316 times one near 1e+316 came out 0, and a factor flushed to zero times one saturated to infinity came out NaN. I didn't bisect the sixth.

The audit covered those four sets only, and that turned out to matter. The first measurement over Rubi 7.1.4/7.1.5 lost five problems (one unsolved, four timeouts) and changed five answers. The pinned answer's derivative carried conditions like -4 - 4 a^2 in RR, which OrderedCondition writes whenever a rule decides comparisons over operands that might be complex (#876). The evaluator had no case for set membership, so the whole value was undefined at every point, and correct answers were rejected. Membership in RR and CC is decided now. With it, 7.1.4/7.1.5 match master answer for answer.

Cost

The first version of the decimal tier made the check much slower than the hundred-digit evaluation on Rubi's hyperbolic files. On 1/(a + b csch(c + d x)^2)^3 the doubles were too wide at all twenty points of its five checks, because the answer's coefficients are large polynomials in the pinned values and cancel. Those checks spent 3.9 s, where the old evaluation spent 0.33 s. Four changes brought it level:

  • the contexts are built once per precision, because the library's constant cache is keyed by the context instance and was recomputing pi and the logarithm tables on every evaluation;
  • one evaluation serves every point of a check and remembers each node's value, so the subtrees free of the point are worked out once;
  • the exponential and logarithm use the library's fixed-point series instead of PeterO's, and monotone functions are evaluated at the two ends;
  • a product or quotient takes two directed operations instead of eight wherever neither operand straddles zero.

Measured

Against master 856831cf on the same machine, each set run on master and then on the branch back to back, with intbench's answers dumped on both sides and compared as text. Time and allocation are over the problems both solve.

problems solved verdicts changed answers changed time allocation
Rubi 4.7.2, 4.3.10, 4.4.10 312 269 → 269 0 0 −18.2% −32.2%
the 1774-problem suite 1774 1719 → 1719 0 0 −0.2% −1.3%
Rubi 7.1.4 – 7.2.5 376 345 → 345 0 0 −2.9% −1.8%
family 1 (6 a file) 228 170 → 169 1 0 +6.5%¹ +0.3%
family 4 (6 a file) 422 332 → 332 0 0 −1.0% −1.9%
family 5 (15 a file) 257 228 → 228 0 0 −0.3% +0.4%
family 6 (20 a file) 417 386 → 386 1 1 −5.6%² +0.2%
family 7 (15 a file) 270 240 → 240 0 0 −0.2% −1.1%

¹ Remeasured on a quiet machine, and it reproduces. Almost all of it is the six problems that follow 1.1.2.4's timeout. They run in the fresh worker the harness starts after a timeout: the first goes from 71 ms to 471 ms, and by the seventh it's 278 ms against 285. Without those six, family 1 is +1%.

² My own probe runs overlapped the master half, so this time isn't a clean measurement.

The three 4.7.2 problems whose checks returned NaN on master are 4–6 times faster: cos(x)^2 sin(x)^2/(a cos(x) + b sin(x))^2 goes from 11.4 s to 1.7 s.

Family 1's verdict is (c + d x^2)^(5/2)/(x^2 (a + b x^2)^2) (1.1.2.4), a timeout at the 5 s budget where master solved it in 5.3 s. The branch's check accepts a candidate that the old one rejected because a factor flushed to zero. The route after that candidate is slower, and it gives a correct answer in 11 s: its derivative matches the integrand at every point in 60-digit intervals. The old evaluation disagrees with that derivative at every point, so the probe marks it wrong. Family 6's verdict is a timeout on master that declines as unsolved on the branch.

The unit suite passes on net10.0: 12767 passed, 14 skipped, none failed, out of 12781, of which 54 are the evaluator's new tests. The allocation gate passes on all 19 gated benchmarks.

Changed answers

One, on Rubi's suites. e^(c (a + b x)) (cosh(a c + b c x)^2)^(3/2) (6.2.5) is now answered by the first route that reaches it, in tanh(a c + b c x). On master the old check turned that route away, because its derivative evaluated to NaN at three of the four points. At the point I bisected, the evaluation divided 2.40e-45 by 2.23e-91, and a divisor below PrecisionErrorZeroRange (1e-16) counts as zero. A later route then answered in exponentials. The intervals agree with the integrand at all four points, the old evaluation agrees at the fourth, and intbench verifies the new answer. 2.5.0 leaves it unevaluated. The row is in BREAKING-CHANGES.md.

One rule deliberately follows the old evaluation even though it isn't interval arithmetic: an exact zero times anything is zero. The old evaluation simplifies 0 x to 0 before evaluating x, and an answer's derivative can carry a term whose coefficient is exactly zero at the pinned values, beside a factor with no value there. Reading that product as undefined lost three answers in Rubi 6.x (see the commit that reverts it).

🤖 Generated with Claude Code

https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura

Rafael-SOWNet and others added 7 commits September 26, 2026 23:12
…ad no setting

HoldsAtSampledPoints and IsZeroAtPinnedSymbols compared hundred-digit decimals from
EvalNumerical under DowncastingEnabled off. They now evaluate in intervals: doubles with
outward rounding first, then decimal intervals of 40 and 80 digits where cancellation
leaves the doubles too wide to tell. The accuracy is an argument and no global setting is
read (#1019, item 10). The evaluator is internal, in AngouriMath.Numerics.

The one scope left is around the simplification of the pinned answer in
DerivativeHoldsAtSampledPoints, which folds number nodes: downcast, a pinned 1.37 is
137/100 and the answer to u^5/((u^2 - 1)^5 (2 a u + b (u^2 + 1))) took 28 s to simplify,
against 3 ms in decimals.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura
…each function's own error

The three are read on a real argument as the decimal evaluation defines them: arctan(1/x)
with pi/2 at zero, arccos(1/x) and arcsin(1/x). Without them, a check over an integrand
holding arccot could compare no point, and would turn every answer away.

In doubles, the sine and cosine widen the function's value before adding the reach, and
the complex logarithm's argument takes the same slack as every other function. In decimals,
the sine, cosine and inverse functions carry the absolute error of the fixed point they are
computed in, the sine's growing with the argument. A whole power is read directly only in
the range of an int: 2^(-2^63) negated its exponent into itself.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura
A value is real where its imaginary part is exactly zero, which real arithmetic keeps it,
and not real where that part is certainly not zero. The pinned answer's derivative carries
conditions of this kind on its constants, `-4 - 4 a^2 in RR`. Left undecided, they left the
whole value undecided at every point, and the check rejected a correct answer:
e^arsinh(a + b x)/x^3 came back unsolved, and four more of Rubi 7.1.4 and 7.1.5 timed out.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura
…the old one did

On Rubi 6.6.7 the decimal tier made the check about ten times slower than the hundred-digit
evaluation it replaced: 1/(a + b csch(c + d x)^2)^3 spent 3.9 s in five checks where the old
evaluation spent 0.33 s. Four things, measured together to 4.9 s against 4.5 s for the whole
integral:

- The contexts are made once per precision. The library's constant cache and guard-digit
  contexts are keyed by the context instance, so pi and the logarithm tables were recomputed
  for every evaluation.
- One evaluation serves every point of a check, remembering each node's value, so the
  subtrees free of the point -- which the substitutions share -- are worked out once.
- The exponential and the logarithm are the library's fixed-point series, and the monotone
  functions are evaluated at the two ends rather than at the middle with a slope, which cost
  the exponential a second evaluation.
- A product or quotient takes two directed operations where neither operand straddles zero,
  rather than eight.

And an undefined operand stops the evaluation, and makes a product undefined even beside an
exact zero: 0 * (1/0) was exactly 0.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura
…zero is one of them

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura
…actor

The decimal evaluation simplifies 0 x to 0 before it evaluates x, and an answer's derivative
can carry a term whose coefficient is exactly zero at the pinned values beside a factor with
no value there. Reading such a product as undefined, as the previous commit did, skipped
points the old check compared: Rubi 6.6.7's 1/(a + b csch(c + d x)^2) came back unsolved,
6.4.7's (a + b coth(x)^2)^(3/2) tanh(x)^2 timed out, and 6.3.7's tanh(c + d x)^2/(a + b
tanh(c + d x)^2) got an answer forty times as long. With the rule reverted all three are
answered as on master. A product's operands are both evaluated again; a sum's,
difference's and quotient's still stop at an undefined first operand.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura
@Happypig375

Copy link
Copy Markdown
Member

Would yield returns fit better than array of lazy here?

@Rafael-SOWNet

Copy link
Copy Markdown
Member Author

In HoldsAtSampledPoints the instances have to outlive a single point. Each one's node memo is what lets the four points of a check share the subtrees that don't depend on the point (the second change under Cost; I measured the four changes together, not that one alone). An iterator would build a fresh evaluation, with an empty memo, each time a point enumerates it. The Lazy array builds each precision at most once per check, and only when the doubles can't tell. A caching iterator would do the same, with a list field that the loop fills as it goes. I'm happy to switch if you'd rather read it that way.

In IsZeroAtPinnedSymbols each array is enumerated once, so an iterator is equivalent there.

This branch has not been deployed

No deployments
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants