diff --git a/Project.toml b/Project.toml index cc59232a..ac7bcb2c 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "AbstractQAtlas" uuid = "dcea2817-62f8-4a75-b498-1b50a9ed1e4d" -version = "0.7.16" +version = "0.7.17" authors = ["sota shimozono "] [deps] diff --git a/examples/critical_entanglement.jl b/examples/critical_entanglement.jl new file mode 100644 index 00000000..52037f1f --- /dev/null +++ b/examples/critical_entanglement.jl @@ -0,0 +1,87 @@ +# Reading a central charge off measured entanglement, and what a convention +# declaration buys. +# +# The numbers are not pasted in from anywhere: a critical free-fermion chain's +# block correlation matrix is a closed form, so the entropies below are computed +# here and the central charge that comes out is a measurement, not a fixture. Its +# distance from the exact `c = 1` is the finite-block correction, which is why +# `test/core/test_examples.jl` checks it to 1e-3 and not to machine precision. +# +# Wrapped in a module because `test/core/test_examples.jl` includes this file, and +# every test file is included into one namespace: `block_entropy` and a couple of +# block sizes at top level would be in every other test file's way. +# +# Run: julia --project=. examples/critical_entanglement.jl + +module CriticalEntanglementExample + +using AbstractQAtlas +using LinearAlgebra: eigvals, Symmetric + +export block_entropy, central_charges + +""" + correlation_matrix(ℓ) -> Symmetric + +The `ℓ x ℓ` block of the half-filled critical free-fermion chain's two-point +function, `C_jk = sin(π(j-k)/2) / (π(j-k))`, `C_jj = 1/2`. +""" +function correlation_matrix(ℓ::Integer) + return Symmetric([ + j == k ? 0.5 : sin(π * (j - k) / 2) / (π * (j - k)) for j in 1:ℓ, k in 1:ℓ + ]) +end + +"Entanglement entropy of a block of `ℓ` sites, in nats." +function block_entropy(ℓ::Integer) + return free_fermion_entanglement_entropy(eigvals(correlation_matrix(ℓ))) +end + +""" + central_charges(small = 16, large = 64) -> NamedTuple + +`c` read off the same measurement entered three ways: in nats, in bits with +nothing declared, and in bits declared. + +Two block sizes, because the charge is not falsifiable from one: the logarithmic +forms each carry a non-universal constant that can absorb any `c`, and the slope +through two blocks carries none. The middle entry is the failure the convention +layer exists for, and it is wrong by exactly the base nobody mentioned. +""" +function central_charges(small::Integer=16, large::Integer=64) + a, b = Region((1:small)...), Region((1:large)...) + span = log(large) - log(small) + slope(bg) = (bg[entanglement_entropy(b)] - bg[entanglement_entropy(a)]) / span + pack(f) = (entanglement_entropy(a) => f(small), entanglement_entropy(b) => f(large)) + # Read out of the bag, so a declared convention reaches the slope through it. + c(bg) = derive_crosschecked(:c; dS_dlogℓ=slope(bg), ncuts=2) + in_bits(ℓ) = block_entropy(ℓ) / log(2) + return ( + nats=c(bag(pack(block_entropy)...)), + bits_undeclared=c(bag(pack(in_bits)...)), + bits_declared=c( + bag(conventions(AbstractEntanglementMeasure => Bits), pack(in_bits)...) + ), + ) +end + +end # module CriticalEntanglementExample + +if abspath(PROGRAM_FILE) == @__FILE__ + let cs = CriticalEntanglementExample.central_charges() + println("c from blocks of 16 and 64 sites, exact value 1:") + println(" nats, as the calculation produced them : ", round(cs.nats; digits=4)) + println( + " the same numbers in bits, undeclared : ", + round(cs.bits_undeclared; digits=4), + ) + println( + " which is 1/ln2 = ", + round(1 / log(2); digits=4), + ", the base it was never told", + ) + println( + " the same numbers in bits, declared : ", round(cs.bits_declared; digits=4) + ) + end +end diff --git a/ext/AbstractQAtlasForwardDiffExt.jl b/ext/AbstractQAtlasForwardDiffExt.jl index b6fb7c79..89d3576b 100644 --- a/ext/AbstractQAtlasForwardDiffExt.jl +++ b/ext/AbstractQAtlasForwardDiffExt.jl @@ -10,15 +10,24 @@ module AbstractQAtlasForwardDiffExt using AbstractQAtlas -import AbstractQAtlas: peierls_current, thermal_derivative # extended below → must import +import AbstractQAtlas: nth_derivative, peierls_current, thermal_derivative using AbstractQAtlas: - _genealogy_derivative, response_order, indices, Susceptibility, SpecificHeat, Energy + _genealogy_derivative, + _susceptibility_derivative, + indices, + Susceptibility, + SpecificHeat, + Energy using ForwardDiff: derivative # n-th derivative of a scalar function by nested ForwardDiff (n small — # response orders are 1–3). _nth(f, x, n::Integer) = n == 0 ? f(x) : _nth(y -> derivative(f, y), x, n - 1) +# `AutoDiff` is a route like the others: this is the method whose absence made +# every route-taking entry point special-case it by name. +nth_derivative(::AbstractQAtlas.AutoDiff, f, x, n::Integer) = _nth(f, x, n) + # ── generic, genealogy-driven response ──────────────────────────────────── # Which order, which field and the net sign all come from `_genealogy_derivative` # (src/derivative_routes.jl), shared with the finite-difference routes so the two @@ -28,19 +37,11 @@ function thermal_derivative(q::AbstractQuantity, F, x::Number) return _genealogy_derivative(q, F, x, _nth) end -# χ⁽ⁿ⁾_{α;β₁…βₙ} = −∂ⁿ⁺¹F/∂h_α∂h_{β₁}…∂h_{βₙ}. With a SINGLE-field function -# F(h) only the DIAGONAL component (all indices equal) is defined — an -# off-diagonal component needs partials w.r.t. distinct field directions, -# so guard against silently returning the diagonal for an off-diagonal ask. +# χ⁽ⁿ⁾_{α;β₁…βₙ} = −∂ⁿ⁺¹F/∂h_α∂h_{β₁}…∂h_{βₙ}, diagonal only from a single-field +# F(h). The guard and the order live in `_susceptibility_derivative`, shared with +# the finite-difference routes; this supplies the AD `nth`. function thermal_derivative(χ::Susceptibility, F, h::Number) - idx = indices(χ) - all(==(idx[1]), idx) || error( - "thermal_derivative(::Susceptibility, F, h::Number) with a single-field " * - "function computes only the DIAGONAL χ⁽ⁿ⁾ (all indices equal); got " * - "off-diagonal $(idx). Pass a multi-field potential F(h⃗) and the field-" * - "component ordering: thermal_derivative(χ, F, h⃗, components).", - ) - return -_nth(F, h, response_order(χ) + 1) + return _susceptibility_derivative(χ, F, h, _nth) end # Multi-field / off-diagonal: F is a function of a field VECTOR `h⃗`, and diff --git a/src/AbstractQAtlas.jl b/src/AbstractQAtlas.jl index 45c6cb90..919cfc8c 100644 --- a/src/AbstractQAtlas.jl +++ b/src/AbstractQAtlas.jl @@ -136,6 +136,7 @@ end module QuantumInformation using ..AbstractQAtlas import ..AbstractQAtlas: fetch # unexported, so bare `fetch` here would be Base's + import ..AbstractQAtlas: _solve # extended for a closed-form inverse (EntanglementSpectrumCorrelation:ζ) using ExperimentalAPI: @experimental include("relations/entanglement.jl") end diff --git a/src/core/conventions.jl b/src/core/conventions.jl index 8944874a..f5869158 100644 --- a/src/core/conventions.jl +++ b/src/core/conventions.jl @@ -203,16 +203,31 @@ export conventions """ declared_convention(cs::ConventionSet, Q::Type) -> Union{Convention,Nothing} -What `cs` says `Q`'s values are written in, walking up to `Q`'s supertypes and -taking the most specific entry; `nothing` when nothing in `cs` covers `Q`. +What `cs` says `Q`'s values are written in, taking the most specific entry that +`Q` is a subtype of; `nothing` when nothing in `cs` covers `Q`. + +Matched by `<:`, not by walking `supertype`, because a parametric quantity's +supertype chain SKIPS its own family: `supertype(Energy{:per_site})` is +`AbstractThermalPotential`, so a walk never reaches the `Energy` a project keyed +its declaration on, and the value goes into the bag unconverted. + +Covers with no unique most specific member are refused rather than resolved by +`Dict` order. Decided after collecting every cover, not folded pairwise: two +unrelated covers can be reconciled by a third that refines both, and a fold that +errors on meeting the first incomparable pair reports a false ambiguity for four +of the six orders that Dict iteration can hand it. """ function declared_convention(cs::ConventionSet, @nospecialize(Q::Type)) - T = Q - while T !== Any - haskey(cs.declared, T) && return cs.declared[T] - T = supertype(T) + covers = [(T, c) for (T, c) in cs.declared if Q <: T] + isempty(covers) && return nothing + for (T, c) in covers + all(U -> T <: U, first.(covers)) && return c end - return nothing + return error( + "declared_convention: $Q is covered by $(join(first.(covers), ", ")), with no " * + "one of them a subtype of all the others, so none is the most specific. Key " * + "the declaration on whichever one the values were actually written in.", + ) end export declared_convention diff --git a/src/derivative_routes.jl b/src/derivative_routes.jl index 39888e29..0a7785e2 100644 --- a/src/derivative_routes.jl +++ b/src/derivative_routes.jl @@ -84,6 +84,10 @@ struct Richardson <: DerivativeRoute end export Richardson +# Whether the named package is loaded at all, which separates "install it" from +# "it is here and the route's own method is missing". +_package_loaded(name::Symbol) = any(k -> k.name == String(name), keys(Base.loaded_modules)) + """ nth_derivative(route::DerivativeRoute, f, x, n::Integer) -> value @@ -93,19 +97,71 @@ The `n`-th derivative of the scalar function `f` at `x`, taken along `route`. This is the one method a new route has to define. """ function nth_derivative(route::DerivativeRoute, f, x, n::Integer) - return error("nth_derivative: no method for $(typeof(route)).") + pkg = backend_package(route) + pkg === nothing && return error("nth_derivative: no method for $(typeof(route)).") + # Reaching the fallback with the backend LOADED means the route's own method is + # missing or its signature does not match, which reinstalling does not fix. + # Saying "not loaded" there points away from the defect. + _package_loaded(pkg) && return error( + "nth_derivative: $(typeof(route)) declares the $pkg backend and $pkg is " * + "loaded, but no method matched. The route's own `nth_derivative` is missing " * + "or its signature differs.", + ) + return throw(MissingRouteBackend(route)) end export nth_derivative -function nth_derivative(::AutoDiff, f, x, n::Integer) - return error( - "nth_derivative(AutoDiff(), ...) needs an automatic-differentiation " * - "backend: run `using ForwardDiff` to load the AbstractQAtlas AD extension, " * - "or take a finite-difference route (CentralDifference / Richardson), which " * - "needs none.", +""" + MissingRouteBackend(route) <: Exception + +Thrown when a [`DerivativeRoute`](@ref) needs a package extension that is not +loaded. + +Its own type, because [`derivative_report`](@ref) has to tell it from every other +way a route can fail. Catching `Exception` there would turn a diagnosed refusal, +an off-diagonal susceptibility or a potential evaluated outside its domain, into +the same `NaN` row as an unloaded backend. +""" +struct MissingRouteBackend <: Exception + route::DerivativeRoute + function MissingRouteBackend(route::DerivativeRoute) + backend_package(route) === nothing && throw( + ArgumentError( + "MissingRouteBackend: $(typeof(route)) declares no `backend_package`, " * + "so there is no extension for it to be missing. Its failure is not a " * + "missing backend.", + ), + ) + return new(route) + end +end +export MissingRouteBackend + +function Base.showerror(io::IO, e::MissingRouteBackend) + return print( + io, + "MissingRouteBackend: $(typeof(e.route)) needs the $(backend_package(e.route)) ", + "extension, which is not loaded. Run `using $(backend_package(e.route))`, or ", + "take a finite-difference route (CentralDifference / Richardson), which needs ", + "no backend.", ) end +""" + backend_package(route::DerivativeRoute) -> Union{Symbol,Nothing} + +The package whose extension supplies `route`'s [`nth_derivative`](@ref), or +`nothing` for a route that needs none. + +A route declaring one and finding no method gets +[`MissingRouteBackend`](@ref) rather than a bare "no method", which is the +difference between "install this" and "this route does not exist". Declared here +and not in the extension: the point is to answer when the extension is ABSENT. +""" +backend_package(::DerivativeRoute) = nothing +backend_package(::AutoDiff) = :ForwardDiff +export backend_package + _central(f, x, h) = (f(x + h) - f(x - h)) / (2h) function nth_derivative(route::CentralDifference, f, x, n::Integer) @@ -129,6 +185,36 @@ function nth_derivative(route::Richardson, f, x, n::Integer) return only(t) end +""" + step_size(route::DerivativeRoute) -> Union{Real,Nothing} + with_step_size(route::DerivativeRoute, h::Real) -> DerivativeRoute + +The step `route` takes, and the same route at a different step. `nothing` means +the route has no step, which is what [`AutoDiff`](@ref) reports. + +Part of the route contract alongside [`nth_derivative`](@ref), and the pair +[`observed_order`](@ref) needs. Each missing half names itself: a route reporting +a `step_size` with no `with_step_size` is told so by `with_step_size`, and one +declaring neither is told it reports no step, which for it is true. A closed +`Union` over the routes that happened to exist told a third route the second +thing whether or not it was true. +""" +step_size(::DerivativeRoute) = nothing +export step_size + +function with_step_size(route::DerivativeRoute, h::Real) + return error( + "with_step_size: $(typeof(route)) defines no `with_step_size`. A route that " * + "reports a `step_size` needs one, so `observed_order` can halve it.", + ) +end +export with_step_size + +step_size(r::CentralDifference) = r.h +step_size(r::Richardson) = r.h +with_step_size(::CentralDifference, h::Real) = CentralDifference(h) +with_step_size(r::Richardson, h::Real) = Richardson(h; levels=r.levels) + """ observed_order(route::DerivativeRoute, f, x, n::Integer) -> Float64 @@ -136,35 +222,34 @@ The convergence order the route actually shows on `f` at `x`, from the values at `h`, `h/2` and `h/4`: `log2(|D(h) - D(h/2)| / |D(h/2) - D(h/4)|)`. The number to look at before trusting a step, rather than a tolerance guessed in -advance. A central difference on a smooth potential returns close to 2; a value -well below that means `h` has reached the roundoff side, and a value near 0 means -`f` is not smooth at `x`. Returns `NaN` when the two differences are both zero, -which is the step being so small that the quotient stopped moving. +advance. A central difference on a smooth potential returns close to 2. Anything +else says the step or the potential is not what the route assumed, and the value +does not identify which: `h` on the roundoff side and a non-smooth `f` both land +off 2, and a kink gives exactly 1 rather than anything near 0. + +`Inf` when only the second difference vanishes and `NaN` when both do, which is +the quotient having stopped moving between halvings. Meaningful only while the successive differences are above roundoff. A route that has already reached machine precision, which [`Richardson`](@ref) does on a smooth potential, is differencing noise and reports a number with no order in it. -Defined for the step-carrying routes; [`AutoDiff`](@ref) has no step to halve. +Defined for any route reporting a [`step_size`](@ref); [`AutoDiff`](@ref) reports +`nothing` and is refused. """ function observed_order(route::DerivativeRoute, f, x, n::Integer) - return error("observed_order: $(typeof(route)) carries no step to halve.") -end -export observed_order - -_with_step(r::CentralDifference, h) = CentralDifference(h) -_with_step(r::Richardson, h) = Richardson(h; levels=r.levels) -_step(r::CentralDifference) = r.h -_step(r::Richardson) = r.h - -function observed_order(route::Union{CentralDifference,Richardson}, f, x, n::Integer) - h = _step(route) - d = [nth_derivative(_with_step(route, h / 2^k), f, x, n) for k in 0:2] + h = step_size(route) + h === nothing && error( + "observed_order: $(typeof(route)) reports no `step_size`, so there is no step " * + "to halve.", + ) + d = [nth_derivative(with_step_size(route, h / 2^k), f, x, n) for k in 0:2] a, b = abs(d[2] - d[1]), abs(d[3] - d[2]) (a == 0 && b == 0) && return NaN b == 0 && return Inf return log2(a / b) end +export observed_order # ── the genealogy, written once ────────────────────────────────────────── # @@ -202,7 +287,6 @@ thermal_derivative(Magnetization(:z), F, 0.3, Richardson(1e-2)) # tanh(0.3) ``` """ function thermal_derivative(q::AbstractQuantity, F, x::Number, route::DerivativeRoute) - route isa AutoDiff && return thermal_derivative(q, F, x) return _genealogy_derivative(q, F, x, (g, y, n) -> nth_derivative(route, g, y, n)) end @@ -210,33 +294,37 @@ end # `U`, and `U = ∂(βF)/∂β` takes `βF`. Both are a plain first derivative of the # function passed, with no sign flip, which is why they cannot go through the # generic path above. -function thermal_derivative(::SpecificHeat, U, T::Number, route::DerivativeRoute) - route isa AutoDiff && return thermal_derivative(SpecificHeat(), U, T) - return nth_derivative(route, U, T, 1) -end -function thermal_derivative(::Energy, βF, β::Number, route::DerivativeRoute) - route isa AutoDiff && return thermal_derivative(Energy(), βF, β) - return nth_derivative(route, βF, β, 1) +function thermal_derivative( + ::Union{SpecificHeat,Energy}, f, x::Number, route::DerivativeRoute +) + return nth_derivative(route, f, x, 1) end -# A single-field potential fixes only the DIAGONAL susceptibility; the same guard -# the AD path carries, so a route change cannot turn a refusal into a wrong number. -function thermal_derivative(χ::Susceptibility, F, h::Number, route::DerivativeRoute) - route isa AutoDiff && return thermal_derivative(χ, F, h) +# A single-field potential fixes only the DIAGONAL susceptibility: an off-diagonal +# component is a mixed partial in distinct field directions. Shared with the +# extension for the same reason `_genealogy_derivative` is, so a route change +# cannot turn a refusal into a wrong number. +function _susceptibility_derivative(χ::Susceptibility, F, h, nth) idx = indices(χ) all(==(idx[1]), idx) || error( "thermal_derivative: with a single-field function only the DIAGONAL χ⁽ⁿ⁾ " * "(all indices equal) is defined; got off-diagonal $(idx). Pass a multi-field " * "potential F(h⃗) and the field-component ordering.", ) - return -nth_derivative(route, F, h, response_order(χ) + 1) + return -nth(F, h, response_order(χ) + 1) +end + +function thermal_derivative(χ::Susceptibility, F, h::Number, route::DerivativeRoute) + return _susceptibility_derivative(χ, F, h, (g, y, n) -> nth_derivative(route, g, y, n)) end # The `n` the route is asked for, so a report on a third-order response halves its -# step against the third derivative and not the first. +# step against the third derivative and not the first. Callers check the edge first: +# there is no order to report for a quantity that is not a derivative of anything. function _route_order(q::AbstractQuantity) e = derivative_edge(q) - e === nothing && return 1 + e === nothing && + error("_route_order: $(typeof(q)) has no derivative_edge, so it has no order.") return derivative_order(q, e.field()) end @@ -267,17 +355,34 @@ rather than aborting the sweep, since the usual reason is a missing backend and the other rows are still the answer. """ function derivative_report(q::AbstractQuantity, F, x::Number, routes) + # Refused up front, so no row is built for a quantity that has no derivative. + # Per row it would surface as `_route_order` throwing inside the order column. + derivative_edge(q) === nothing && error( + "derivative_report: $(typeof(q)) is not a response function (no " * + "derivative_edge), so there is no derivative for a route to take.", + ) + n = _route_order(q) out = DerivativeRouteRow[] for r in routes v = try Float64(thermal_derivative(q, F, x, r)) - catch + catch e + e isa MissingRouteBackend || rethrow() NaN end - o = try - Float64(observed_order(r, F, x, _route_order(q))) - catch + # Asked whether the route has a step BEFORE calling, so the catch can stay as + # narrow as the value's. Absorbing `ErrorException` here instead would turn a + # route that reports a `step_size` and defines no `with_step_size` into the + # same NaN as one that legitimately has no step. + o = if step_size(r) === nothing NaN + else + try + Float64(observed_order(r, F, x, n)) + catch e + e isa MissingRouteBackend || rethrow() + NaN + end end push!(out, DerivativeRouteRow(r, v, o)) end diff --git a/src/relations/derivation.jl b/src/relations/derivation.jl index 5888699f..3496f560 100644 --- a/src/relations/derivation.jl +++ b/src/relations/derivation.jl @@ -117,12 +117,17 @@ end # Forward-chaining closure: keep firing any step whose inputs are all known # until nothing new is produced. Records the ordered steps actually used. # Stops early once `stop` (if given) becomes known. -function _forward_chain(known::Dict{Symbol,Any}, stop::Union{Symbol,Nothing}) +function _forward_chain( + known::Dict{Symbol,Any}, + stop::Union{Symbol,Nothing}; + avoid::Union{Symbol,Nothing}=nothing, +) used = DerivationStep[] progress = true while progress && !(stop !== nothing && haskey(known, stop)) progress = false for step in derivation_steps() + step.output === avoid && continue haskey(known, step.output) && continue v = _try_step(step, known) v === nothing && continue @@ -334,12 +339,18 @@ end # Forward-chaining closure over VariableKey nodes: fire any step whose inputs are all # known until nothing new is produced (stopping early once `stop` is known). All # "known" tests are aliasing-aware, so `known` never accumulates both β and T. -function _typed_chain!(known::Bag, extras, stop::Union{VariableKey,Nothing}) +function _typed_chain!( + known::Bag, + extras, + stop::Union{VariableKey,Nothing}; + avoid::Union{VariableKey,Nothing}=nothing, +) used = TypedStep[] progress = true while progress && !(stop !== nothing && _known(stop.type, known)) progress = false for step in typed_derivation_steps() + step.output == avoid && continue _known(step.output.type, known) && continue v = _try_typed_step(step, known, extras) v === nothing && continue @@ -778,3 +789,274 @@ function consistent(b::Bag; kwargs...) _refuse_vacuous(rows) return all(r -> r.agree, rows) end + +# ─── Every route to one target, not the first one found ────────────────── +# +# `derive` stops at the first relation that produces the target, and 66 of the +# 219 symbol-keyed outputs (38 of 84 typed ones) have more than one producing +# relation, up to nineteen. Which one runs is registry iteration order, and if +# two disagree the caller is handed a number and told nothing. +# +# The trap in checking them is circular confirmation: derive the target first +# and an intermediate built FROM it will confirm a second route trivially. So +# the closure here is built with the target held out, both as a given and as a +# derivable node, and every route is then run against data that does not contain +# it. + +""" + DerivationRouteRow + +One route to a target: the `relation`, the `inputs` it consumed, and either the +`value` it returned or the `error` it raised on the way. + +A route that raised because of the DATA is a row rather than an absence: an +impossible input makes a relation throw where it would otherwise have DISAGREED, +and dropping it silently turns the strongest evidence the data is wrong into one +fewer route to compare. A route the framework declines (`_route_declined`) is +still an absence, since it was never applicable here. +""" +struct DerivationRouteRow + relation::AbstractRelation + inputs::Vector{Any} + value::Any + error::Union{String,Nothing} + function DerivationRouteRow(rel, inputs, value, error) + # Exactly one of the two, or every consumer's `r.error === nothing` branch is + # wrong about what `value` holds. Neither set crashes `_disagreement` with a + # `MethodError` on `nothing - nothing` instead of any diagnosis. + (value === nothing) == (error === nothing) && throw( + ArgumentError( + "DerivationRouteRow: a row carries either a value or an error, not " * + "both and not neither; got value=$(repr(value)), error=$(repr(error)).", + ), + ) + return new(rel, inputs, value, error) + end +end +DerivationRouteRow(rel, inputs, value) = DerivationRouteRow(rel, inputs, value, nothing) +export DerivationRouteRow + +function Base.show(io::IO, r::DerivationRouteRow) + print(io, nameof(typeof(r.relation)), ": {", join(r.inputs, ", "), "} → ") + return print(io, r.error === nothing ? r.value : "THREW $(r.error)") +end + +# Telling "this relation cannot be applied here" from "it applied and the data +# broke it". The framework declines in exactly two shapes, both raised by +# relations/interface.jl and nowhere else: `solve:` for the affine, parametric and +# abstract-group refusals, and the untyped-slot message for a supplied value the +# caller did not give. `solve:` is therefore reserved vocabulary, pinned by a test, +# because a relation guard that borrows it disappears from the reported set. +# Everything else is the DATA, including a relation's OWN physics guard, +# which raises an `ErrorException` like `CFTEntanglementSlope: ncuts = 0 ...` and +# is a statement about the inputs. Matching the type alone would drop those, and +# `test_derivation_routes.jl` sweeps every target to pin that neither shape leaks +# into the reported set. +_row_ok(r::DerivationRouteRow) = r.error === nothing + +# The push is the same on both doors; only the presence test and the `solve` call +# above it are door-specific. Mirrors `_finite_size_scaling_row!` in finite_size.jl. +function _route_row!(rows, step, v) + return push!(rows, DerivationRouteRow(step.relation, Any[step.inputs...], v)) +end +function _route_row!(rows, step, e::Exception) + return push!( + rows, + DerivationRouteRow( + step.relation, Any[step.inputs...], nothing, sprint(showerror, e) + ), + ) +end + +function _route_declined(e) + e isa ErrorException || return false + return startswith(e.msg, "solve:") || occursin("(untyped slot)", e.msg) +end + +""" + derivation_routes(target::Symbol; knowns...) -> Vector{DerivationRouteRow} + derivation_routes(Q::Type, bag::Bag; extras...) -> Vector{DerivationRouteRow} + +EVERY relation that can produce `target` from the knowns, each with the value it +gives, where [`derive`](@ref) runs whichever one the registry reaches first. + +The target is held out of the data the routes are run against, as a given and as +a derivable node both, so a route cannot read a value that was itself derived +from the target and confirm itself. A route needing something only reachable +through the target therefore does not appear, which is the correct answer for it. + +Supplying the target is the useful case: the rows are then what the rest of the +data predicts for a number already measured. + +```julia +derivation_routes(:c; dS_dlogℓ = 0.1667, ncuts = 2) +``` +""" +function derivation_routes(target::Symbol; knowns...) + known = Dict{Symbol,Any}(pairs(knowns)) + pop!(known, target, nothing) + _forward_chain(known, nothing; avoid=target) + rows = DerivationRouteRow[] + for step in derivation_steps() + step.output === target || continue + all(v -> haskey(known, v), step.inputs) || continue + try + _route_row!( + rows, + step, + solve( + step.relation, Val(step.output); (v => known[v] for v in step.inputs)... + ), + ) + catch e + _route_declined(e) && continue + _route_row!(rows, step, e) + end + end + return rows +end + +# β and T are one quantity under two names, so holding out the target means +# holding out whichever of the pair the bag carries. +function _holdout!(known::Bag, @nospecialize(Q::Type)) + delete!(known, VariableKey(Q)) + Q === Temperature && delete!(known, VariableKey(InverseTemperature)) + Q === InverseTemperature && delete!(known, VariableKey(Temperature)) + return known +end + +function derivation_routes(@nospecialize(Q::Type), bag::Bag; extras...) + known = copy(bag) + _check_one_temperature(known) + _holdout!(known, Q) + target = VariableKey(Q) + _typed_chain!(known, values(extras), nothing; avoid=target) + rows = DerivationRouteRow[] + for step in typed_derivation_steps() + step.output == target || continue + all(k -> _known(k.type, known), step.inputs) || continue + try + _route_row!( + rows, step, solve(step.relation, step.output.type, known; extras...) + ) + catch e + _route_declined(e) && continue + _route_row!(rows, step, e) + end + end + return rows +end +export derivation_routes + +# The largest pairwise difference in a set of values, and their largest magnitude, +# to be compared as `d <= atol + rtol*m`. `nothing` for fewer than two, which is +# not agreement and must not share a sentinel with what a degenerate route returns. +function _disagreement(vs) + length(vs) < 2 && return nothing + return (maximum(abs(a - b) for a in vs, b in vs), maximum(abs, vs)) +end + +""" + derive_crosschecked(target::Symbol; atol=0, rtol=1e-8, min_routes=1, knowns...) + derive_crosschecked(Q::Type, bag::Bag; atol=0, rtol=1e-8, min_routes=1, extras...) + +[`derive`](@ref), refusing when the data reaches the target two ways that differ +by more than `atol + rtol * max|value|`, which is `isapprox`'s rule. + +`atol` defaults to `0`, so small values are judged relatively. A quantity whose +routes are genuinely noise-dominated near zero needs an `atol` saying so, rather +than a floor built into the comparison. + +One route returning a number is not evidence the data is consistent about it: +several relations can produce one target, and which one [`derive`](@ref) ran was +chosen by registry order. + +Supplying the target is the case to reach for. [`derive`](@ref) hands it straight +back without looking at anything else, and this compares it against every route +the rest of the data affords, which is the question a measured number raises. + +`min_routes` defaults to `1`: data affording no independent route is refused, +because returning a number from it is what this verb's name would otherwise be +claiming it had checked. `min_routes = 0` is the opt-out, and `2` or more is how +to demand a genuine cross-check rather than a single unopposed route. +""" +function derive_crosschecked( + target::Symbol; atol::Real=0, rtol::Real=1e-8, min_routes::Int=1, knowns... +) + rows = derivation_routes(target; knowns...) + supplied = get(Dict{Symbol,Any}(pairs(knowns)), target, nothing) + _refuse_broken_routes(":$target", rows) + _require_routes(":$target", rows, min_routes) + isempty(rows) && return supplied === nothing ? derive(target; knowns...) : supplied + _refuse_disagreement(":$target", rows, supplied, atol, rtol) + return supplied === nothing ? first(rows).value : supplied +end + +function derive_crosschecked( + @nospecialize(Q::Type), + bag::Bag; + atol::Real=0, + rtol::Real=1e-8, + min_routes::Int=1, + extras..., +) + rows = derivation_routes(Q, bag; extras...) + sup = _slot_value(Q, bag) + supplied = sup === nothing ? nothing : something(sup) + _refuse_broken_routes(string(nameof(Q)), rows) + _require_routes(string(nameof(Q)), rows, min_routes) + isempty(rows) && return supplied === nothing ? derive(Q, bag; extras...) : supplied + _refuse_disagreement(string(nameof(Q)), rows, supplied, atol, rtol) + return supplied === nothing ? first(rows).value : supplied +end +export derive_crosschecked + +function _refuse_disagreement(what, rows, supplied, atol, rtol) + vs = Any[r.value for r in rows if _row_ok(r)] + supplied === nothing || push!(vs, supplied) + # A route that returned NaN makes every difference NaN, and `isnan` as a + # "nothing to compare" sentinel would then read that as agreement. It is the + # opposite: a route degenerated and the others were never compared to it. + bad = findall(v -> v isa Number && isnan(v), vs) + isempty(bad) || error( + "derive_crosschecked: $(length(bad)) of the $(length(vs)) values for $what is " * + "NaN, so nothing was compared. A route degenerated on this data:\n " * + join(string.(rows), "\n "), + ) + dm = _disagreement(vs) + dm === nothing && return nothing + d, m = dm + d <= atol + rtol * m && return nothing + lines = string.(rows) + supplied === nothing || push!(lines, "supplied: $supplied") + return error( + "derive_crosschecked: the data reaches $what $(length(vs)) ways that differ by " * + "$d, past the $(atol + rtol * m) allowed (atol = $atol, rtol = $rtol). One of " * + "the inputs is wrong, or they are not all describing the same system:\n " * + join(lines, "\n "), + ) +end + +# A route that raised is reported before any comparison: it is a relation that +# would have disagreed, prevented from doing so by the data itself. +function _refuse_broken_routes(what, rows) + broken = [r for r in rows if !_row_ok(r)] + isempty(broken) && return nothing + return error( + "derive_crosschecked: $(length(broken)) route(s) to $what raised instead of " * + "returning a value, so they never got to disagree. Either an input is outside " * + "the relation's domain, or the relation guards the point `solve` probed the " * + "target at and needs a specialized `_solve`:\n " * + join(string.(broken), "\n "), + ) +end + +function _require_routes(what, rows, min_routes) + n = count(_row_ok, rows) + n >= min_routes && return nothing + return error( + "derive_crosschecked: $what is reached by $n independent " * + "route(s), fewer than the $min_routes asked for. The data affords no " * + "cross-check here, and a value returned from it would not have had one.", + ) +end diff --git a/src/relations/entanglement.jl b/src/relations/entanglement.jl index 5bed9f2b..33622475 100644 --- a/src/relations/entanglement.jl +++ b/src/relations/entanglement.jl @@ -103,6 +103,22 @@ is out of domain and [`CFTEntanglementChordSlope`](@ref) is the one to use. `nc dS_dlogℓ - ncuts * c / 6 end +# `_require_cuts` fires at the solver's probe of `ncuts = 0` before the slope is +# ever read, so without these the count is unreachable and the refusal names the +# probe as though it were the caller's. The guard belongs on the ANSWER. +function _solve(::CFTEntanglementSlope, ::Val{:ncuts}; dS_dlogℓ, c, _extra...) + return _solved_cuts(:CFTEntanglementSlope, dS_dlogℓ, c) +end +function _solve(::CFTEntanglementChordSlope, ::Val{:ncuts}; dS_dlogchord, c, _extra...) + return _solved_cuts(:CFTEntanglementChordSlope, dS_dlogchord, c) +end +function _solved_cuts(what::Symbol, slope, c) + iszero(c) && error("$what: c = 0 carries no slope, so it fixes no cut count.") + n = 6 * slope / c + _require_cuts(what, n) + return n +end + """ InfiniteRandomnessEntanglementSlope <: AbstractRelation @@ -665,6 +681,18 @@ Variables: `ε`, `ζ`. ε::EntanglementSpectrumLevel, ζ::CorrelationMatrixEigenvalue ) = ε - log((1 - ζ) / ζ) +# The docstring's own inverse, handed to the solver. Without it the generic +# three-point prober evaluates the kernel at ζ = 2, where `(1-ζ)/ζ = -0.5` and +# `log` raises, so `derive` called this unreachable and `derivation_routes` filed +# a broken route, for EVERY ε. The failure is a property of the probe points, not +# of the caller's data, which no classification downstream can tell apart. +function _solve(::EntanglementSpectrumCorrelation, ::Val{:ζ}; ε, _extra...) + isfinite(ε) || error( + "EntanglementSpectrumCorrelation: ζ = 1/(exp(ε) + 1) needs a finite ε; got $ε." + ) + return 1 / (exp(ε) + 1) +end + """ free_fermion_entanglement_entropy(ζ) -> Float64 diff --git a/src/relations/scaling.jl b/src/relations/scaling.jl index 5cf70949..4b3b1335 100644 --- a/src/relations/scaling.jl +++ b/src/relations/scaling.jl @@ -762,7 +762,7 @@ Variables: `dlogδTc_dlogL` (caller-computed), `ν`. # a flat slope returns an infinity whose SIGN comes from the caller's zero. function _solve(::PseudocriticalWidthScaling, ::Val{:ν}; dlogδTc_dlogL, _extra...) iszero(dlogδTc_dlogL) && error( - "solve: PseudocriticalWidthScaling has no ν at dlogδTc_dlogL = 0. A width " * + "PseudocriticalWidthScaling: no ν at dlogδTc_dlogL = 0. A width " * "that does not shift with L does not identify a correlation-length exponent.", ) return -1 / dlogδTc_dlogL diff --git a/test/core/test_conventions.jl b/test/core/test_conventions.jl index b3a90e37..8531553a 100644 --- a/test/core/test_conventions.jl +++ b/test/core/test_conventions.jl @@ -17,6 +17,27 @@ struct HalfUnits <: Convention end AbstractQAtlas.canonical_convention(::Type{ConventionProbeQuantity}) = WholeUnits() AbstractQAtlas.convert_convention(::WholeUnits, ::HalfUnits, ::Type, v) = 2v +# A parametric quantity, which is where a supertype WALK loses the declaration. +struct ParametricProbeQuantity{I} <: AbstractQuantity end +# A default parameter, as the package's own parametric quantities carry: without +# one, `test/core/test_invariants.jl`'s reflection sweep over every concrete +# `AbstractQuantity` leaf cannot build this and goes red, but only when the two +# files land in the same shard. +ParametricProbeQuantity() = ParametricProbeQuantity{:probe}() +AbstractQAtlas.canonical_convention(::Type{<:ParametricProbeQuantity}) = WholeUnits() + +# Three covers of one quantity where two are unrelated to each other and the third +# refines both. This is the shape a pairwise fold gets wrong. +abstract type AmbProbeParent <: AbstractQuantity end +struct AmbProbeSideA <: AbstractQuantity end +struct AmbProbeSideB <: AbstractQuantity end +struct AmbProbeQuantity <: AmbProbeParent end +struct UnitsA <: Convention end +struct UnitsB <: Convention end +struct UnitsC <: Convention end +const AMB_WIDE_A = Union{AmbProbeParent,AmbProbeSideA} +const AMB_WIDE_B = Union{AmbProbeParent,AmbProbeSideB} + @testset "an axis is declared per quantity, never by supertype" begin # The twelve whose ABQ definition contains a logarithm, or is an additive # combination of ones that do. @@ -45,18 +66,43 @@ AbstractQAtlas.convert_convention(::WholeUnits, ::HalfUnits, ::Type, v) = 2v @test canonical_convention(Temperature) === nothing end +why(f) = + try + f() + "" + catch e + sprint(showerror, e) + end + @testset "a declaration is refused when it claims something it cannot mean" begin - @test_throws ErrorException conventions(Float64 => Bits) - @test_throws ErrorException conventions(VonNeumannEntropy => 2) - @test_throws ErrorException conventions( - VonNeumannEntropy => Bits, VonNeumannEntropy => Nats + # Each branch has to DIAGNOSE, not merely throw: swapping the four messages + # between the four conditions leaves every `@test_throws ErrorException` green + # while handing the caller the wrong reason for their mistake. + @test occursin("not a relation variable", why(() -> conventions(Float64 => Bits))) + @test occursin("not a Convention", why(() -> conventions(VonNeumannEntropy => 2))) + @test occursin( + "duplicate key", + why(() -> conventions(VonNeumannEntropy => Bits, VonNeumannEntropy => Nats)), ) # Naming a concrete type is a claim about that type, so a type with no axis # is an error; naming its supertype is a sweep, and skips it silently. - @test_throws ErrorException conventions(TsallisEntropy => Bits) + @test occursin( + "declares no convention axis", why(() -> conventions(TsallisEntropy => Bits)) + ) @test conventions(AbstractEntanglementMeasure => Bits) isa ConventionSet end +@testset "conversion is not restricted to scalars" begin + # A bag holds whatever the calculation produced, and a spectrum or a sweep of + # region entropies is the normal shape. A conversion narrowed to `Float64` + # would ship green against every scalar fixture in this file. + cs = conventions(AbstractEntanglementMeasure => Bits) + v = bag(cs, VonNeumannEntropy() => [1.0, 2.0, 3.0])[VariableKey(VonNeumannEntropy)] + @test v ≈ [1.0, 2.0, 3.0] .* log(2) + @test convert_convention(Nats, Bits, VonNeumannEntropy, [1.0 2.0; 3.0 4.0]) ≈ + [1.0 2.0; 3.0 4.0] .* log(2) +end + @testset "lookup is most specific first" begin cs = conventions(AbstractEntanglementMeasure => Bits, VonNeumannEntropy => Nats) @test declared_convention(cs, VonNeumannEntropy) === Nats @@ -64,6 +110,50 @@ end @test declared_convention(cs, Temperature) === nothing end +@testset "a parametric quantity finds the declaration keyed on its family" begin + # The language fact the matching has to survive: a parametric type's supertype + # chain SKIPS its own family, so walking `supertype` never reaches the name a + # project keyed its declaration on. `Energy{:per_site}` is a live bag key here + # (FreeEnergyLegendre takes it), which is what makes this more than academic. + @test Energy{:per_site} <: Energy + @test supertype(Energy{:per_site}) !== Energy + cs = conventions(ParametricProbeQuantity => HalfUnits()) + @test declared_convention(cs, ParametricProbeQuantity{:a}) === HalfUnits() + @test bag(cs, ParametricProbeQuantity{:a}() => 2.5)[VariableKey( + ParametricProbeQuantity{:a} + )] == 5.0 +end + +@testset "the most specific cover is found whatever order the Dict yields" begin + # The earlier spelling of this testset paired the query with a type it is not a + # subtype of, so the ambiguity branch was never reached and the assertion held + # with the second entry deleted. These two DO both cover it. + @test AmbProbeQuantity <: AMB_WIDE_A + @test AmbProbeQuantity <: AMB_WIDE_B + @test !(AMB_WIDE_A <: AMB_WIDE_B) && !(AMB_WIDE_B <: AMB_WIDE_A) + @test AmbProbeParent <: AMB_WIDE_A && AmbProbeParent <: AMB_WIDE_B + + # Every insertion order must give the one cover that refines both. A fold that + # errors on meeting the first incomparable pair gets this right for four of the + # six orders and reports a false ambiguity for two. + entries = [AMB_WIDE_A => UnitsA(), AMB_WIDE_B => UnitsB(), AmbProbeParent => UnitsC()] + for o in + [[a, b, c] for a in 1:3 for b in 1:3 for c in 1:3 if length(unique([a, b, c])) == 3] + cs = ConventionSet(Dict{Type,Convention}(entries[i] for i in o)) + @test declared_convention(cs, AmbProbeQuantity) === UnitsC() + end + + # With no common refinement there IS no most specific cover, and guessing one by + # Dict order is the thing being refused. + genuine = ConventionSet( + Dict{Type,Convention}(AMB_WIDE_A => UnitsA(), AMB_WIDE_B => UnitsB()) + ) + @test occursin( + "no one of them a subtype of all the others", + why(() -> declared_convention(genuine, AmbProbeQuantity)), + ) +end + @testset "conversion" begin @test convert_convention(Nats, Bits, VonNeumannEntropy, 1.0) ≈ log(2) @test convert_convention(Bits, Nats, VonNeumannEntropy, log(2)) ≈ 1.0 diff --git a/test/core/test_derivative_routes.jl b/test/core/test_derivative_routes.jl index 1189e412..60f26139 100644 --- a/test/core/test_derivative_routes.jl +++ b/test/core/test_derivative_routes.jl @@ -7,9 +7,40 @@ using AbstractQAtlas using AbstractQAtlas: _route_order, _genealogy_derivative using Test: @test, @test_throws, @testset +why(f) = + try + f() + "" + catch e + sprint(showerror, e) + end + F(h) = -log(2cosh(h)) # M = -F'(h) = tanh(h) Φ(T) = -T * log(2cosh(1 / T)) # S = -Φ'(T) -kinked(x) = x < 0 ? x^2 : 2x^2 # derivative discontinuous at 0 +kinked(x) = x < 0 ? x^2 : 2x^2 # f and f' continuous at 0, f'' jumps 2 -> 4 + +# A genealogy that roots somewhere other than a thermodynamic potential, so +# `_genealogy_derivative`'s root guard has something that can actually fire it. +struct RootProbeParent <: AbstractQuantity end +struct RootProbeQuantity <: AbstractQuantity end +function AbstractQAtlas.derivative_edge(::Type{RootProbeQuantity}) + return DerivativeEdge(RootProbeParent, Temperature) +end + +# A route that breaks the contract: it reports a step and gives no way to change +# one. `observed_order` needs both, so this must surface rather than become a NaN. +struct BrokenStepRoute <: DerivativeRoute + h::Float64 +end +function AbstractQAtlas.nth_derivative(r::BrokenStepRoute, f, x, n::Integer) + return nth_derivative(CentralDifference(r.h), f, x, n) +end +AbstractQAtlas.step_size(r::BrokenStepRoute) = r.h + +# A route naming a backend that IS loaded, so the fallback must not blame the +# package for a method the route itself never defined. +struct LoadedBackendRoute <: DerivativeRoute end +AbstractQAtlas.backend_package(::LoadedBackendRoute) = :LinearAlgebra @testset "a finite-difference route reaches the closed form" begin x = 0.3 @@ -43,15 +74,73 @@ end # In the asymptotic regime a central difference shows its nominal order. @test observed_order(CentralDifference(1e-2), F, 0.3, 1) ≈ 2 atol = 0.05 # A step small enough to be dominated by cancellation does not, even though - # its error happens to be smaller. Stated as "does not show the nominal - # order" rather than a number, because the value there is roundoff and would - # be a different number on another machine (NaN included, when the successive - # differences both vanish). - @test !(observed_order(CentralDifference(1e-9), F, 0.3, 1) > 1.9) - # The control the diagnostic needs: a function it should FAIL on. A kink at - # the evaluation point leaves the quotient first-order, exactly. + # its error happens to be smaller. Written as "outside the band around 2" + # rather than "below 2": within one decade of this step the quotient also + # returns `Inf` (only the second difference vanishes) and `NaN` (both do), and + # `!(o > 1.9)` is false for `Inf`, so that spelling would fail on a step 12% + # away in log space with no platform difference needed. + @test !(1.9 <= observed_order(CentralDifference(1e-9), F, 0.3, 1) <= 2.1) + @test !(1.9 <= observed_order(CentralDifference(1e-11), F, 0.3, 1) <= 2.1) + # The control the diagnostic needs: a function it should FAIL on. The second + # derivative jumps at 0, which leaves the quotient first-order exactly: the + # h^2 term of the central difference does not cancel, so D(h) = h/2. @test observed_order(CentralDifference(1e-3), kinked, 0.0, 1) ≈ 1 atol = 1e-9 - @test_throws ErrorException observed_order(AutoDiff(), F, 0.3, 1) + # AutoDiff reports no step, so there is nothing to halve. The message has to + # say that, not the generic "no method". + msg = try + observed_order(AutoDiff(), F, 0.3, 1) + "" + catch e + sprint(showerror, e) + end + @test occursin("step_size", msg) + @test step_size(AutoDiff()) === nothing + @test step_size(CentralDifference(1e-3)) == 1e-3 + @test step_size(Richardson(1e-2)) == 1e-2 + # Richardson carries a step too, so it must not fall out of `observed_order` + # the way a closed Union over the routes that happened to exist would drop a + # third one. Its value is noise once the route is at machine precision, so the + # claim is that it RUNS and returns a number, not what the number is. + @test isfinite(observed_order(Richardson(1e-1), F, 0.3, 1)) || + isnan(observed_order(Richardson(1e-1), F, 0.3, 1)) +end + +@testset "a route that carries a step must say so, not be named in a Union" begin + # The contract is `step_size` + `with_step_size`, so a future route gets the + # honest refusal instead of the false claim that it carries no step. + @test with_step_size(CentralDifference(1e-2), 1e-3) == CentralDifference(1e-3) + @test with_step_size(Richardson(1e-2; levels=4), 1e-3) == Richardson(1e-3; levels=4) + @test_throws ErrorException with_step_size(AutoDiff(), 1e-3) +end + +@testset "Richardson's levels is a knob, not a decoration" begin + x, exact = 0.3, tanh(0.3) + errs = [ + abs( + thermal_derivative(Magnetization(:z), F, x, Richardson(1e-1; levels=L)) - exact + ) for L in 2:5 + ] + # STRICTLY falling: a `levels` that is ignored gives four equal errors, and + # `issorted` counts ties as sorted, so it alone would pass that. + @test all(errs[i] > errs[i + 1] for i in 1:(length(errs) - 1)) + @test errs[1] / errs[end] > 1e3 + # And a bad `levels` is refused rather than silently behaving as the default. + @test_throws ArgumentError Richardson(1e-2; levels=1) + @test_throws ArgumentError Richardson(1e-2; levels=0) +end + +@testset "the root guard can fire" begin + # Every shipped derivative_edge chains to FreeEnergy or GrandPotential, so + # without a quantity rooted elsewhere this guard is unreachable and deleting it + # changes nothing. + @test potential_root(RootProbeQuantity()) === RootProbeParent + msg = try + thermal_derivative(RootProbeQuantity(), F, 0.3, CentralDifference(1e-3)) + "" + catch e + sprint(showerror, e) + end + @test occursin("RootProbeParent", msg) end @testset "a route change cannot turn a refusal into a number" begin @@ -81,6 +170,17 @@ end for r in (CentralDifference(1e-3), Richardson(1e-2)) @test nth_derivative(r, F, 0.3, 0) == F(0.3) end + # `nth_derivative` is exported, so a caller can reach the AutoDiff route + # directly rather than through `thermal_derivative`. Without a backend it has + # to name the routes that need none, not fall through to the generic + # "no method" of an unrecognised route. + msg = try + nth_derivative(AutoDiff(), F, 0.3, 2) + "" + catch e + sprint(showerror, e) + end + @test isempty(msg) || occursin("CentralDifference", msg) end @testset "a report compares routes instead of trusting one" begin @@ -95,4 +195,83 @@ end @test isnan(ad.value) || isapprox(ad.value, tanh(0.3); atol=1e-12) @test isnan(ad.order) @test all(r -> isapprox(r.value, tanh(0.3); atol=1e-3), rows[1:2]) + # The row has to name the route it ran, or a report that always stored the + # first route would read the same. + @test [r.route for r in rows] == [CentralDifference(1e-2), Richardson(1e-2), AutoDiff()] + # Only a missing backend is absorbed into a NaN row. A quantity with no + # genealogy edge, and a guard the route itself raises, both propagate: a NaN + # there would read as "install a package" for a mistake no package fixes. + # Pinned by MESSAGE: `_route_order` has its own guard one line later, so a + # type-only assertion passes whichever of the two fired. + @test occursin( + "is not a response function", + why( + () -> derivative_report(PartitionFunction(), F, 0.3, (CentralDifference(1e-2),)) + ), + ) + @test_throws ErrorException derivative_report( + Susceptibility(:x, :y), F, 0.3, (CentralDifference(1e-2),) + ) +end + +@testset "the backend route's method belongs to the extension alone" begin + # Defining `nth_derivative(::AutoDiff, ...)` in BOTH the package and the + # extension is a method overwrite, which makes the extension fail to + # precompile while every test here stays green, because the fallback path + # still loads. So the package must own no method for that signature, and the + # "which package" answer lives in a trait instead. + owned = [ + m for m in methods(nth_derivative) if + m.module === AbstractQAtlas && m.sig.parameters[2] === AutoDiff + ] + @test isempty(owned) + @test backend_package(AutoDiff()) === :ForwardDiff + @test backend_package(CentralDifference(1e-3)) === nothing + @test backend_package(Richardson(1e-2)) === nothing + # Without the extension the fallback has to say WHICH package, not "no method". + msg = try + nth_derivative(AutoDiff(), F, 0.3, 1) + "" + catch e + sprint(showerror, e) + end + @test isempty(msg) || occursin("ForwardDiff", msg) +end + +@testset "a route that breaks the step contract is not absorbed as a NaN" begin + # The order column used to catch `ErrorException` as well as + # `MissingRouteBackend`, so a route reporting a `step_size` with no + # `with_step_size` produced the same NaN as one that legitimately has no step. + # Different mistakes, and only one of them is the caller's. + @test occursin( + "with_step_size", why(() -> observed_order(BrokenStepRoute(1e-2), F, 0.3, 1)) + ) + @test_throws ErrorException derivative_report( + Magnetization(:z), F, 0.3, (BrokenStepRoute(1e-2),) + ) + # A route declaring no step is still a quiet NaN, the one case the column may + # absorb, and its message says which half is missing. + @test isnan(only(derivative_report(Magnetization(:z), F, 0.3, (AutoDiff(),))).order) + @test occursin( + "reports no `step_size`", why(() -> observed_order(AutoDiff(), F, 0.3, 1)) + ) +end + +@testset "a missing method is not blamed on a package that is loaded" begin + # LinearAlgebra is loaded by this package, so reaching the fallback means the + # ROUTE is incomplete. Telling its author to reinstall points away from that. + msg = why(() -> nth_derivative(LoadedBackendRoute(), F, 0.3, 1)) + @test occursin("is loaded, but no method matched", msg) + @test !occursin("which is not loaded", msg) +end + +@testset "MissingRouteBackend cannot name an extension that does not exist" begin + # It is exported, so an extension author can construct it. Built for a route + # with no `backend_package` it used to render "needs the nothing extension". + @test_throws ArgumentError MissingRouteBackend(CentralDifference(1e-3)) + @test_throws ArgumentError MissingRouteBackend(Richardson(1e-2)) + e = MissingRouteBackend(AutoDiff()) + @test e isa Exception + @test occursin("MissingRouteBackend:", sprint(showerror, e)) + @test occursin("ForwardDiff", sprint(showerror, e)) end diff --git a/test/core/test_examples.jl b/test/core/test_examples.jl new file mode 100644 index 00000000..21a02b6c --- /dev/null +++ b/test/core/test_examples.jl @@ -0,0 +1,31 @@ +# The shipped examples, run and checked. +# +# An example nobody runs is worse than none: it reads as a promise and rots +# silently. This includes the script rather than restating it, so the two cannot +# disagree, and asserts against the EXACT values the physics fixes rather than +# against the digits the script happens to print. + +using AbstractQAtlas +using Test: @test, @testset + +include(joinpath(@__DIR__, "..", "..", "examples", "critical_entanglement.jl")) +using .CriticalEntanglementExample: block_entropy, central_charges + +@testset "examples/critical_entanglement.jl" begin + cs = central_charges() + # The oracle is `c = 1` for a free-fermion chain, not a recorded output. The + # gap is the finite-block correction at ℓ = 16, 64, which is why this is 1e-3. + @test cs.nats ≈ 1 atol = 1e-3 + @test cs.bits_declared ≈ 1 atol = 1e-3 + # Undeclared bits are wrong by exactly the base, and the test says which + # number it is rather than just "not 1": a wrong answer that happens to miss + # by something else would otherwise pass for the same reason. + @test cs.bits_undeclared ≈ 1 / log(2) atol = 1e-3 + @test !isapprox(cs.bits_undeclared, 1; atol=1e-3) + # Declaring the convention has to change the answer, not merely be accepted. + @test cs.bits_declared ≈ cs.nats atol = 1e-12 + + # The entropies themselves are a measurement, so pin that they grow with the + # block rather than pinning digits a BLAS version can move. + @test block_entropy(64) > block_entropy(16) > block_entropy(8) > 0 +end diff --git a/test/ext/test_derivative_routes_ad.jl b/test/ext/test_derivative_routes_ad.jl index 21e38589..1ce70b46 100644 --- a/test/ext/test_derivative_routes_ad.jl +++ b/test/ext/test_derivative_routes_ad.jl @@ -49,3 +49,24 @@ end @test rows[3].order ≈ 2 atol = 0.05 @test maximum(r -> abs(r.value - rows[1].value), rows) < 1e-4 end + +@testset "the off-diagonal guard is one guard, not two copies" begin + # It used to be pasted in both paths, and the two copies had already drifted + # to different wording. Same message from both is what says they are shared. + off = Susceptibility(:x, :y) + ad = try + thermal_derivative(off, Fad, 0.3, AutoDiff()) + "" + catch e + sprint(showerror, e) + end + fd = try + thermal_derivative(off, Fad, 0.3, Richardson(1e-2)) + "" + catch e + sprint(showerror, e) + end + @test !isempty(ad) + @test ad == fd + @test occursin("DIAGONAL", ad) +end diff --git a/test/relations/test_derivation_routes.jl b/test/relations/test_derivation_routes.jl new file mode 100644 index 00000000..9d17b899 --- /dev/null +++ b/test/relations/test_derivation_routes.jl @@ -0,0 +1,332 @@ +# Every route to a target, rather than the first one the registry reaches. +# +# 66 of the 219 symbol-keyed outputs have more than one producing relation (38 +# of 84 typed ones), up to nineteen, and `derive` runs whichever comes first. +# +# The trap this file exists to pin is not the disagreement, it is the vacuous +# agreement: derive the target first and the closure manufactures its own inputs +# from it, so several routes hand the supplied number back and read as +# confirmation of something nothing independent touched. + +using AbstractQAtlas +using AbstractQAtlas: _disagreement, derivation_steps, typed_derivation_steps + +why(f) = + try + f() + "" + catch e + sprint(showerror, e) + end +using Test: @test, @test_throws, @testset + +slope(c) = 2 * c / 6 # CFTEntanglementSlope: dS/dlnℓ = ncuts·c/6, ncuts = 2 + +@testset "every INDEPENDENT route appears, and only those" begin + rows = derivation_routes(:c; dS_dlogℓ=slope(0.5), ncuts=2) + # One measurement, one route. Letting the chain derive `c` first and then + # rebuild the chord slope and the halved-chain difference FROM it returns + # three rows that all read 0.5, which is this measurement counted three + # times and would pass any agreement test put to it. + @test length(rows) == 1 + @test nameof(typeof(first(rows).relation)) === :CFTEntanglementSlope + @test first(rows).value ≈ 0.5 + # So asking for two independent routes here must fail: there is only one. + @test_throws ErrorException derive_crosschecked( + :c; dS_dlogℓ=slope(0.5), ncuts=2, min_routes=2 + ) +end + +@testset "the target is held out, so a route cannot confirm itself" begin + # With only `c` and `ncuts`, the closure CAN manufacture the slopes: they are + # derivable from the supplied `c`. + reach = derivable(; c=0.9, ncuts=2) + @test :dS_dlogℓ in reach + @test :dS_dlogchord in reach + # And yet no route is reported, because every one of them would be reading a + # value built from the target. Three rows all returning 0.9 is what this + # emptiness replaces. + @test isempty(derivation_routes(:c; c=0.9, ncuts=2)) + # The same graph does reach `c` once something independent is supplied, so the + # emptiness above is the holdout and not an unreachable target. + @test !isempty(derivation_routes(:c; c=0.9, dS_dlogℓ=slope(0.5), ncuts=2)) +end + +@testset "a contradiction the first-route solver returns anyway" begin + # `derive` hands back the supplied value without consulting anything. + @test derive(:c; c=0.9, dS_dlogℓ=slope(0.5), ncuts=2) == 0.9 + msg = try + derive_crosschecked(:c; c=0.9, dS_dlogℓ=slope(0.5), ncuts=2) + "" + catch e + sprint(showerror, e) + end + @test occursin("differ by 0.4", msg) # the size, not just that it threw + @test occursin("CFTEntanglementSlope", msg) # the route is named + @test occursin("supplied: 0.9", msg) # and so is the value it contradicts + # Consistent data passes, and returns the supplied value. + @test derive_crosschecked(:c; c=0.5, dS_dlogℓ=slope(0.5), ncuts=2) == 0.5 + # No supplied target: the routes are the answer. + @test derive_crosschecked(:c; dS_dlogℓ=slope(0.5), ncuts=2) ≈ 0.5 +end + +@testset "min_routes defaults to demanding that a cross-check happened" begin + # Data affording no independent route is refused by default: returning 0.9 from + # it is what the verb's name would otherwise be claiming it had checked. + @test_throws ErrorException derive_crosschecked(:c; c=0.9, ncuts=2) + @test derive_crosschecked(:c; c=0.9, ncuts=2, min_routes=0) == 0.9 # opt-out + # One route is enough to compare a supplied value against. + @test derive_crosschecked(:c; c=0.5, dS_dlogℓ=slope(0.5), ncuts=2) == 0.5 + # And two can be demanded where one is not evidence. + @test_throws ErrorException derive_crosschecked( + :c; c=0.5, dS_dlogℓ=slope(0.5), ncuts=2, min_routes=2 + ) +end + +@testset "a route that raised is reported, not dropped" begin + # Z must be positive: it is a sum of Boltzmann weights. `FreeEnergyFromZ` needs + # log(Z) and throws, which is the relation being PREVENTED from disagreeing. + # Dropping it silently leaves one route and a clean "cross-checked" answer. + impossible = bag( + PartitionFunction => -2.0, + InverseTemperature => 1.0, + Energy(:per_site) => -0.4, + ThermalEntropy => 0.3, + ) + rows = derivation_routes(FreeEnergy, impossible) + @test length(rows) == 2 + threw = only(r for r in rows if r.error !== nothing) + @test nameof(typeof(threw.relation)) === :FreeEnergyFromZ + @test occursin("DomainError", threw.error) + @test threw.value === nothing + # The row DISPLAYS as a failure. Printing `r.value` unconditionally would show + # `nothing`, and the refusal message is built from `string.(rows)`. + @test occursin("THREW", string(threw)) + @test occursin("DomainError", string(threw)) + # Neither state and both states are unrepresentable, so no consumer's + # `error === nothing` branch can be wrong about what `value` holds. + @test_throws ArgumentError DerivationRouteRow(threw.relation, Any[], nothing, nothing) + @test_throws ArgumentError DerivationRouteRow(threw.relation, Any[], 1.0, "boom") + msg = try + derive_crosschecked(FreeEnergy, impossible) + "" + catch e + sprint(showerror, e) + end + @test occursin("never got to disagree", msg) + @test occursin("FreeEnergyFromZ", msg) + # A relation merely declining to be solved for a slot is NOT a broken route, or + # every ordinary call would refuse. The consistent bag still passes. + β, Z, U = 0.8, 3.0, 0.4 + F = -log(Z) / β + ok = bag( + PartitionFunction => Z, + InverseTemperature => β, + Energy(:per_site) => U, + ThermalEntropy => β * (U - F), + ) + @test all(r -> r.error === nothing, derivation_routes(FreeEnergy, ok)) + @test derive_crosschecked(FreeEnergy, ok; min_routes=2) ≈ F +end + +@testset "a NaN route is a degeneration, not an agreement" begin + # Every difference against NaN is NaN, so an `isnan` short-circuit meant for + # "fewer than two values" would read a degenerate route as agreement. + msg = try + derive_crosschecked(:c; c=NaN, dS_dlogℓ=slope(0.5), ncuts=2) + "" + catch e + sprint(showerror, e) + end + @test occursin("NaN", msg) + @test occursin("nothing was compared", msg) +end + +@testset "the typed door holds the target out the same way" begin + β, Z, U = 0.8, 3.0, 0.4 + F = -log(Z) / β + S = β * (U - F) # FreeEnergyLegendre: F = U - S/β + # Two genuinely independent routes to F: through Z, and through the Legendre + # transform. They agree here. + good = bag( + PartitionFunction => Z, + InverseTemperature => β, + Energy(:per_site) => U, + ThermalEntropy => S, + ) + rows = derivation_routes(FreeEnergy, good) + @test length(rows) == 2 + @test Set(nameof(typeof(r.relation)) for r in rows) == + Set([:FreeEnergyFromZ, :FreeEnergyLegendre]) + @test all(r -> r.value ≈ F, rows) + @test derive_crosschecked(FreeEnergy, good; min_routes=2) ≈ F + + # Break the entropy. `derive` still returns the RIGHT number, off the other + # route, and says nothing about the input that contradicts it. + bad = bag( + PartitionFunction => Z, + InverseTemperature => β, + Energy(:per_site) => U, + ThermalEntropy => S + 0.5, + ) + @test derive(FreeEnergy, bad) ≈ F + @test_throws ErrorException derive_crosschecked(FreeEnergy, bad) + + # The holdout, on the typed door. `FreeEnergyLegendre` produces BOTH `F` and + # `S`, so with F supplied the closure can manufacture the S that gives F back. + circular = bag(FreeEnergy => F, Energy(:per_site) => U, InverseTemperature => β) + @test VariableKey(ThermalEntropy) in derivable(circular) # it could be built + @test isempty(derivation_routes(FreeEnergy, circular)) # and it is not used + + # β and T are one quantity under two names, so a bag carrying β must not let + # `KelvinRelation` hand T back through a Peltier coefficient built from it. + alias = bag(InverseTemperature => 0.5, Thermopower(:x, :x) => 3.0) + @test VariableKey(PeltierCoefficient{(:x, :x)}) in derivable(alias) + @test isempty(derivation_routes(Temperature, alias)) +end + +@testset "the tolerance is isapprox's, not a floor that goes absolute below one" begin + # Two routes returning +1e-12 and -1e-12 is a sign flip. Dividing the + # difference by `max(maximum(abs, vs), 1)` reports it as 2e-12 and passes it + # at any sane rtol, which is why the comparison is `d <= atol + rtol*m`. + d, m = _disagreement([1e-12, -1e-12]) + @test (d, m) == (2e-12, 1e-12) + @test !(d <= 0 + 1e-8 * m) + # A genuine agreement to 1e-9 relative still passes. + d2, m2 = _disagreement([1.0, 1.0 + 1e-9]) + @test d2 <= 0 + 1e-8 * m2 + # Fewer than two values is not agreement. It gets `nothing`, not a NaN that a + # degenerate route would also produce. + @test _disagreement([1.0]) === nothing + @test _disagreement(Any[]) === nothing + + # `atol` is how a caller says their routes are noise-dominated, and it is the + # only thing that lets the broken-entropy bag through. + β, Z, U = 0.8, 3.0, 0.4 + F = -log(Z) / β + bad = bag( + PartitionFunction => Z, + InverseTemperature => β, + Energy(:per_site) => U, + ThermalEntropy => β * (U - F) + 0.5, + ) + @test_throws ErrorException derive_crosschecked(FreeEnergy, bad) + @test derive_crosschecked(FreeEnergy, bad; atol=1.0) ≈ F +end + +@testset "multiple producers is the common case the section comment claims" begin + # The comment above `derivation_routes` carries measured counts. Pinned as + # floors rather than exact numbers: adding a relation that produces an + # already-produced output raises them, and a floor does not rot for that. + function multi(steps) + d = Dict{Any,Set{Any}}() + for st in steps + push!(get!(d, st.output, Set{Any}()), nameof(typeof(st.relation))) + end + return count(v -> length(v) > 1, values(d)), maximum(length, values(d)) + end + n_sym, max_sym = multi(derivation_steps()) + n_typed, max_typed = multi(typed_derivation_steps()) + @test n_sym >= 60 + @test max_sym >= 15 + @test n_typed >= 35 + @test max_typed >= 12 +end + +@testset "the declining/broken split holds across the whole registry" begin + # `_route_declined` separates "this relation cannot be applied here" (skip) + # from "it applied and the data broke it" (report). Getting it wrong in either + # direction is invisible: too narrow and a broken input is silently dropped, + # too wide and ordinary calls refuse. + # + # Both framework declination shapes, on fixtures that actually produce them. + @test isempty(derivation_routes(:β; C=1.0, var_E=1.0)) # solve: not affine + @test isempty( # untyped slot absent + derivation_routes( + Temperature, bag(InverseTemperature => 0.5, Thermopower(:x, :x) => 3.0) + ), + ) + # A relation's own physics guard is about the DATA and must be reported, not + # skipped. `ncuts = 0` makes the residual independent of `c`. + ncz = derivation_routes(:c; dS_dlogℓ=0.1667, ncuts=0) + @test length(ncz) == 1 + @test occursin("ncuts = 0", only(ncz).error) + + # The sweep: hand every target's inputs the same nonsense value and check that + # nothing the framework MEANT as a declination ends up in the reported set. A + # seventh declination site added to interface.jl and not registered here fails + # this, which is the whole point of the assertion. + leaked = String[] + total = 0 + for t in unique(s.output for s in derivation_steps()) + ins = unique( + vcat([collect(s.inputs) for s in derivation_steps() if s.output === t]...) + ) + rows = try + derivation_routes(t; (i => 0.7 for i in ins)...) + catch + continue + end + total += length(rows) + for r in rows + r.error === nothing && continue + (startswith(r.error, "solve:") || occursin("(untyped slot)", r.error)) && + push!(leaked, r.error) + end + end + @test total > 300 # the sweep reached the registry, not two rows + @test isempty(leaked) +end + +@testset "a guard on the solver's probe must not make a slot unreachable" begin + # `solve` probes its target at 0, 1, 2. A relation guarding one of those points + # then refuses for EVERY input, and the refusal names a value the caller never + # supplied. Both shapes below did that, and the fix is the closed form each + # relation already had in prose. + # + # `EntanglementSpectrumCorrelation` is `ε - log((1-ζ)/ζ)`: at the probe ζ = 2 + # the argument of `log` is -0.5. Its docstring carried `ζ = 1/(e^ε + 1)` already. + @test derive(:ζ; ε=0.7) ≈ 1 / (exp(0.7) + 1) + @test derive_crosschecked(:ζ; ε=0.7) ≈ 1 / (exp(0.7) + 1) + @test isempty([r for r in derivation_routes(:ζ; ε=0.7) if r.error !== nothing]) + + # The slope relations guard `ncuts = 0`, which is exactly the first probe. + @test derive(:ncuts; dS_dlogℓ=slope(0.5), c=0.5) ≈ 2 + @test derive_crosschecked(:ncuts; dS_dlogℓ=slope(0.5), c=0.5) ≈ 2 + @test isempty([ + r for + r in derivation_routes(:ncuts; dS_dlogℓ=slope(0.5), c=0.5) if r.error !== nothing + ]) + # The guard still holds, now read off the ANSWER rather than the probe, and it + # reaches the caller. Through `derive` it does not: that verb drops the route + # and reports the target unreachable, which is the difference this layer makes. + @test occursin("ncuts = 0", only(derivation_routes(:ncuts; dS_dlogℓ=0.0, c=0.5)).error) + @test occursin("ncuts = 0", why(() -> derive_crosschecked(:ncuts; dS_dlogℓ=0.0, c=0.5))) + @test occursin("c = 0", why(() -> derive_crosschecked(:ncuts; dS_dlogℓ=0.3, c=0.0))) + @test occursin("not reachable", why(() -> derive(:ncuts; dS_dlogℓ=0.0, c=0.5))) +end + +@testset "`solve:` is the framework's vocabulary, and only the framework's" begin + # `_route_declined` classifies by message prefix, so a relation guard that + # borrows `solve:` is silently SKIPPED instead of reported. The sweep above + # cannot see that direction: a skipped route leaves no row to inspect. One + # guard in `scaling.jl` had borrowed it, and only a reader caught that. + # + # Read off the source, because the claim is about what future authors write. + dir = joinpath(@__DIR__, "..", "..", "src", "relations") + offenders = String[] + for f in sort(readdir(dir)) + endswith(f, ".jl") && f != "interface.jl" || continue + text = replace(read(joinpath(dir, f), String), "\r\n" => "\n") + for (i, line) in enumerate(split(text, "\n")) + occursin("\"solve: ", line) && push!(offenders, "$f:$i") + end + end + @test isempty(offenders) + # The positive control: the vocabulary IS used, in the one file that owns it. + iface = replace( + read(joinpath(@__DIR__, "..", "..", "src", "relations", "interface.jl"), String), + "\r\n" => "\n", + ) + @test count(!isempty, [m.match for m in eachmatch(r"\"solve: ", iface)]) >= 5 +end