Skip to content

Revised Belief Propagation #26

Draft
jack-dunham wants to merge 41 commits intoITensor:mainfrom
jack-dunham:bp
Draft

Revised Belief Propagation #26
jack-dunham wants to merge 41 commits intoITensor:mainfrom
jack-dunham:bp

Conversation

@jack-dunham
Copy link
Contributor

@jack-dunham jack-dunham commented Nov 25, 2025

This PR express belief propagation in terms of the new interface based on AlgorithmsInterface.jl and the included AlgorithmsInterfaceExtensions.jl library.

Note, the tolerance for the second belief propagation test has been increased from 1.0e-14 to 1.0e-12, however this test was not passing reliable before the changes in this PR anyway.

JoeyT1994 and others added 29 commits January 6, 2026 09:55
Introduce `BeliefPropagationProblem` wrapper to hold the cache and the
error `diff` field.

Also simplifies some kwargs wrangling.
Also includes some fixes to the way `TensorNetwork` types are
constructed based on index structure.
…instead of trying to operate on existing graphs

The reason for this is:
- One only cares about the edges of the input graph
- A simple graph cannot be used as it "forgets" its edge names resulting
in recursion
- As shown with `TensorNetwork`, removing edges may not always be
defined.
This was caused by the change to the `cache` being backed by a directed
graph.
@github-actions
Copy link
Contributor

github-actions bot commented Jan 9, 2026

Your PR no longer requires formatting changes. Thank you for your contribution!

@jack-dunham jack-dunham force-pushed the bp branch 2 times, most recently from a20ddcb to ac74a7f Compare February 11, 2026 14:44
@codecov
Copy link

codecov bot commented Feb 13, 2026

Codecov Report

❌ Patch coverage is 50.27778% with 179 lines in your changes missing coverage. Please review.
✅ Project coverage is 57.88%. Comparing base (6f3521f) to head (fef588d).

Files with missing lines Patch % Lines
...eliefpropagation/abstractbeliefpropagationcache.jl 36.78% 55 Missing ⚠️
src/tensornetwork.jl 14.06% 55 Missing ⚠️
src/beliefpropagation/beliefpropagationcache.jl 52.50% 38 Missing ⚠️
src/abstracttensornetwork.jl 22.72% 17 Missing ⚠️
src/beliefpropagation/beliefpropagationproblem.jl 89.21% 11 Missing ⚠️
src/LazyNamedDimsArrays/symbolicnameddimsarray.jl 0.00% 2 Missing ⚠️
src/LazyNamedDimsArrays/lazynameddimsarray.jl 0.00% 1 Missing ⚠️

❗ There is a different number of reports uploaded between BASE (6f3521f) and HEAD (fef588d). Click for more details.

HEAD has 5 uploads less than BASE
Flag BASE (6f3521f) HEAD (fef588d)
6 1
Additional details and impacted files
@@            Coverage Diff             @@
##             main      #26      +/-   ##
==========================================
- Coverage   63.58%   57.88%   -5.71%     
==========================================
  Files          21       21              
  Lines         736     1040     +304     
==========================================
+ Hits          468      602     +134     
- Misses        268      438     +170     
Flag Coverage Δ
docs 0.00% <0.00%> (ø)

Flags with carried forward coverage won't be shown. Click here to find out more.

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:
  • ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

messages = incoming_messages(bp_cache, vertex)
state = factors(bp_cache, vertex)

return (reduce(*, messages) * reduce(*, state))[]
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This should be done with a contraction sequence finder. Probably with a default like "greedy" but the ability to overload the sequence kwargs would be ideal.

message_type(bpc::AbstractBeliefPropagationCache) = message_type(typeof(bpc))
message_type(::Type{<:AbstractBeliefPropagationCache{<:Any, <:Any, ED}}) where {ED} = ED

function free_energy(bp_cache::AbstractBeliefPropagationCache)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe this should be called logscalar ? Technically it returns log(scalar(bp_cache)) where scalar(bp_cache) is the BP approximation of the underlying network. The free energy is really -kBTlog(Z) so free_energy might be a misnomer.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Agreed, math terminology would be preferred over physics terminology.

return rv
end

function induced_subgraph_bpcache(graph, subvertices)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is there a reason this can't ne called induced_subgraph and dispatch based on the type?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It is a convention we are using where induced_subgraph_bpcache is a function defining the implementation of induced_subgraph for bpcache (perhaps this suffix should be renamed to something more explicit). Then if one defines a type NotAnAbstractBeliefPropagationCache that doesn't subtype AbstractBeliefPropagationCache but should behave like an AbstractBeliefPropagationCache (at least for subgraph purposes) one can define:

function induced_subgraph_from_vertices(cache::NotAnAbstractBeliefPropagationCache, subvertices)
    return induced_subgraph_bpcache(cache, subvertices)
end

Note the appropriate function to overload is usually induced_subgraph_from_vertices, which assumes the subvertices argument as already been canonized via to_vertices. In any case, overloading induced_subgraph directly can be lead to annoying method ambiguities as induced_subgraph is a function in Graphs.jl that we have no control over. Matt and I have discussed refining this aspect of the code to move away from induced_subgraph entirely and just having the call stack be:

subgraph(graph, vertices) = subgraph_<implementation>(graph, to_vertices(graph, vertices))

which would simplify things greatly.

using NamedGraphs.PartitionedGraphs: quotientvertices

@kwdef struct StopWhenConverged <: AI.StoppingCriterion
tol::Float64 = 0.0
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can we use Number instead of Float64. It will be good to be able to use Float32 (and other precisions) arithmetic if we wish.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We can use AbstractFloat.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Just a note that probably we should parametrize the type (we can constrain the parameter to AbstractFloat or Real).


if algorithm.normalize
# TODO: use `sum` not `norm`
message_norm = LinearAlgebra.norm(state.iterate)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

message_norm = sum(state.iterate) .

This is important for stability in the presence of complex numbers (avoids blow up of non-PSDness of messages).

maxiter = is_tree(cache) ? 1 : nothing,
tol = -Inf,
message_diff_function = if tol > -Inf
(m1, m2) -> norm(m1 / norm(m1) - m2 / norm(m2))
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This message diff is currently dependent on phase, i.e. m1 = -m2 will count as non-converged. 1 - abs2(dot(m1, m2)) might be better, with m1 and m2 normalized.

@assert U * LinearAlgebra.diagm(S) * V ≈ W
id = [1.0 0.0; 0.0 1.0]
set!(sqrt_Ws, e, id)
set!(sqrt_Ws, reverse(e), U * LinearAlgebra.diagm(S) * V)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Currently the SVD is not really playing a role here.
set!(sqrt_Ws, reverse(e), W) would work. Or, probably best, take the square root via the SVD and set U*sqrt(S) on e and sqrt(S)*V on reverse(e).

end

bpc = BeliefPropagationCache(tn)
bpc = ITensorNetworksNext.beliefpropagation(bpc; maxiter = 1)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It would be nice if the maxiter default could be set to 1 if graph(bpc) ==1 and unset otherwise.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@JoeyT1994 what do you mean by graph(bpc) == 1? Shouldn't the check be that it is a tree?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes sorry I mean is_tree(graph(bpc))

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

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants