smart_integrate() clips Piecewise integrals to their breakpoints
Summary
smart_integrate() replaces the caller's integration interval with the minimum and
maximum Piecewise breakpoints. A Piecewise expression with one breakpoint therefore
gets an empty interval and always integrates to zero.
Explanation
After get_bounds() finds discontinuities, lines 190-194 of
src/jaff/common/_integrators.py use those points as new integration limits instead
of using them only as hints to scipy.integrate.quad.
This minimal reproduction returns 0.0, although the exact integral is 3.0:
import sympy as sp
from jaff.common._integrators import smart_integrate
x = sp.Symbol("x")
expr = sp.Piecewise((1, x < 1), (2, True))
assert smart_integrate(expr, x, (0, 2)) == 3.0
The same helper is used to partition radiation source terms across bands, so a
Piecewise dRad can silently produce zero or omit the portions outside the outermost
breakpoints.
Proposed patch
Keep a, b = t_low, t_high for the entire requested interval. Treat values returned
by get_bounds() only as internal quad(..., points=...) breakpoints after filtering
them with a < point < b. Add a regression test for both a single-breakpoint
Piecewise expression and a multi-breakpoint expression whose requested bounds extend
beyond the first and last breakpoint.
smart_integrate()clips Piecewise integrals to their breakpointsSummary
smart_integrate()replaces the caller's integration interval with the minimum andmaximum Piecewise breakpoints. A Piecewise expression with one breakpoint therefore
gets an empty interval and always integrates to zero.
Explanation
After
get_bounds()finds discontinuities, lines 190-194 ofsrc/jaff/common/_integrators.pyuse those points as new integration limits insteadof using them only as hints to
scipy.integrate.quad.This minimal reproduction returns
0.0, although the exact integral is3.0:The same helper is used to partition radiation source terms across bands, so a
Piecewise
dRadcan silently produce zero or omit the portions outside the outermostbreakpoints.
Proposed patch
Keep
a, b = t_low, t_highfor the entire requested interval. Treat values returnedby
get_bounds()only as internalquad(..., points=...)breakpoints after filteringthem with
a < point < b. Add a regression test for both a single-breakpointPiecewise expression and a multi-breakpoint expression whose requested bounds extend
beyond the first and last breakpoint.