From 31f78de9529c352e970e27a096d081348e7ff42b Mon Sep 17 00:00:00 2001 From: KingArth0r Date: Sun, 4 Oct 2026 23:13:15 -0500 Subject: [PATCH] feat: integrate piecewise functions over finite real intervals --- src/compute-engine/library/calculus.ts | 23 +- .../symbolic/piecewise-integrate.ts | 483 ++++++++++++++++++ test/compute-engine/compile-integrate.test.ts | 18 +- .../integrate-piecewise.test.ts | 431 ++++++++++++++++ 4 files changed, 940 insertions(+), 15 deletions(-) create mode 100644 src/compute-engine/symbolic/piecewise-integrate.ts create mode 100644 test/compute-engine/integrate-piecewise.test.ts diff --git a/src/compute-engine/library/calculus.ts b/src/compute-engine/library/calculus.ts index c79eafca1..efd672a8c 100644 --- a/src/compute-engine/library/calculus.ts +++ b/src/compute-engine/library/calculus.ts @@ -79,6 +79,7 @@ import { // Self-registers the `expr.explain('D')` driver (see explain.ts) import '../symbolic/explain-derivative.js'; import { antiderivative } from '../symbolic/antiderivative.js'; +import { integratePiecewise } from '../symbolic/piecewise-integrate.js'; import { definiteIntegralByResidues, divergentIntegralValue, @@ -3318,7 +3319,9 @@ volumes if ( isSymbol(op) && opDef !== undefined && - (!opDef.lazy || op.symbol === 'Add' || op.symbol === 'Multiply') && + (!opDef.lazy || + op.symbol === 'Add' || + op.symbol === 'Multiply') && orders.reduce((sum, n) => sum + n, 0) <= MAX_NESTED_PARTIAL_ORDER ) { const names = orders.map((_, i) => `_${i + 1}`); @@ -4537,14 +4540,16 @@ volumes continue; } if (isDefinite) { - const split = integrateAcrossKinks( - ce, - integrand, - variable, - lower, - upper, - numericApproximation ?? false - ); + const split = + integratePiecewise(ce, integrand, variable, lower, upper) ?? + integrateAcrossKinks( + ce, + integrand, + variable, + lower, + upper, + numericApproximation ?? false + ); if (split !== undefined) { isIndefinite = false; expr = diff --git a/src/compute-engine/symbolic/piecewise-integrate.ts b/src/compute-engine/symbolic/piecewise-integrate.ts new file mode 100644 index 000000000..5a5b230f1 --- /dev/null +++ b/src/compute-engine/symbolic/piecewise-integrate.ts @@ -0,0 +1,483 @@ +import type { Expression, IComputeEngine } from '../global-types.js'; +import { isFunction, isNumber, sym } from '../boxed-expression/type-guards.js'; +import { boundVariableNamesInOperand } from '../boxed-expression/binders.js'; +import { getPolynomialCoefficients } from '../boxed-expression/polynomials.js'; +import { solveOverDomain } from '../boxed-expression/solve-domain.js'; +import { isExactConstantExpression } from '../boxed-expression/compare.js'; +import { differentiate } from './derivative.js'; +import { exactSign } from './contour.js'; +import { checkDeadline } from '../../common/interruptible.js'; +import { shadowsLibraryName } from '../library-shadowing.js'; + +const PIECEWISE = new Set([ + 'If', + 'Which', + 'Min', + 'Max', + 'Floor', + 'Ceil', + 'Fract', +]); +const COMPARISONS = new Set([ + 'Less', + 'LessEqual', + 'Greater', + 'GreaterEqual', + 'Equal', + 'NotEqual', +]); +const MAX_CUTS = 128; +const MAX_NODES = 4096; +const MAX_SOLVES = 128; + +function finiteReal(e: Expression): boolean { + return ( + e.isValid && + e.isNaN !== true && + e.isFinite !== false && + e.type.matches('real') + ); +} + +function sign(e: Expression): -1 | 0 | 1 | undefined { + // The shared exact resolver must separate close constants, not identify + // overlapping enclosures as equal at its precision limit. + if (isExactConstantExpression(e) || (isNumber(e) && e.isExact)) + return exactSign(e); + const value = e.evaluate(); + if (value.isSame(0)) return 0; + if (value.isPositive === true) return 1; + if (value.isNegative === true) return -1; + return undefined; +} + +function binds(e: Expression, index: number, variable: string): boolean { + return boundVariableNamesInOperand(e, index).includes(variable); +} + +/** Integrate finite real intervals by resolving branches on open cells. + * `undefined` means no relevant piecewise expression; `inert` means the complete + * partition or a cell integral could not be established. Point values at cuts + * do not affect these ordinary integrals, but poles still reach Integrate's + * one-sided endpoint checks when each resolved cell is evaluated. */ +export function integratePiecewise( + ce: IComputeEngine, + integrand: Expression, + variable: string, + lower: Expression, + upper: Expression +): Expression | 'inert' | undefined { + if (!integrand.has([...PIECEWISE])) return undefined; + let nodes = 0; + let exhausted = false; + const contains = (e: Expression): boolean => { + if (++nodes > MAX_NODES) { + exhausted = true; + return false; + } + if (!isFunction(e)) return false; + if (PIECEWISE.has(e.operator) && e.has(variable)) return true; + return e.ops.some((op, i) => !binds(e, i, variable) && contains(op)); + }; + if (!contains(integrand)) return exhausted ? 'inert' : undefined; + // Partitioning and resolving must observe the same function and bounds. + if (!integrand.isPure || !lower.isPure || !upper.isPure) return 'inert'; + const a = lower.evaluate(); + const b = upper.evaluate(); + if (!finiteReal(a) || !finiteReal(b)) return 'inert'; + const direction = sign(b.sub(a)); + if (direction === undefined) return 'inert'; + if (direction === 0) return ce.Zero; + const lo = direction > 0 ? a : b; + const hi = direction > 0 ? b : a; + const cuts: Expression[] = []; + const domainChecks = new Set(); + const rootCache = new Map(); + let solves = 0; + + const affine = (e: Expression): [Expression, Expression] | undefined => { + if (!e.has(variable)) { + const constant = e.evaluate(); + return finiteReal(constant) ? [constant, ce.Zero] : undefined; + } + const coefficients = getPolynomialCoefficients(e, variable); + if (!coefficients || coefficients.length > 2) return undefined; + const m = (coefficients[0] ?? ce.Zero).evaluate(); + const k = (coefficients[1] ?? ce.Zero).evaluate(); + return finiteReal(m) && finiteReal(k) ? [m, k] : undefined; + }; + + const addCut = (cut: Expression): boolean => { + const left = sign(cut.sub(lo)); + const right = sign(hi.sub(cut)); + if (left === undefined || right === undefined) return false; + if (left <= 0 || right <= 0) return true; + for (let i = 0; i < cuts.length; i++) { + const order = sign(cut.sub(cuts[i])); + if (order === undefined) return false; + if (order === 0) return true; + if (order < 0) { + if (cuts.length >= MAX_CUTS) return false; + cuts.splice(i, 0, cut); + return true; + } + } + if (cuts.length >= MAX_CUTS) return false; + cuts.push(cut); + return true; + }; + + const roots = (e: Expression): readonly Expression[] | undefined => { + if (++nodes > MAX_NODES) return undefined; + const coefficients = affine(e); + if (coefficients) { + const [m, k] = coefficients; + const slope = sign(k); + if (slope === undefined) return undefined; + return slope === 0 ? [] : [m.neg().div(k).evaluate()]; + } + if (isFunction(e, 'Divide')) return roots(e.op1); + if (isFunction(e, 'Power') && !e.op2.has(variable)) { + const exponent = sign(e.op2); + if (exponent === -1) return []; + if (exponent === 1) return roots(e.op1); + } + if ( + isFunction(e, 'Exp') || + (isFunction(e, 'Power') && !e.op1.has(variable) && sign(e.op1) === 1) + ) + return []; + if (isFunction(e, 'Multiply')) { + const result: Expression[] = []; + for (const op of e.ops) { + const factor = roots(op); + if (factor === undefined) return undefined; + result.push(...factor); + } + return result; + } + if (rootCache.has(e)) return rootCache.get(e); + if (++solves > MAX_SOLVES) return undefined; + checkDeadline(ce._deadlineFrame); + // A bounded domain expands trig root families and rejects partial answers. + // Keep the exact returned expressions, never numerical root approximations. + const found = solveOverDomain(ce, e.canonical, { + unknown: variable, + domain: ce.function('Interval', [lo, hi]), + }); + const result = found?.every( + (r) => + finiteReal(r) && + (isNumber(r) ? r.isExact : isExactConstantExpression(r)) + ) + ? found + : undefined; + rootCache.set(e, result); + return result; + }; + + const zero = (e: Expression): boolean => { + const found = roots(e); + return found !== undefined && found.every(addCut); + }; + + const integerPart = ( + value: Expression, + ceiling = false + ): Expression | undefined => { + const n = ce.function(ceiling ? 'Ceil' : 'Floor', [value]).evaluate(); + if (!isNumber(n) || !n.isExact || n.isInteger !== true) return undefined; + // Floor/Ceil may resolve an enclosure at their precision limit as a tie. + // Partitioning requires strict separation from the adjacent integer. + const lower = sign(value.sub(ceiling ? n.sub(ce.One) : n)); + const upper = sign((ceiling ? n : n.add(1)).sub(value)); + if (lower === undefined || upper === undefined) return undefined; + return (ceiling ? lower > 0 && upper >= 0 : lower >= 0 && upper > 0) + ? n + : undefined; + }; + + // On each open cell these operations are continuous and real, provided their + // domain checks hold. Poles and branch-domain endpoints are cuts too: zeros + // alone do not suffice to classify a rational inequality or a logarithm. + const continuous = (e: Expression): boolean => { + if (++nodes > MAX_NODES) return false; + if (sym(e) === variable) return true; + if (!e.has(variable)) return finiteReal(e.evaluate()); + if (!isFunction(e)) return false; + if (shadowsLibraryName(ce, e.operator)) return false; + if (!e.ops.every(continuous)) return false; + const op = e.operator; + if ( + [ + 'Add', + 'Subtract', + 'Multiply', + 'Negate', + 'Square', + 'Sin', + 'Cos', + 'Exp', + 'Sinh', + 'Cosh', + 'Tanh', + 'Arctan', + ].includes(op) + ) + return true; + domainChecks.add(e); + if (op === 'Divide') return zero(e.op2); + if (op === 'Ln' || op === 'Sqrt') return zero(e.op1); + if (op === 'Tan' || op === 'Sec') return zero(ce.function('Cos', [e.op1])); + if (op === 'Cot' || op === 'Csc') return zero(ce.function('Sin', [e.op1])); + if (op === 'Power') { + const exponent = e.op2.evaluate(); + if ( + !e.op2.has(variable) && + isNumber(exponent) && + exponent.isInteger && + exponent.isPositive + ) + return true; + // A fixed positive base admits real exponents everywhere. Other powers + // can change their real domain only where the base vanishes. + if (!e.op1.has(variable) && sign(e.op1) === 1) return true; + if (!e.op2.has(variable)) return zero(e.op1); + } + return false; + }; + + // Sampling chooses a branch only after all roots and domain boundaries of + // its supported continuous operands have been accounted for. + const condition = (e: Expression): boolean => { + if (++nodes > MAX_NODES) return false; + if (!e.has(variable)) { + const value = sym(e.evaluate()); + return value === 'True' || value === 'False'; + } + if (!isFunction(e)) return false; + if (shadowsLibraryName(ce, e.operator)) return false; + if (e.operator === 'Not' && e.nops === 1) return condition(e.op1); + if (e.operator === 'And' || e.operator === 'Or') + return e.ops.every(condition); + return ( + COMPARISONS.has(e.operator) && + e.nops === 2 && + continuous(e.op1) && + continuous(e.op2) && + zero(e.op1.sub(e.op2)) + ); + }; + + const collect = (e: Expression): boolean => { + if (++nodes > MAX_NODES) return false; + if (!isFunction(e) || !e.has(variable)) return true; + if (shadowsLibraryName(ce, e.operator)) return false; + if (e.operator === 'If' || e.operator === 'Which') { + if (e.operator === 'If') { + if (e.nops < 2 || e.nops > 3 || !condition(e.op1)) return false; + const truth = !e.op1.has(variable) ? sym(e.op1.evaluate()) : undefined; + if (truth === 'True') return collect(e.op2); + if (truth === 'False') return e.nops === 3 && collect(e.op3); + return e.ops.slice(1).every(collect); + } + if (e.nops % 2 !== 0) return false; + for (let i = 0; i < e.nops; i += 2) { + const predicate = e.ops[i]; + if (!condition(predicate)) return false; + const truth = !predicate.has(variable) + ? sym(predicate.evaluate()) + : undefined; + if (truth === 'False') continue; + if (!collect(e.ops[i + 1])) return false; + if (truth === 'True') break; + } + return true; + } + if (e.operator === 'Min' || e.operator === 'Max') { + if (e.nops === 0 || !e.ops.every(continuous)) return false; + for (let i = 0; i < e.nops; i++) + for (let j = 0; j < i; j++) { + if (++nodes > MAX_NODES || !zero(e.ops[i].sub(e.ops[j]))) + return false; + } + return true; + } + if ( + e.operator === 'Floor' || + e.operator === 'Ceil' || + e.operator === 'Fract' + ) { + if (e.nops !== 1) return false; + const coefficients = affine(e.op1); + if (!coefficients) { + if (!continuous(e.op1)) return false; + const d = differentiate(e.op1.canonical, variable); + if (d === undefined) return false; + const critical = roots(d.evaluate()); + if (critical === undefined) return false; + // All extrema of a smooth branch occur at endpoints or critical + // points. Domain cuts are included; a pole there makes the range + // unbounded and prevents a finite integer-level partition. + const values = [lo, hi, ...cuts, ...critical].map((p) => + e.op1.subs({ [variable]: p }).evaluate() + ); + if (!values.every(finiteReal)) return false; + const levels = values.map((v) => integerPart(v)); + if (!levels.every((v) => v !== undefined && Number.isSafeInteger(v.re))) + return false; + const start = Math.min(...levels.map((v) => v!.re)); + const end = Math.max(...levels.map((v) => v!.re)); + if (end - start > MAX_CUTS) return false; + // Include the lowest integer too: an attained minimum may be a jump + // of Ceil, even though Floor takes its usual value there. + for (let n = start; n <= end; n++) + if (!zero(e.op1.sub(ce.number(n)))) return false; + return true; + } + const [m, k] = coefficients; + const slope = sign(k); + if (slope === undefined) return false; + if (slope === 0) return true; + const first = integerPart(m.add(k.mul(lo))); + const last = integerPart(m.add(k.mul(hi))); + if ( + first === undefined || + last === undefined || + !Number.isSafeInteger(first.re) || + !Number.isSafeInteger(last.re) + ) + return false; + const start = Math.min(first.re, last.re); + const end = Math.max(first.re, last.re); + if (end - start > MAX_CUTS) return false; + for (let n = start + 1; n <= end; n++) + if (!addCut(ce.number(n).sub(m).div(k).evaluate())) return false; + return true; + } + if ((e.operator === 'Abs' || e.operator === 'Sign') && e.nops === 1) + return continuous(e.op1) && zero(e.op1); + return e.ops.every((op, i) => binds(e, i, variable) || collect(op)); + }; + if (!collect(integrand)) return 'inert'; + + const truthAt = (e: Expression, sample: Expression): boolean | undefined => { + if (!e.has(variable)) { + const value = sym(e.evaluate()); + return value === 'True' ? true : value === 'False' ? false : undefined; + } + if (!isFunction(e)) return undefined; + if (e.operator === 'Not') { + const truth = truthAt(e.op1, sample); + return truth === undefined ? undefined : !truth; + } + if (e.operator === 'And' || e.operator === 'Or') { + const values = e.ops.map((op) => truthAt(op, sample)); + const short = e.operator === 'Or'; + if (values.includes(short)) return short; + return values.includes(undefined) ? undefined : !short; + } + const order = sign(e.op1.sub(e.op2).subs({ [variable]: sample })); + if (order === undefined) return undefined; + switch (e.operator) { + case 'Less': + return order < 0; + case 'LessEqual': + return order <= 0; + case 'Greater': + return order > 0; + case 'GreaterEqual': + return order >= 0; + case 'Equal': + return order === 0; + case 'NotEqual': + return order !== 0; + default: + return undefined; + } + }; + + const resolve = ( + e: Expression, + sample: Expression + ): Expression | undefined => { + if (!isFunction(e) || !e.has(variable)) return e; + const at = (op: Expression) => op.subs({ [variable]: sample }).evaluate(); + if (e.operator === 'If') { + const truth = truthAt(e.op1, sample); + if (truth === true) return resolve(e.op2, sample); + if (truth === false && e.nops === 3) return resolve(e.op3, sample); + return undefined; + } + if (e.operator === 'Which') { + for (let i = 0; i < e.nops; i += 2) { + const truth = truthAt(e.ops[i], sample); + if (truth === true) return resolve(e.ops[i + 1], sample); + if (truth !== false) return undefined; + } + return undefined; + } + if (e.operator === 'Min' || e.operator === 'Max') { + let best = e.op1; + for (const op of e.ops.slice(1)) { + const order = sign(at(op).sub(at(best))); + if (order === undefined) return undefined; + if ( + (e.operator === 'Min' && order < 0) || + (e.operator === 'Max' && order > 0) + ) + best = op; + } + return best; + } + if ( + e.operator === 'Floor' || + e.operator === 'Ceil' || + e.operator === 'Fract' + ) { + const integer = integerPart(at(e.op1), e.operator === 'Ceil'); + if (integer === undefined) return undefined; + return e.operator === 'Fract' ? e.op1.sub(integer) : integer; + } + if (e.operator === 'Abs' || e.operator === 'Sign') { + const s = sign(at(e.op1)); + if (s === undefined) return undefined; + return e.operator === 'Sign' ? ce.number(s) : s < 0 ? e.op1.neg() : e.op1; + } + const operands: Expression[] = []; + for (let i = 0; i < e.nops; i++) { + const op = binds(e, i, variable) ? e.ops[i] : resolve(e.ops[i], sample); + if (op === undefined) return undefined; + operands.push(op); + } + return operands.every((op, i) => op === e.ops[i]) + ? e + : ce.function(e.operator, operands); + }; + + const points = [lo, ...cuts, hi]; + const values: Expression[] = []; + for (let i = 0; i + 1 < points.length; i++) { + checkDeadline(ce._deadlineFrame); + const [left, right] = [points[i], points[i + 1]]; + const sample = left.add(right).div(2); + if ( + ![...domainChecks].every((e) => + finiteReal(e.subs({ [variable]: sample }).evaluate()) + ) + ) + return 'inert'; + const body = resolve(integrand, sample); + if (body === undefined) return 'inert'; + const value = ce + .function('Integrate', [ + body, + ce.function('Limits', [ce.symbol(variable), left, right]), + ]) + .evaluate(); + if (value.has('Integrate')) return 'inert'; + values.push(value); + } + const total = ce.function('Add', values).evaluate(); + return direction > 0 ? total : total.neg().evaluate(); +} diff --git a/test/compute-engine/compile-integrate.test.ts b/test/compute-engine/compile-integrate.test.ts index 183248590..07b24fba7 100644 --- a/test/compute-engine/compile-integrate.test.ts +++ b/test/compute-engine/compile-integrate.test.ts @@ -152,7 +152,7 @@ describe('COMPILE Integrate — adaptive Gauss–Kronrod', () => { expect(Math.abs(got - 3.7408935992e-6)).toBeLessThan(1e-6); }); - test('piecewise integrand ∫_0^2 f, f = 1 for t<1 else 2 → 3', () => { + test('affine piecewise integral resolves exactly before compilation', () => { // Built via Which (jump discontinuity at t = 1). const expr = ce.box([ 'Integrate', @@ -161,18 +161,24 @@ describe('COMPILE Integrate — adaptive Gauss–Kronrod', () => { ]); const r = compile(expr); expect(r.success).toBe(true); - expect(r.code).toContain('_SYS.integrate('); + expect(r.code).not.toContain('_SYS.integrate'); + expect(r.code).not.toContain('integrateMC'); expect(r.run() as number).toBeCloseTo(3, 6); }); describe('quadrature option', () => { - // A piecewise (Which) integrand has no elementary antiderivative, so the - // antiderivative-first path declines and the quadrature emitter is exercised - // (∫_0^2 of {1 for t<1, else 2} = 3). + // The non-periodic trig argument is outside the symbolic splitter, so this + // exercises quadrature. On [0, 2], sin(x²-1) < 0 selects [0, 1). const piecewise = () => ce.box([ 'Integrate', - ['Which', ['Less', 'x', 1], 1, 'True', 2], + [ + 'Which', + ['Less', ['Sin', ['Subtract', ['Square', 'x'], 1]], 0], + 1, + 'True', + 2, + ], ['Limits', 'x', 0, 2], ]); diff --git a/test/compute-engine/integrate-piecewise.test.ts b/test/compute-engine/integrate-piecewise.test.ts new file mode 100644 index 000000000..7a89832dc --- /dev/null +++ b/test/compute-engine/integrate-piecewise.test.ts @@ -0,0 +1,431 @@ +import { ComputeEngine } from '../../src/compute-engine'; +import type { ExpressionInput } from '../../src/compute-engine/global-types'; + +function integrate( + body: ExpressionInput, + lower: ExpressionInput = 0, + upper: ExpressionInput = 1, + ce = new ComputeEngine() +) { + return ce.box(['Integrate', body, ['Limits', 'x', lower, upper]]).evaluate(); +} + +describe('FINITE REAL PIECEWISE INTEGRALS', () => { + test.each([ + [['Min', 'x', ['Subtract', 1, 'x']], 0, 1, 1 / 4], + [['Max', 'x', ['Subtract', 1, 'x']], 0, 1, 3 / 4], + [['Min', 'x', ['Negate', 'x'], 1], -1, 1, -1], + [['Floor', 'x'], 0, 3, 3], + [['Ceil', 'x'], 0, 3, 6], + [['Fract', 'x'], 0, 3, 3 / 2], + [['Floor', 'x'], -1.5, 1.5, -1.5], + [['Ceil', 'x'], -1.5, 1.5, 1.5], + [['Fract', 'x'], -1.5, 1.5, 1.5], + [['Floor', ['Subtract', 1, ['Multiply', 2, 'x']]], 0, 1, -1 / 2], + [['Floor', ['Add', ['Multiply', 2, 'x'], ['Rational', 1, 2]]], 0, 1, 1], + [['Multiply', 'x', ['Floor', 'x']], 0, 3, 13 / 2], + [['Floor', 'x'], 3, 0, -3], + [['Max', 'x', 0], 1, -1, -1 / 2], + ] as [ExpressionInput, number, number, number][])( + '%j on [%s, %s]', + (body, a, b, expected) => { + const result = integrate(body, a, b); + expect(result.has('Integrate')).toBe(false); + expect(result.N().re).toBeCloseTo(expected, 12); + } + ); + + test('exact results and exact rational boundaries', () => { + expect(integrate(['Min', 'x', ['Subtract', 1, 'x']]).json).toEqual([ + 'Rational', + 1, + 4, + ]); + expect( + integrate(['Floor', 'x'], ['Rational', -3, 2], ['Rational', 3, 2]).json + ).toEqual(['Rational', -3, 2]); + }); + + test('If resolves open cells, including reversed bounds', () => { + const f: ExpressionInput = ['If', ['Less', 'x', 0], ['Negate', 'x'], 'x']; + expect(integrate(f, -1, 1).json).toBe(1); + expect(integrate(f, 1, -1).json).toBe(-1); + }); + + test('Which keeps first-true ordering for overlapping conditions', () => { + const f: ExpressionInput = [ + 'Which', + ['Less', 'x', 1], + 2, + ['Less', 'x', 2], + 5, + 'True', + 9, + ]; + expect(integrate(f, 0, 3).json).toBe(16); + }); + + test('Boolean combinations of affine conditions', () => { + const f: ExpressionInput = [ + 'If', + ['And', ['GreaterEqual', 'x', 0], ['Less', 'x', 1]], + 'x', + 0, + ]; + expect(integrate(f, -2, 2).json).toEqual(['Rational', 1, 2]); + expect( + integrate(['If', ['Not', ['Equal', 'x', 0]], 1, 100], -1, 1).json + ).toBe(2); + }); + + test('isolated exceptional values do not change the integral', () => { + expect(integrate(['If', ['Equal', 'x', 0], 999, 1], -1, 1).json).toBe(2); + expect( + integrate(['Which', ['Less', 'x', 0], -1, ['Greater', 'x', 0], 1], -1, 1) + .json + ).toBe(0); + }); + + test('missing branches on an open interval remain unevaluated', () => { + expect( + integrate(['Which', ['Less', 'x', 0], 1], -1, 1).has('Integrate') + ).toBe(true); + expect(integrate(['If', ['Less', 'x', 0], 1], -1, 1).has('Integrate')).toBe( + true + ); + }); + + test('nested selectors and weighted pieces', () => { + const f: ExpressionInput = [ + 'Multiply', + 'x', + ['If', ['Less', 'x', 1], ['If', ['Less', 'x', 0], -1, 1], 2], + ]; + expect(integrate(f, -1, 2).json).toBe(4); + }); + + test('LaTeX cases use the same branch semantics', () => { + const ce = new ComputeEngine(); + const result = ce.parse( + '\\int_{-1}^{1} \\begin{cases} x & x>0 \\\\ -x & \\text{otherwise} \\end{cases} \\, dx' + ); + expect(result.isValid).toBe(true); + expect(result.evaluate().json).toBe(1); + }); + + test('affine absolute values and signs share the partition', () => { + expect( + integrate(['Multiply', ['Floor', 'x'], ['Abs', 'x']], -1, 1).json + ).toEqual(['Rational', -1, 2]); + expect( + integrate(['Multiply', ['Ceil', 'x'], ['Sign', 'x']], -1, 1).json + ).toBe(1); + }); + + test('random branch boundaries are not sampled to construct a partition', () => { + expect( + integrate(['If', ['Less', 'x', ['Random']], 1, 2]).has('Integrate') + ).toBe(true); + }); + + test('a pole in an active branch is not hidden by splitting', () => { + const f: ExpressionInput = ['If', ['Less', 'x', 0], 0, ['Power', 'x', -2]]; + expect(integrate(f, -1, 1).isInfinity).toBe(true); + }); + + test('an inactive singular branch is not integrated', () => { + expect( + integrate(['If', ['Less', 'x', 0], ['Power', 'x', -2], 1], 0, 1).json + ).toBe(1); + }); + + test('unknown root positions and unsupported root families remain unevaluated', () => { + expect(integrate(['If', ['Less', 'x', 'a'], 1, 2]).has('Integrate')).toBe( + true + ); + expect( + integrate(['If', ['Less', ['Sin', ['Square', 'x']], 0], 1, 2]).has( + 'Integrate' + ) + ).toBe(true); + expect( + integrate(['Floor', ['Sin', ['Divide', 1, 'x']]], 0, 2).has('Integrate') + ).toBe(true); + }); + + test.each([ + [['If', ['Less', ['Square', 'x'], 2], 1, 0], -2, 2, 2 * Math.sqrt(2)], + [['Max', ['Square', 'x'], 'x'], -1, 2, 19 / 6], + [['Min', ['Square', 'x'], 'x'], -1, 2, 4 / 3], + [['If', ['Greater', ['Square', 'x'], 0], 1, 99], -1, 1, 2], + [['If', ['Less', ['Add', ['Square', 'x'], 1], 0], 1, 2], -1, 1, 4], + [['Floor', ['Square', 'x']], 0, 2, 5 - Math.sqrt(2) - Math.sqrt(3)], + [['Ceil', ['Square', 'x']], 0, 2, 7 - Math.sqrt(2) - Math.sqrt(3)], + [['Fract', ['Square', 'x']], 0, 2, Math.sqrt(2) + Math.sqrt(3) - 7 / 3], + [ + ['Floor', ['Square', 'x']], + -2, + 2, + 10 - 2 * Math.sqrt(2) - 2 * Math.sqrt(3), + ], + [['Floor', ['Subtract', 2, ['Square', 'x']]], -1, 1, 2], + [['If', ['Greater', ['Divide', 1, 'x'], 0], 1, 2], -1, 1, 3], + [['If', ['Less', ['Sqrt', 'x'], 1], 1, 2], 0, 4, 7], + [['If', ['Less', ['Exp', 'x'], 2], 1, 0], 0, 2, Math.log(2)], + [['If', ['Less', ['Ln', 'x'], 1], 1, 0], 1, 4, Math.E - 1], + ] as [ExpressionInput, number, number, number][])( + 'nonlinear switches: %j on [%s, %s]', + (body, a, b, expected) => { + const result = integrate(body, a, b); + expect(result.has('Integrate')).toBe(false); + expect(result.N().re).toBeCloseTo(expected, 10); + } + ); + + test('irrational switch positions remain exact', () => { + const ce = new ComputeEngine(); + const result = integrate( + ['If', ['Less', ['Square', 'x'], 2], 1, 0], + -2, + 2, + ce + ); + expect(result.isSame(ce.box(['Multiply', 2, ['Sqrt', 2]]).evaluate())).toBe( + true + ); + }); + + test('periodic sine switches include every period and retain Pi', () => { + const result = integrate(['If', ['Greater', ['Sin', 'x'], 0], 1, 0], 0, [ + 'Multiply', + 4, + 'Pi', + ]); + expect(result.has('Integrate')).toBe(false); + expect(result.has('Pi')).toBe(true); + expect(result.N().re).toBeCloseTo(2 * Math.PI, 10); + }); + + test('tangent switches include poles, not only zeros', () => { + const result = integrate(['If', ['Greater', ['Tan', 'x'], 0], 1, 2], 0, [ + 'Multiply', + 2, + 'Pi', + ]); + expect(result.has('Integrate')).toBe(false); + expect(result.N().re).toBeCloseTo(3 * Math.PI, 10); + }); + + test('periodic switches respect the engine angular unit', () => { + const ce = new ComputeEngine(); + ce.angularUnit = 'deg'; + expect( + integrate(['If', ['Greater', ['Sin', 'x'], 0], 1, 0], 0, 360, ce).json + ).toBe(180); + }); + + test('a user-defined Sin is not classified using sine root families', () => { + const ce = new ComputeEngine(); + ce.declare('Sin', { + signature: '(number) -> number', + evaluate: ([t], { engine }) => t.add(engine.One), + }); + expect( + integrate(['If', ['Greater', ['Sin', 'x'], 0], 1, 0], -2, 2, ce).has( + 'Integrate' + ) + ).toBe(true); + }); + + test('shifted periodic roots and tangent contacts are all included', () => { + const result = integrate( + ['If', ['Greater', ['Sin', 'x'], ['Rational', 1, 2]], 1, 0], + 0, + ['Multiply', 4, 'Pi'] + ); + expect(result.has('Integrate')).toBe(false); + expect(result.N().re).toBeCloseTo((4 * Math.PI) / 3, 10); + expect( + integrate( + ['If', ['LessEqual', ['Square', ['Subtract', 'x', 1]], 0], 99, 1], + 0, + 2 + ).json + ).toBe(2); + }); + + test('nonlinear floor partitions include interior extrema', () => { + expect( + integrate(['Floor', ['Subtract', 2, ['Square', 'x']]], -2, 2).N().re + ).toBeCloseTo(-6 + 2 * Math.sqrt(2) + 2 * Math.sqrt(3), 10); + const wave = integrate(['Floor', ['Multiply', 2, ['Sin', 'x']]], 0, [ + 'Multiply', + 2, + 'Pi', + ]); + expect(wave.has('Integrate')).toBe(false); + expect(wave.N().re).toBeCloseTo(-Math.PI, 10); + }); + + test('nearby polynomial crossings are kept distinct', () => { + const epsilon: ExpressionInput = ['Rational', 1, 1000000000000]; + const f: ExpressionInput = [ + 'Multiply', + ['Subtract', 'x', 1], + ['Subtract', 'x', ['Add', 1, epsilon]], + ]; + expect(integrate(['If', ['Less', f, 0], 1, 0], 0, 2).json).toEqual(epsilon); + }); + + test('exact irrational and transcendental constants retain tiny cells', () => { + const ce = new ComputeEngine({ precision: 15 }); + const epsilon: ExpressionInput = ['Power', 10, -30]; + for (const center of [['Sqrt', 2], 'Pi', ['Ln', 2]] as ExpressionInput[]) { + const f: ExpressionInput = [ + 'Which', + ['Less', 'x', center], + 0, + ['Less', 'x', ['Add', center, epsilon]], + 1, + 'True', + 0, + ]; + const result = integrate(f, 0, 4, ce); + expect(result.isSame(ce.box(epsilon).evaluate())).toBe(true); + } + expect(ce.precision).toBe(15); + }); + + test('real-domain holes and unbounded nonlinear integer levels are not hidden', () => { + expect( + integrate( + ['If', ['Greater', ['Ln', ['Square', 'x']], 0], 1, 0], + -2, + 2 + ).has('Integrate') + ).toBe(false); + expect(integrate(['Floor', ['Tan', 'x']], 0, 'Pi').has('Integrate')).toBe( + true + ); + expect( + integrate(['Floor', ['Divide', 1, ['Square', 'x']]], -1, 1).has( + 'Integrate' + ) + ).toBe(true); + }); + + test('non-real condition regions and incomplete root sets stay unevaluated', () => { + expect( + integrate(['If', ['Less', ['Sqrt', 'x'], 1], 1, 2], -1, 1).has( + 'Integrate' + ) + ).toBe(true); + const partial: ExpressionInput = [ + 'Multiply', + ['Subtract', 'x', 1], + ['Add', 'x', ['Exp', 'x']], + ]; + expect( + integrate(['If', ['Greater', partial, 0], 1, 0], -2, 2).has('Integrate') + ).toBe(true); + }); + + test('unbounded and excessive integer-cell partitions stay unevaluated', () => { + expect( + integrate(['Floor', 'x'], 0, 'PositiveInfinity').has('Integrate') + ).toBe(true); + expect(integrate(['Fract', 'x'], 0, 1000000).has('Integrate')).toBe(true); + expect(integrate(['Floor', ['Divide', 1, 'x']]).has('Integrate')).toBe( + true + ); + }); + + test('an assigned outer variable does not replace the bound variable', () => { + const ce = new ComputeEngine(); + ce.assign('x', 10); + expect(integrate(['Max', 'x', 0], -1, 1, ce).json).toEqual([ + 'Rational', + 1, + 2, + ]); + expect(ce.box('x').evaluate().json).toBe(10); + }); + + test('tiny distinct cuts are not merged by a numerical tolerance', () => { + const epsilon: ExpressionInput = ['Rational', 1, 1000000000000]; + const f: ExpressionInput = [ + 'Which', + ['Less', 'x', 0], + 0, + ['Less', 'x', epsilon], + 1, + 'True', + 0, + ]; + expect(integrate(f, -1, 1).json).toEqual(epsilon); + }); + + test('rounding at integer endpoints does not create a spurious cell', () => { + expect(integrate(['Floor', ['Negate', 'x']], 0, 1).json).toBe(-1); + expect(integrate(['Ceil', ['Negate', 'x']], 0, 1).json).toBe(0); + expect(integrate(['Fract', ['Negate', 'x']], 0, 1).json).toEqual([ + 'Rational', + 1, + 2, + ]); + }); + + test('open branches are evaluated even when their endpoints select another value', () => { + const f: ExpressionInput = [ + 'If', + ['Or', ['Equal', 'x', 0], ['Equal', 'x', 1]], + 99, + 'x', + ]; + expect(integrate(f).json).toEqual(['Rational', 1, 2]); + }); + + test('a pole at a branch boundary diverges with the correct orientation', () => { + const f: ExpressionInput = [ + 'Which', + ['Less', 'x', 1], + ['Power', ['Subtract', 'x', 1], -2], + 'True', + 0, + ]; + expect(integrate(f, 0, 2).json).toBe('PositiveInfinity'); + expect(integrate(f, 2, 0).json).toBe('NegativeInfinity'); + }); + + test('smooth cells are sent through the registered integration provider', () => { + const ce = new ComputeEngine(); + const visited: string[] = []; + ce._integrationProvider = (f) => { + visited.push(f.operator); + if (f.operator === 'Sin') return ce.box(['Negate', ['Cos', 'x']]); + if (f.operator === 'Cos') return ce.box(['Sin', 'x']); + return null; + }; + const value = integrate( + ['If', ['Less', 'x', 0], ['Sin', 'x'], ['Cos', 'x']], + -1, + 1, + ce + ); + expect(visited).toEqual(['Sin', 'Cos']); + expect(value.N().re).toBeCloseTo(Math.cos(1) - 1 + Math.sin(1), 12); + }); + + test('nested integration binds its own variable', () => { + const inner: ExpressionInput = [ + 'Integrate', + ['Max', 'x', 0], + ['Limits', 'x', -1, 1], + ]; + expect(integrate(['Multiply', 'x', inner]).json).toEqual([ + 'Rational', + 1, + 4, + ]); + }); +});