From f60ed31c3fd33834d79370d332e6793e691e22a6 Mon Sep 17 00:00:00 2001 From: "David P. Sanders" Date: Mon, 20 Apr 2026 09:50:00 -0500 Subject: [PATCH 1/3] Adopt ReversePropagation 0.8 SSAFunction + RuntimeGeneratedFunctions API Contractor now lowers its expression to an SSAFunction once via cse_equations and passes that SSA to forward_backward_contractor. The SSA is stored on the contractor so it can be reused across multiple constraints over the same expression. Separator shares the SSA between its forward evaluator and its HC4Revise contractor, so each constructor does the symbolic CSE pass only once. make_function switches from eval to @RuntimeGeneratedFunction, matching RP's own shift: generated functions can now be built and invoked within the same dynamic extent without hitting Julia world-age errors, with much faster build time. Bumps compat to ReversePropagation = "0.8" and adds RuntimeGeneratedFunctions as a direct dep. ICP package version bumped to 0.16.0. Drive-by fix to a garbled .gitignore entry. Co-Authored-By: Claude Opus 4.7 (1M context) --- .gitignore | 3 ++- CLAUDE.md | 25 +++++++++++++++++++++++++ Project.toml | 6 ++++-- src/IntervalConstraintProgramming.jl | 7 ++++++- src/contractor.jl | 24 +++++++++++++++++++++--- src/utils.jl | 3 ++- 6 files changed, 60 insertions(+), 8 deletions(-) diff --git a/.gitignore b/.gitignore index 44a110c..fa8b6af 100644 --- a/.gitignore +++ b/.gitignore @@ -2,4 +2,5 @@ *.jl.*.cov *.jl.mem docs/build/ -Manifest.toml \ No newline at end of file +Manifest.toml +papers/ diff --git a/CLAUDE.md b/CLAUDE.md index 395e1fc..51dc309 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -43,3 +43,28 @@ This file tracks decisions and context discovered while working on this package. - Chained comparisons like `0 <= x^2+y^2 <= 1` don't work (Julia lowers to `&&` which requires `Bool`). Users should use `x^2+y^2 ∈ interval(0, 1)` instead. A `@constraint` macro could fix this — see `future.md`. - ReversePropagation emits many "Method definition overwritten" warnings with Symbolics v7. Harmless but noisy — see `future.md`. + +## ReversePropagation SSAFunction / RuntimeGeneratedFunctions migration (2026-04-20) + +### Versions + +Bumped `ReversePropagation` compat from `"0.3 - 0.4"` to `"0.8"`. Added a direct dependency on `RuntimeGeneratedFunctions = "0.5"` (transitively available through RP, but worth owning directly since we now call `@RuntimeGeneratedFunction` ourselves). Bumped the package version to `0.16.0`. + +### Why + +RP 0.5+ reshapes `cse_equations`, `forward_backward_code`, `forward_backward_contractor`, and `gradient_code` around a new `SSAFunction` value: SSA `code`, an `output` variable, and a `variables` NamedTuple carrying extras like `constraint`. It also replaces `eval` with `@RuntimeGeneratedFunction` so a generated function can be built and called within the same dynamic extent (no world-age error) and with a 4–15× faster build time. + +### Code changes + +**`src/IntervalConstraintProgramming.jl`** — bring in `SSAFunction` and `cse_equations` from RP; call `RuntimeGeneratedFunctions.init(@__MODULE__)` so `@RuntimeGeneratedFunction` works inside this module; export `SSAFunction`. + +**`src/utils.jl`** — `make_function` now wraps `build_function(...)` in `@RuntimeGeneratedFunction` instead of `eval`. Same reasons as upstream RP: no world-age hazard, faster build. + +**`src/contractor.jl`** — `Contractor` now stores the `SSAFunction` it was built from: +- `Contractor(ex, vars)` calls `cse_equations(ex)` once and forwards to +- `Contractor(ssa::SSAFunction, ex, vars)` which calls `forward_backward_contractor(ssa, vars)` (the SSA-taking method). +- `Separator(ex, vars, constraint::Interval)` runs `cse_equations(ex)` once and shares the SSA between `make_function(ex, vars)` (forward evaluator) and `Contractor(ssa, ex, vars)` (HC4Revise), so each constructor does CSE on `ex` only once rather than twice. + +### Known issues + +- This change requires ReversePropagation `0.8.0`, which is not registered yet at the time of writing — to develop/test against it, `]dev` the local RP checkout. diff --git a/Project.toml b/Project.toml index 56e51d3..2107ee1 100644 --- a/Project.toml +++ b/Project.toml @@ -1,12 +1,13 @@ name = "IntervalConstraintProgramming" uuid = "138f1668-1576-5ad7-91b9-7425abbf3153" -version = "0.15.0" +version = "0.16.0" [deps] IntervalArithmetic = "d1acc4aa-44c8-5952-acd4-ba5d80a2a253" IntervalBoxes = "43d83c95-ebbb-40ec-8188-24586a1458ed" IntervalContractors = "15111844-de3b-5229-b4ba-526f2f385dc9" ReversePropagation = "527681c1-8309-4d3f-8790-caf822a419ba" +RuntimeGeneratedFunctions = "7e49a35a-f44a-4d26-94aa-eba1b4ca6b47" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7" @@ -14,7 +15,8 @@ Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7" IntervalArithmetic = "1" IntervalBoxes = "0.3" IntervalContractors = "0.6" -ReversePropagation = "0.4.1" +ReversePropagation = "0.8" +RuntimeGeneratedFunctions = "0.5" StaticArrays = "1" Symbolics = "6 - 7" julia = "1.10" diff --git a/src/IntervalConstraintProgramming.jl b/src/IntervalConstraintProgramming.jl index 12b02be..7672701 100644 --- a/src/IntervalConstraintProgramming.jl +++ b/src/IntervalConstraintProgramming.jl @@ -14,6 +14,10 @@ using Symbolics: @register_symbolic using StaticArrays using ReversePropagation +using ReversePropagation: SSAFunction, cse_equations + +using RuntimeGeneratedFunctions +RuntimeGeneratedFunctions.init(@__MODULE__) # using MacroTools @@ -50,7 +54,8 @@ export add_constraint!, variables, add_variables!, - ConstraintProblem + ConstraintProblem, + SSAFunction export ¬ diff --git a/src/contractor.jl b/src/contractor.jl index 3181bbd..a482a07 100644 --- a/src/contractor.jl +++ b/src/contractor.jl @@ -1,13 +1,28 @@ abstract type AbstractContractor end -struct Contractor{V, E, CC} <: AbstractContractor +""" + Contractor(ex, vars) + Contractor(ssa::SSAFunction, ex, vars) + +HC4Revise forward--backward contractor for the constraint `ex(vars...) ∈ constraint`. + +The symbolic expression is first lowered to SSA form via +`ReversePropagation.cse_equations`; the resulting `SSAFunction` is stored on the +contractor so it can be reused (e.g. shared between a `Separator` and other +contractors over the same expression). +""" +struct Contractor{V, E, S, CC} <: AbstractContractor vars::V ex::E + ssa::S contractor::CC end -Contractor(ex, vars) = Contractor(vars, ex, forward_backward_contractor(ex, vars)) +Contractor(ex, vars) = Contractor(cse_equations(ex), ex, vars) + +Contractor(ssa::SSAFunction, ex, vars) = + Contractor(vars, ex, ssa, forward_backward_contractor(ssa, vars)) (CC::Contractor)(X, constraint=interval(0.0)) = IntervalBox(CC.contractor(X, constraint)[1]) @@ -39,7 +54,10 @@ function Separator(orig_expr, vars) return Separator(ex, vars, constraint) end -Separator(ex, vars, constraint::Interval) = Separator(vars, ex, constraint, make_function(ex, vars), Contractor(ex, vars)) +function Separator(ex, vars, constraint::Interval) + ssa = cse_equations(ex) + return Separator(vars, ex, constraint, make_function(ex, vars), Contractor(ssa, ex, vars)) +end function separate_infinite_box(S::Separator, X::IntervalBox) # for an box that extends to infinity we cannot evaluate at a corner diff --git a/src/utils.jl b/src/utils.jl index 64b76cd..bd17adc 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -1,6 +1,7 @@ -make_function(ex, vars) = eval(build_function(ex, vars, nanmath=false)) +make_function(ex, vars) = + @RuntimeGeneratedFunction(build_function(ex, vars, nanmath=false)) """ From cd1a1bd6e379a8b14eefba44e4c6bf9b76ec5f3b Mon Sep 17 00:00:00 2001 From: "David P. Sanders" Date: Mon, 20 Apr 2026 10:11:55 -0500 Subject: [PATCH 2/3] Add exact=true kwarg to suppress NG flag on constraint output MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit IntervalArithmetic v1 taints any interval that touches a plain (non- ExactReal) Real with the "not guaranteed" (NG) flag. ICP's generated code is full of those — the literals from the user's expression (the 1 in x^2+y^2==1) and the integer exponents CSE inserts (the 2 in x^2) — so C.f([3..4, 5..6]) came back as [33.0, 51.0]_com_NG. Wrapping the literals inside the symbolic expression doesn't survive: Symbolics normalises ExactReal back to a plain number during its own simplification. The fix has to happen on the generated Expr. Add exactify(ex), which walks an Expr and wraps every non-Bool Real literal in IntervalArithmetic.exact. Thread an exact::Bool = false keyword through Contractor, Separator, and separator/constraint: - make_function applies exactify to the build_function output before @RuntimeGeneratedFunction when exact=true. - build_fb_contractor replicates RP's outer wrapper locally so it can exactify the body before RGF. When exact=false it forwards to forward_backward_contractor unchanged (no overhead in the default). Default behaviour is unchanged; opting in gives julia> constraint(x^2+y^2==1, vars; exact=true).f([3..4, 5..6]) [33.0, 51.0]_com Co-Authored-By: Claude Opus 4.7 (1M context) --- CLAUDE.md | 30 ++++++++++++++++++++++++++ src/contractor.jl | 55 ++++++++++++++++++++++++++++++++++++++--------- src/utils.jl | 50 +++++++++++++++++++++++++++++------------- 3 files changed, 110 insertions(+), 25 deletions(-) diff --git a/CLAUDE.md b/CLAUDE.md index 51dc309..ffa7048 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -68,3 +68,33 @@ RP 0.5+ reshapes `cse_equations`, `forward_backward_code`, `forward_backward_con ### Known issues - This change requires ReversePropagation `0.8.0`, which is not registered yet at the time of writing — to develop/test against it, `]dev` the local RP checkout. + +## Suppressing the "NG" flag via `exact=true` (2026-04-20) + +### Problem + +`IntervalArithmetic` v1 sets a "not guaranteed" (NG) flag on any interval that has ever touched a plain (non-`ExactReal`) `Real`. ICP's generated code is full of those — both the literals from the user's expression (the `1` in `x^2 + y^2 == 1`) and the integer exponents that CSE inserts (the `2` in `x^2`). So even `C.f([3..4, 5..6])` came back as `[33.0, 51.0]_com_NG`. + +### Why not on the symbolic side + +Wrapping literals in the user-facing symbolic expression (`x^2 + y^2 - exact(1)`) does not survive — `Symbolics` normalises `ExactReal` back to a plain `Int`/`Float` during its own simplification. Any fix has to happen on the generated `Expr`, after `Symbolics.build_function` / `ReversePropagation.forward_backward_expr`. + +### Design + +Added an `exactify(ex)` helper in `src/utils.jl` that walks an `Expr` and wraps every non-`Bool` `Real` literal with `IntervalArithmetic.exact(...)`. Threaded an `exact::Bool = false` keyword argument through `Contractor`, `Separator`, and `separator`/`constraint`: + +- `make_function(ex, vars; exact=false)` — runs `exactify` on the `build_function` output when `exact=true`, before `@RuntimeGeneratedFunction`. +- `Contractor(ssa, ex, vars; exact=false)` — when `exact=true`, delegates to `build_fb_contractor`, which replicates RP's outer wrapper locally so it can `exactify` the body before RGF. When `exact=false`, still forwards to `forward_backward_contractor` unchanged (no overhead in the default path). +- `Separator(...; exact=false)` and `separator(...; exact=false)` — threads the keyword through. + +### Usage + +```julia +julia> C = constraint(x^2 + y^2 == 1, vars; exact=true); +julia> C.f([3..4, 5..6]) +[33.0, 51.0]_com # no NG +julia> C.contractor((-10..10, -10..10)) +[-1.0, 1.0]_trv² # no NG +``` + +Default behaviour (`exact=false`) is unchanged. diff --git a/src/contractor.jl b/src/contractor.jl index a482a07..66161b5 100644 --- a/src/contractor.jl +++ b/src/contractor.jl @@ -2,8 +2,8 @@ abstract type AbstractContractor end """ - Contractor(ex, vars) - Contractor(ssa::SSAFunction, ex, vars) + Contractor(ex, vars; exact=false) + Contractor(ssa::SSAFunction, ex, vars; exact=false) HC4Revise forward--backward contractor for the constraint `ex(vars...) ∈ constraint`. @@ -11,6 +11,10 @@ The symbolic expression is first lowered to SSA form via `ReversePropagation.cse_equations`; the resulting `SSAFunction` is stored on the contractor so it can be reused (e.g. shared between a `Separator` and other contractors over the same expression). + +If `exact=true`, every numeric literal in the generated forward--backward code +is wrapped in `IntervalArithmetic.exact` (see [`exactify`](@ref)), suppressing +the "NG" flag on results. """ struct Contractor{V, E, S, CC} <: AbstractContractor vars::V @@ -19,10 +23,39 @@ struct Contractor{V, E, S, CC} <: AbstractContractor contractor::CC end -Contractor(ex, vars) = Contractor(cse_equations(ex), ex, vars) - -Contractor(ssa::SSAFunction, ex, vars) = - Contractor(vars, ex, ssa, forward_backward_contractor(ssa, vars)) +Contractor(ex, vars; exact::Bool = false) = + Contractor(cse_equations(ex), ex, vars; exact = exact) + +Contractor(ssa::SSAFunction, ex, vars; exact::Bool = false) = + Contractor(vars, ex, ssa, build_fb_contractor(ssa, vars; exact = exact)) + +# Build the forward--backward contractor function. When `exact=false`, just +# forward to RP's `forward_backward_contractor`. When `exact=true`, re-assemble +# the outer wrapper so we can walk the generated code and mark every Real +# literal as `IntervalArithmetic.exact` before handing it to RGF — otherwise +# the literals from the user's constraint (and the integer exponents inserted +# by CSE) drag the NG flag into the result. +function build_fb_contractor(ssa::SSAFunction, vars; exact::Bool = false) + exact || return forward_backward_contractor(ssa, vars) + + code, final_var, constraint_var = + ReversePropagation.forward_backward_expr(ssa, vars) + + input_vars = toexpr(Symbolics.MakeTuple(vars)) + final = toexpr(final_var) + constraint = toexpr(constraint_var) + + function_code = :( + (__inputs, __constraint) -> begin + $input_vars = __inputs + $constraint = __constraint + $(exactify(code)) + return $input_vars, $(final) + end + ) + + return @RuntimeGeneratedFunction(function_code) +end (CC::Contractor)(X, constraint=interval(0.0)) = IntervalBox(CC.contractor(X, constraint)[1]) @@ -48,15 +81,17 @@ function Base.show(io::IO, S::AbstractSeparator) end end -function Separator(orig_expr, vars) +function Separator(orig_expr, vars; exact::Bool = false) ex, constraint = analyse(orig_expr) - return Separator(ex, vars, constraint) + return Separator(ex, vars, constraint; exact = exact) end -function Separator(ex, vars, constraint::Interval) +function Separator(ex, vars, constraint::Interval; exact::Bool = false) ssa = cse_equations(ex) - return Separator(vars, ex, constraint, make_function(ex, vars), Contractor(ssa, ex, vars)) + return Separator(vars, ex, constraint, + make_function(ex, vars; exact = exact), + Contractor(ssa, ex, vars; exact = exact)) end function separate_infinite_box(S::Separator, X::IntervalBox) diff --git a/src/utils.jl b/src/utils.jl index bd17adc..ae79775 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -1,7 +1,27 @@ -make_function(ex, vars) = - @RuntimeGeneratedFunction(build_function(ex, vars, nanmath=false)) +""" + exactify(ex) + +Walk `ex` and replace every non-`Bool` `Real` literal with `IntervalArithmetic.exact(literal)`. + +Used on the `Expr` produced by symbolic code generation (`Symbolics.build_function` +or `ReversePropagation.forward_backward_expr`) so that, when the generated function +is called on `Interval` arguments, numeric literals from the user's constraint +(like the `1` in `x^2 + y^2 == 1`) are treated as exact and do not introduce the +"NG" (not-guaranteed) flag in the result. +""" +exactify(ex::Bool) = ex +exactify(ex::Real) = :($(IntervalArithmetic.exact)($ex)) +exactify(ex::Expr) = Expr(ex.head, exactify.(ex.args)...) +exactify(ex) = ex + + +make_function(ex, vars; exact::Bool = false) = begin + code = build_function(ex, vars, nanmath=false) + exact && (code = exactify(code)) + @RuntimeGeneratedFunction(code) +end """ @@ -42,53 +62,53 @@ end -function separator(ex, vars) - ex2 = ex +function separator(ex, vars; exact::Bool = false) + ex2 = ex if ex isa Num ex2 = value(ex) end - + op = operation(ex2) if op == ¬ arg = arguments(ex2)[1] - return ¬(separator(arg, vars)) + return ¬(separator(arg, vars; exact = exact)) end lhs, rhs = arguments(ex2) if op == & - return separator(lhs, vars) ⊓ separator(rhs, vars) + return separator(lhs, vars; exact = exact) ⊓ separator(rhs, vars; exact = exact) elseif op == | - return separator(lhs, vars) ⊔ separator(rhs, vars) + return separator(lhs, vars; exact = exact) ⊔ separator(rhs, vars; exact = exact) elseif op ∈ (≤, <) constraint = interval(-Inf, 0) - Separator(Num(lhs - rhs), vars, constraint) + Separator(Num(lhs - rhs), vars, constraint; exact = exact) elseif op ∈ (≥, >) constraint = interval(0, +Inf) - Separator(Num(lhs - rhs), vars, constraint) + Separator(Num(lhs - rhs), vars, constraint; exact = exact) elseif op == (==) constraint = interval(0, 0) - Separator(Num(lhs - rhs), vars, constraint) + Separator(Num(lhs - rhs), vars, constraint; exact = exact) else - return Separator(ex, vars, interval(0, 0)) # implicit "== 0" + return Separator(ex, vars, interval(0, 0); exact = exact) # implicit "== 0" end - + end -function separator(ex) +function separator(ex; exact::Bool = false) vars = Symbolics.get_variables(ex) - return separator(ex, vars) + return separator(ex, vars; exact = exact) end From 1beb7f25b10317bb321437a84bac10e26c844160 Mon Sep 17 00:00:00 2001 From: "David P. Sanders" Date: Mon, 20 Apr 2026 10:27:57 -0500 Subject: [PATCH 3/3] Preserve decorations through HC4Revise contraction under exact=true MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit IEEE 1788 weakens intersect_interval(a, b) to the trv decoration because, in general, the result may not be a subset of the range of any function. Every reverse op in IntervalContractors is built on that same intersect, so both the outer intersect_interval(_value, _constraint) that RP emits and every narrowing step in the reverse sweep drag trv through the contractor output. In HC4Revise, however, each intersection/reverse-op is narrowing *within* an already-evaluated function range, so min(decoration(inputs)...) is the valid decoration — not trv. Add preserve_decorations(ex) in src/utils.jl: a walker that rewrites the body returned by ReversePropagation.forward_backward_expr so that every IntervalArithmetic.intersect_interval call and every ReversePropagation._rev_* reverse-op call has its output re-decorated with min(decoration(inputs)...). Because Symbolics' toexpr emits the call head as the actual function value (not Expr(:., Module, :name)), the matchers compare ex.args[1] === IntervalArithmetic.intersect_interval and parentmodule(f) === ReversePropagation && startswith(nameof(f), "_rev_") directly. _redecorate(x, d) unconditionally replaces the decoration (not min with the op's own output decoration, which is precisely the weakened trv we're overriding). Empty intervals stay trv (empty semantics). build_fb_contractor now runs preserve_decorations alongside exactify under exact=true. With both on: julia> constraint(x^2+y^2==1, vars; exact=true).contractor((-10..10, -10..10)) [-1.0, 1.0]_com² # was [-1.0, 1.0]_trv_NG² (The outer of the Separator triple can still end up _trv because Separator computes it via boundary ⊔ corner and hull intrinsically weakens under IEEE 1788 — intended semantics.) Co-Authored-By: Claude Opus 4.7 (1M context) --- CLAUDE.md | 16 ++++++-- src/contractor.jl | 20 ++++++---- src/utils.jl | 94 +++++++++++++++++++++++++++++++++++++++++++++++ 3 files changed, 119 insertions(+), 11 deletions(-) diff --git a/CLAUDE.md b/CLAUDE.md index ffa7048..062f074 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -84,17 +84,25 @@ Wrapping literals in the user-facing symbolic expression (`x^2 + y^2 - exact(1)` Added an `exactify(ex)` helper in `src/utils.jl` that walks an `Expr` and wraps every non-`Bool` `Real` literal with `IntervalArithmetic.exact(...)`. Threaded an `exact::Bool = false` keyword argument through `Contractor`, `Separator`, and `separator`/`constraint`: - `make_function(ex, vars; exact=false)` — runs `exactify` on the `build_function` output when `exact=true`, before `@RuntimeGeneratedFunction`. -- `Contractor(ssa, ex, vars; exact=false)` — when `exact=true`, delegates to `build_fb_contractor`, which replicates RP's outer wrapper locally so it can `exactify` the body before RGF. When `exact=false`, still forwards to `forward_backward_contractor` unchanged (no overhead in the default path). +- `Contractor(ssa, ex, vars; exact=false)` — when `exact=true`, delegates to `build_fb_contractor`, which replicates RP's outer wrapper locally so it can `exactify` and run `preserve_decorations` on the body before RGF. When `exact=false`, still forwards to `forward_backward_contractor` unchanged (no overhead in the default path). - `Separator(...; exact=false)` and `separator(...; exact=false)` — threads the keyword through. +### Decoration preservation through contraction + +IEEE 1788 weakens `intersect_interval(a, b)` to the `trv` decoration because, in general, the result may not be a subset of the range of any function. Every reverse op in `IntervalContractors` (`plus_rev`, `power_rev`, …) is built on top of that same intersect, so both the outer `intersect_interval(_value, _constraint)` that RP emits *and* every narrowing step in the reverse sweep drag `trv` through the result. HC4Revise however is only ever narrowing *within* an already-evaluated function range, so the right answer is `min(decoration(inputs)...)` — not `trv`. + +`preserve_decorations(ex)` walks the body `Expr` returned by `ReversePropagation.forward_backward_expr` and wraps each `IntervalArithmetic.intersect_interval` call and each `ReversePropagation._rev_*` reverse-op call so its output interval(s) are re-decorated with `min(decoration(inputs)...)`. Empty intervals stay `trv` (empty semantics). Because Symbolics' `toexpr` emits the call head as the actual function *value* (not `Expr(:., Module, :name)`), the matchers compare `ex.args[1] === IntervalArithmetic.intersect_interval` and `parentmodule(f) === ReversePropagation && startswith(nameof(f), "_rev_")` directly. + +Re-decoration strategy: `_redecorate(x, d) = interval(inf(x), sup(x), d)` unconditionally replaces the decoration; it does *not* `min` with `decoration(x)`, because the op's own output decoration is exactly the weakened `trv` we're trying to override. Narrowing is a subset, so the pre-op input decoration is the valid upper bound regardless of what the op stamped on. Empty intervals keep their existing (`trv`) decoration. + ### Usage ```julia julia> C = constraint(x^2 + y^2 == 1, vars; exact=true); julia> C.f([3..4, 5..6]) -[33.0, 51.0]_com # no NG +[33.0, 51.0]_com # no NG julia> C.contractor((-10..10, -10..10)) -[-1.0, 1.0]_trv² # no NG +[-1.0, 1.0]_com² # no NG, no trv ``` -Default behaviour (`exact=false`) is unchanged. +Default behaviour (`exact=false`) is unchanged. Note: the outer of the full `Separator` triple can still end up `_trv` because it is computed via `boundary ⊔ corner` in `(SS::Separator)(X)` — hull (`⊔`) intrinsically weakens in IEEE 1788, and that is the intended semantics. diff --git a/src/contractor.jl b/src/contractor.jl index 66161b5..84bc724 100644 --- a/src/contractor.jl +++ b/src/contractor.jl @@ -13,8 +13,11 @@ contractor so it can be reused (e.g. shared between a `Separator` and other contractors over the same expression). If `exact=true`, every numeric literal in the generated forward--backward code -is wrapped in `IntervalArithmetic.exact` (see [`exactify`](@ref)), suppressing -the "NG" flag on results. +is wrapped in `IntervalArithmetic.exact` (see [`exactify`](@ref)) and every +`intersect_interval` / `_rev_*` call site is wrapped so its output decoration +is `min(decoration(inputs)...)` instead of the IEEE 1788 default `trv` (see +[`preserve_decorations`](@ref)). Together these suppress the "NG" flag and +preserve `com` / `dac` decorations through contraction. """ struct Contractor{V, E, S, CC} <: AbstractContractor vars::V @@ -31,10 +34,11 @@ Contractor(ssa::SSAFunction, ex, vars; exact::Bool = false) = # Build the forward--backward contractor function. When `exact=false`, just # forward to RP's `forward_backward_contractor`. When `exact=true`, re-assemble -# the outer wrapper so we can walk the generated code and mark every Real -# literal as `IntervalArithmetic.exact` before handing it to RGF — otherwise -# the literals from the user's constraint (and the integer exponents inserted -# by CSE) drag the NG flag into the result. +# the outer wrapper so we can (a) mark every Real literal as +# `IntervalArithmetic.exact` (so constraint literals and CSE-inserted exponents +# don't drag the NG flag through the result) and (b) clamp decorations through +# intersect/reverse-op call sites to `min(input decorations)` instead of the +# IEEE 1788 default `trv` weakening. function build_fb_contractor(ssa::SSAFunction, vars; exact::Bool = false) exact || return forward_backward_contractor(ssa, vars) @@ -45,11 +49,13 @@ function build_fb_contractor(ssa::SSAFunction, vars; exact::Bool = false) final = toexpr(final_var) constraint = toexpr(constraint_var) + body = preserve_decorations(exactify(code)) + function_code = :( (__inputs, __constraint) -> begin $input_vars = __inputs $constraint = __constraint - $(exactify(code)) + $body return $input_vars, $(final) end ) diff --git a/src/utils.jl b/src/utils.jl index ae79775..1b0785a 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -17,6 +17,100 @@ exactify(ex::Expr) = Expr(ex.head, exactify.(ex.args)...) exactify(ex) = ex +# --- decoration preservation for HC4Revise --------------------------------- +# +# IEEE 1788 weakens `intersect_interval(a, b)` to the `trv` decoration because, +# in general, the intersection may not be the range of any function. In +# HC4Revise however each intersection/reverse op is narrowing *within* an +# already-evaluated function range, so `min(decoration(inputs)...)` is the +# right answer. The helpers below re-decorate the output accordingly; the +# Expr walker `preserve_decorations!` rewrites the call sites that RP's +# `forward_backward_expr` emits so the weakening never takes hold. + +_dec_of(x::Interval) = decoration(x) +_dec_of(x) = IntervalArithmetic.com + +_min_input_decoration(inputs::Tuple) = minimum(_dec_of, inputs) + +# Re-decorate an interval with `d`. Used to undo the `trv` weakening that +# `intersect_interval` and reverse ops perform: in HC4 the output is always a +# subset of one of the inputs, so `d = min(decoration(inputs)...)` is the valid +# decoration, regardless of what the op itself stamped on. Empty intervals +# stay as-is (`trv` / empty semantics). +_redecorate(x::Interval, d::IntervalArithmetic.Decoration) = + isempty_interval(x) ? x : interval(inf(x), sup(x), d) +_redecorate(x, ::IntervalArithmetic.Decoration) = x + +"Wrap an `intersect_interval`-style call: single interval result." +function preserve_dec_intersect(result::Interval, inputs::Tuple) + return _redecorate(result, _min_input_decoration(inputs)) +end + +"Wrap a reverse-op call: tuple result, each interval element re-decorated." +function preserve_dec_reverse(result::Tuple, inputs::Tuple) + d = _min_input_decoration(inputs) + return map(x -> _redecorate(x, d), result) +end + + +# ---- Expr walker ---------------------------------------------------------- +# +# Matches the two forms RP's generated body emits: +# _lhs = (IntervalArithmetic.intersect_interval)(a, b) +# (t1, t2, …) = (ReversePropagation.var"_rev_OP")(arg1, arg2, …) +# and rewrites the rhs so that output decorations are clamped to +# `min(decoration(arg_i))`. + +# Symbolics' `toexpr` resolves call heads to actual function *values* (not +# `Expr(:., Module, QuoteNode(name))`), so these checks compare against the +# function objects directly. + +_is_intersect_call(ex) = false +function _is_intersect_call(ex::Expr) + ex.head === :call || return false + return ex.args[1] === IntervalArithmetic.intersect_interval +end + +_is_reverse_op_call(ex) = false +function _is_reverse_op_call(ex::Expr) + ex.head === :call || return false + f = ex.args[1] + f isa Function || return false + name = string(nameof(f)) + return startswith(name, "_rev_") && + parentmodule(f) === ReversePropagation +end + +""" + preserve_decorations(ex) + +Walk `ex` and wrap calls to `IntervalArithmetic.intersect_interval` and +`ReversePropagation._rev_*` so that each output interval's decoration is +clamped to `min(decoration(inputs)...)` instead of being weakened to `trv`. + +Called on the body `Expr` returned by +`ReversePropagation.forward_backward_expr` before `@RuntimeGeneratedFunction`. +""" +preserve_decorations(ex) = ex +function preserve_decorations(ex::Expr) + if ex.head === :(=) && ex.args[2] isa Expr && _is_intersect_call(ex.args[2]) + lhs = ex.args[1] + call = ex.args[2] + inputs = Expr(:tuple, call.args[2:end]...) + new_rhs = :($(preserve_dec_intersect)($call, $inputs)) + return Expr(:(=), lhs, new_rhs) + elseif ex.head === :(=) && ex.args[2] isa Expr && _is_reverse_op_call(ex.args[2]) + lhs = ex.args[1] + call = ex.args[2] + inputs = Expr(:tuple, call.args[2:end]...) + new_rhs = :($(preserve_dec_reverse)($call, $inputs)) + return Expr(:(=), lhs, new_rhs) + else + return Expr(ex.head, map(preserve_decorations, ex.args)...) + end +end + + make_function(ex, vars; exact::Bool = false) = begin code = build_function(ex, vars, nanmath=false) exact && (code = exactify(code))