From d0dda4702ebf2ec336591e359876e827f4c533d5 Mon Sep 17 00:00:00 2001 From: Joey Date: Mon, 2 Oct 2023 16:34:26 -0400 Subject: [PATCH 01/37] Added option to specify order of BP updates. Also set of edges for doing BP updates. Updated BP gauging example to reflect change --- examples/gauging/gauging_itns.jl | 36 ++++++++++++++----- src/beliefpropagation/beliefpropagation.jl | 19 +++++++--- .../sqrt_beliefpropagation.jl | 10 +++++- 3 files changed, 52 insertions(+), 13 deletions(-) diff --git a/examples/gauging/gauging_itns.jl b/examples/gauging/gauging_itns.jl index 6c30cfd4..c390e9ac 100644 --- a/examples/gauging/gauging_itns.jl +++ b/examples/gauging/gauging_itns.jl @@ -45,7 +45,10 @@ end """Bring an ITN into the Vidal gauge, various methods possible. Result is timed""" function benchmark_state_gauging( - ψ::ITensorNetwork; mode="BeliefPropagation", no_iterations=50 + ψ::ITensorNetwork; + mode="BeliefPropagation", + no_iterations=50, + BP_update_order::String="Parallel", ) s = siteinds(ψ) @@ -67,8 +70,9 @@ function benchmark_state_gauging( if mode == "BeliefPropagation" times_iters[i] = @elapsed mts, _ = belief_propagation_iteration( - ψψ, mts; contract_kwargs=(; alg="exact") + ψψ, mts; contract_kwargs=(; alg="exact"), update_order=BP_update_order ) + times_gauging[i] = @elapsed ψ, bond_tensors = vidal_gauge(ψinit, mts) elseif mode == "Eager" times_iters[i] = @elapsed ψ, bond_tensors, mts = eager_gauging(ψ, bond_tensors, mts) @@ -95,22 +99,38 @@ s = siteinds("S=1/2", g) no_iterations = 30 BPG_simulation_times, BPG_Cs = benchmark_state_gauging(ψ; no_iterations) +BPG_sequential_simulation_times, BPG_sequential_Cs = benchmark_state_gauging( + ψ; no_iterations, BP_update_order="Sequential" +) Eager_simulation_times, Eager_Cs = benchmark_state_gauging(ψ; mode="Eager", no_iterations) SU_simulation_times, SU_Cs = benchmark_state_gauging(ψ; mode="SU", no_iterations) epsilon = 1e-6 + println( - "Time for BPG to reach C < epsilon was " * + "Time for BPG (with parallel updates) to reach C < epsilon was " * string(BPG_simulation_times[findfirst(x -> x < 0, BPG_Cs .- epsilon)]) * - " seconds", + " seconds. No iters was " * + string(findfirst(x -> x < 0, BPG_Cs .- epsilon)), ) println( - "Time for Eager to reach C < epsilon was " * + "Time for BPG (with sequential updates) to reach C < epsilon was " * + string( + BPG_sequential_simulation_times[findfirst(x -> x < 0, BPG_sequential_Cs .- epsilon)] + ) * + " seconds. No iters was " * + string(findfirst(x -> x < 0, BPG_sequential_Cs .- epsilon)), +) + +println( + "Time for Eager Gauging to reach C < epsilon was " * string(Eager_simulation_times[findfirst(x -> x < 0, Eager_Cs .- epsilon)]) * - " seconds", + " seconds. No iters was " * + string(findfirst(x -> x < 0, Eager_Cs .- epsilon)), ) println( - "Time for SU to reach C < epsilon was " * + "Time for SU Gauging to reach C < epsilon was " * string(SU_simulation_times[findfirst(x -> x < 0, SU_Cs .- epsilon)]) * - " seconds", + " seconds. No iters was " * + string(findfirst(x -> x < 0, SU_Cs .- epsilon)), ) diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index dc23d4a5..f0b5fa2b 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -81,17 +81,25 @@ function belief_propagation_iteration( mts::DataGraph; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), compute_norm=false, + update_order::String="Parallel", + es=edges(mts), ) new_mts = copy(mts) + if update_order != "Parallel" && update_order != "Sequential" + error( + "Specified update order is not currently implemented. Choose Parallel or Sequential." + ) + end + incoming_mts = update_order == "Parallel" ? mts : new_mts c = 0 - es = edges(mts) for e in es environment_tensornetworks = ITensorNetwork[ - mts[e_in] for e_in in setdiff(boundary_edges(mts, [src(e)]; dir=:in), [reverse(e)]) + incoming_mts[e_in] for + e_in in setdiff(boundary_edges(incoming_mts, [src(e)]; dir=:in), [reverse(e)]) ] new_mts[src(e) => dst(e)] = update_message_tensor( - tn, mts[src(e)], environment_tensornetworks; contract_kwargs + tn, incoming_mts[src(e)], environment_tensornetworks; contract_kwargs ) if compute_norm @@ -111,10 +119,13 @@ function belief_propagation( contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), niters=20, target_precision::Union{Float64,Nothing}=nothing, + update_order::String="Parallel", ) compute_norm = target_precision == nothing ? false : true for i in 1:niters - mts, c = belief_propagation_iteration(tn, mts; contract_kwargs, compute_norm) + mts, c = belief_propagation_iteration( + tn, mts; contract_kwargs, compute_norm, update_order + ) if compute_norm && c <= target_precision println( "Belief Propagation finished. Reached a canonicalness of " * diff --git a/src/beliefpropagation/sqrt_beliefpropagation.jl b/src/beliefpropagation/sqrt_beliefpropagation.jl index 15662367..c7036138 100644 --- a/src/beliefpropagation/sqrt_beliefpropagation.jl +++ b/src/beliefpropagation/sqrt_beliefpropagation.jl @@ -26,11 +26,19 @@ end function sqrt_belief_propagation_iteration( tn::ITensorNetwork, sqrt_mts::DataGraph; + update_order::String="Parallel", + es=edges(sqrt_mts), + # compute_norm=false, ) new_sqrt_mts = copy(sqrt_mts) + if update_order != "Parallel" && update_order != "Sequential" + error( + "Specified update order is not currently implemented. Choose Parallel or Sequential." + ) + end + incoming_mts = update_order == "Parallel" ? mts : new_mts c = 0.0 - es = edges(sqrt_mts) for e in es environment_tensornetworks = ITensorNetwork[ sqrt_mts[e_in] for From 54604ac49de3bc95676ca68dacda7b192322b8e8 Mon Sep 17 00:00:00 2001 From: Joey Date: Mon, 2 Oct 2023 16:39:33 -0400 Subject: [PATCH 02/37] Changed BP ordering to sequential for a test --- test/test_belief_propagation.jl | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/test/test_belief_propagation.jl b/test/test_belief_propagation.jl index 6bc1eb01..bd4720cf 100644 --- a/test/test_belief_propagation.jl +++ b/test/test_belief_propagation.jl @@ -40,7 +40,13 @@ ITensors.disable_warn_order() Z = partition(ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ))))) mts = message_tensors(Z) - mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact"), target_precision=1e-8) + mts = belief_propagation( + ψψ, + mts; + contract_kwargs=(; alg="exact"), + target_precision=1e-8; + update_order="Sequential", + ) numerator_network = approx_network_region( ψψ, mts, [(v, 1)]; verts_tn=ITensorNetwork(ITensor[apply(op("Sz", s[v]), ψ[v])]) From ccc79462ea8a47598d08bef7c3c8cb23a5ea3636 Mon Sep 17 00:00:00 2001 From: Joey Date: Tue, 3 Oct 2023 09:56:26 -0400 Subject: [PATCH 03/37] No capitalisation. Fixed bug in test_belief_propagation.jl --- examples/gauging/gauging_itns.jl | 6 +++--- src/beliefpropagation/beliefpropagation.jl | 10 +++++----- src/beliefpropagation/sqrt_beliefpropagation.jl | 11 ++++++----- test/test_belief_propagation.jl | 4 ++-- 4 files changed, 16 insertions(+), 15 deletions(-) diff --git a/examples/gauging/gauging_itns.jl b/examples/gauging/gauging_itns.jl index c390e9ac..d561295c 100644 --- a/examples/gauging/gauging_itns.jl +++ b/examples/gauging/gauging_itns.jl @@ -48,7 +48,7 @@ function benchmark_state_gauging( ψ::ITensorNetwork; mode="BeliefPropagation", no_iterations=50, - BP_update_order::String="Parallel", + BP_update_order::String="parallel", ) s = siteinds(ψ) @@ -100,7 +100,7 @@ no_iterations = 30 BPG_simulation_times, BPG_Cs = benchmark_state_gauging(ψ; no_iterations) BPG_sequential_simulation_times, BPG_sequential_Cs = benchmark_state_gauging( - ψ; no_iterations, BP_update_order="Sequential" + ψ; no_iterations, BP_update_order="sequential" ) Eager_simulation_times, Eager_Cs = benchmark_state_gauging(ψ; mode="Eager", no_iterations) SU_simulation_times, SU_Cs = benchmark_state_gauging(ψ; mode="SU", no_iterations) @@ -129,7 +129,7 @@ println( string(findfirst(x -> x < 0, Eager_Cs .- epsilon)), ) println( - "Time for SU Gauging to reach C < epsilon was " * + "Time for SU Gauging (with sequential updates) to reach C < epsilon was " * string(SU_simulation_times[findfirst(x -> x < 0, SU_Cs .- epsilon)]) * " seconds. No iters was " * string(findfirst(x -> x < 0, SU_Cs .- epsilon)), diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index f0b5fa2b..5a5b85b0 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -81,16 +81,16 @@ function belief_propagation_iteration( mts::DataGraph; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), compute_norm=false, - update_order::String="Parallel", + update_order::String="parallel", es=edges(mts), ) new_mts = copy(mts) - if update_order != "Parallel" && update_order != "Sequential" + if update_order != "parallel" && update_order != "sequential" error( - "Specified update order is not currently implemented. Choose Parallel or Sequential." + "Specified update order is not currently implemented. Choose parallel or sequential." ) end - incoming_mts = update_order == "Parallel" ? mts : new_mts + incoming_mts = update_order == "parallel" ? mts : new_mts c = 0 for e in es environment_tensornetworks = ITensorNetwork[ @@ -119,7 +119,7 @@ function belief_propagation( contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), niters=20, target_precision::Union{Float64,Nothing}=nothing, - update_order::String="Parallel", + update_order::String="parallel", ) compute_norm = target_precision == nothing ? false : true for i in 1:niters diff --git a/src/beliefpropagation/sqrt_beliefpropagation.jl b/src/beliefpropagation/sqrt_beliefpropagation.jl index c7036138..5091b506 100644 --- a/src/beliefpropagation/sqrt_beliefpropagation.jl +++ b/src/beliefpropagation/sqrt_beliefpropagation.jl @@ -5,12 +5,13 @@ function sqrt_belief_propagation( tn::ITensorNetwork, mts::DataGraph; niters=20, + update_order::String="parallel", # target_precision::Union{Float64,Nothing}=nothing, ) # compute_norm = target_precision == nothing ? false : true sqrt_mts = sqrt_message_tensors(tn, mts) for i in 1:niters - sqrt_mts, c = sqrt_belief_propagation_iteration(tn, sqrt_mts) #; compute_norm) + sqrt_mts, c = sqrt_belief_propagation_iteration(tn, sqrt_mts; update_order) #; compute_norm) # if compute_norm && c <= target_precision # println( # "Belief Propagation finished. Reached a canonicalness of " * @@ -26,18 +27,18 @@ end function sqrt_belief_propagation_iteration( tn::ITensorNetwork, sqrt_mts::DataGraph; - update_order::String="Parallel", + update_order::String="parallel", es=edges(sqrt_mts), # compute_norm=false, ) new_sqrt_mts = copy(sqrt_mts) - if update_order != "Parallel" && update_order != "Sequential" + if update_order != "parallel" && update_order != "sequential" error( - "Specified update order is not currently implemented. Choose Parallel or Sequential." + "Specified update order is not currently implemented. Choose parallel or sequential." ) end - incoming_mts = update_order == "Parallel" ? mts : new_mts + incoming_mts = update_order == "parallel" ? mts : new_mts c = 0.0 for e in es environment_tensornetworks = ITensorNetwork[ diff --git a/test/test_belief_propagation.jl b/test/test_belief_propagation.jl index bd4720cf..d42b204d 100644 --- a/test/test_belief_propagation.jl +++ b/test/test_belief_propagation.jl @@ -44,8 +44,8 @@ ITensors.disable_warn_order() ψψ, mts; contract_kwargs=(; alg="exact"), - target_precision=1e-8; - update_order="Sequential", + target_precision=1e-8, + update_order="sequential", ) numerator_network = approx_network_region( From 6fdd8d85b2ed423cdcf80d4029f6af8b57b3a9ee Mon Sep 17 00:00:00 2001 From: Joey Date: Tue, 3 Oct 2023 10:05:05 -0400 Subject: [PATCH 04/37] es -> edges --- src/beliefpropagation/sqrt_beliefpropagation.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/beliefpropagation/sqrt_beliefpropagation.jl b/src/beliefpropagation/sqrt_beliefpropagation.jl index 5091b506..3f3b31d2 100644 --- a/src/beliefpropagation/sqrt_beliefpropagation.jl +++ b/src/beliefpropagation/sqrt_beliefpropagation.jl @@ -28,7 +28,7 @@ function sqrt_belief_propagation_iteration( tn::ITensorNetwork, sqrt_mts::DataGraph; update_order::String="parallel", - es=edges(sqrt_mts), + edges=Graphs.edges(sqrt_mts), # compute_norm=false, ) @@ -40,7 +40,7 @@ function sqrt_belief_propagation_iteration( end incoming_mts = update_order == "parallel" ? mts : new_mts c = 0.0 - for e in es + for e in edges environment_tensornetworks = ITensorNetwork[ sqrt_mts[e_in] for e_in in setdiff(boundary_edges(sqrt_mts, [src(e)]; dir=:in), [reverse(e)]) @@ -58,7 +58,7 @@ function sqrt_belief_propagation_iteration( # c += 0.5 * norm(LHS - RHS) # end end - return new_sqrt_mts, c / (length(es)) + return new_sqrt_mts, c / (length(edges)) end function update_sqrt_message_tensor( From 6e0f1529d799c7ad4ea41af384ca9e3473b8e9fd Mon Sep 17 00:00:00 2001 From: Joey Date: Tue, 3 Oct 2023 12:39:59 -0400 Subject: [PATCH 05/37] Bug Fix --- src/beliefpropagation/sqrt_beliefpropagation.jl | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/beliefpropagation/sqrt_beliefpropagation.jl b/src/beliefpropagation/sqrt_beliefpropagation.jl index 3f3b31d2..e0a2ce56 100644 --- a/src/beliefpropagation/sqrt_beliefpropagation.jl +++ b/src/beliefpropagation/sqrt_beliefpropagation.jl @@ -38,16 +38,16 @@ function sqrt_belief_propagation_iteration( "Specified update order is not currently implemented. Choose parallel or sequential." ) end - incoming_mts = update_order == "parallel" ? mts : new_mts + incoming_sqrt_mts = update_order == "parallel" ? sqrt_mts : new_sqrt_mts c = 0.0 for e in edges environment_tensornetworks = ITensorNetwork[ - sqrt_mts[e_in] for - e_in in setdiff(boundary_edges(sqrt_mts, [src(e)]; dir=:in), [reverse(e)]) + incoming_sqrt_mts[e_in] for + e_in in setdiff(boundary_edges(incoming_sqrt_mts, [src(e)]; dir=:in), [reverse(e)]) ] new_sqrt_mts[src(e) => dst(e)] = update_sqrt_message_tensor( - tn, sqrt_mts[src(e)], environment_tensornetworks; + tn, incoming_sqrt_mts[src(e)], environment_tensornetworks; ) # if compute_norm From ff699d4b0c5040f6138114dcf43235f1e29de4a5 Mon Sep 17 00:00:00 2001 From: Joey Date: Fri, 6 Oct 2023 15:42:08 -0400 Subject: [PATCH 06/37] Update_order -> update_sequence --- src/beliefpropagation/beliefpropagation.jl | 16 ++++++++-------- src/beliefpropagation/sqrt_beliefpropagation.jl | 10 +++++----- 2 files changed, 13 insertions(+), 13 deletions(-) diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index 5a5b85b0..c20402a0 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -81,18 +81,18 @@ function belief_propagation_iteration( mts::DataGraph; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), compute_norm=false, - update_order::String="parallel", - es=edges(mts), + update_sequence::String="parallel", + edges=Graphs.edges(mts), ) new_mts = copy(mts) - if update_order != "parallel" && update_order != "sequential" + if update_sequence != "parallel" && update_sequence != "sequential" error( "Specified update order is not currently implemented. Choose parallel or sequential." ) end - incoming_mts = update_order == "parallel" ? mts : new_mts + incoming_mts = update_sequence == "parallel" ? mts : new_mts c = 0 - for e in es + for e in edges environment_tensornetworks = ITensorNetwork[ incoming_mts[e_in] for e_in in setdiff(boundary_edges(incoming_mts, [src(e)]; dir=:in), [reverse(e)]) @@ -110,7 +110,7 @@ function belief_propagation_iteration( c += 0.5 * norm(denseblocks(LHS) - denseblocks(RHS)) end end - return new_mts, c / (length(es)) + return new_mts, c / (length(edges)) end function belief_propagation( @@ -119,12 +119,12 @@ function belief_propagation( contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), niters=20, target_precision::Union{Float64,Nothing}=nothing, - update_order::String="parallel", + update_sequence::String="parallel", ) compute_norm = target_precision == nothing ? false : true for i in 1:niters mts, c = belief_propagation_iteration( - tn, mts; contract_kwargs, compute_norm, update_order + tn, mts; contract_kwargs, compute_norm, update_sequence ) if compute_norm && c <= target_precision println( diff --git a/src/beliefpropagation/sqrt_beliefpropagation.jl b/src/beliefpropagation/sqrt_beliefpropagation.jl index e0a2ce56..07b074da 100644 --- a/src/beliefpropagation/sqrt_beliefpropagation.jl +++ b/src/beliefpropagation/sqrt_beliefpropagation.jl @@ -5,13 +5,13 @@ function sqrt_belief_propagation( tn::ITensorNetwork, mts::DataGraph; niters=20, - update_order::String="parallel", + update_sequence::String="parallel", # target_precision::Union{Float64,Nothing}=nothing, ) # compute_norm = target_precision == nothing ? false : true sqrt_mts = sqrt_message_tensors(tn, mts) for i in 1:niters - sqrt_mts, c = sqrt_belief_propagation_iteration(tn, sqrt_mts; update_order) #; compute_norm) + sqrt_mts, c = sqrt_belief_propagation_iteration(tn, sqrt_mts; update_sequence) #; compute_norm) # if compute_norm && c <= target_precision # println( # "Belief Propagation finished. Reached a canonicalness of " * @@ -27,18 +27,18 @@ end function sqrt_belief_propagation_iteration( tn::ITensorNetwork, sqrt_mts::DataGraph; - update_order::String="parallel", + update_sequence::String="parallel", edges=Graphs.edges(sqrt_mts), # compute_norm=false, ) new_sqrt_mts = copy(sqrt_mts) - if update_order != "parallel" && update_order != "sequential" + if update_sequence != "parallel" && update_sequence != "sequential" error( "Specified update order is not currently implemented. Choose parallel or sequential." ) end - incoming_sqrt_mts = update_order == "parallel" ? sqrt_mts : new_sqrt_mts + incoming_sqrt_mts = update_sequence == "parallel" ? sqrt_mts : new_sqrt_mts c = 0.0 for e in edges environment_tensornetworks = ITensorNetwork[ From 9ad1672e410553bfdc52c49ffbcd3b5ea9b7d490 Mon Sep 17 00:00:00 2001 From: Joey Date: Fri, 6 Oct 2023 16:24:41 -0400 Subject: [PATCH 07/37] Bug Fix --- test/test_belief_propagation.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/test_belief_propagation.jl b/test/test_belief_propagation.jl index d42b204d..a5e54556 100644 --- a/test/test_belief_propagation.jl +++ b/test/test_belief_propagation.jl @@ -45,7 +45,7 @@ ITensors.disable_warn_order() mts; contract_kwargs=(; alg="exact"), target_precision=1e-8, - update_order="sequential", + update_sequence="sequential", ) numerator_network = approx_network_region( From 8efe3011f209c29f8b003dd07801bafd98be9b00 Mon Sep 17 00:00:00 2001 From: Joey Date: Mon, 23 Oct 2023 10:59:40 -0400 Subject: [PATCH 08/37] Working on optimal orders --- examples/belief_propagation/bpexample.jl | 2 + src/beliefpropagation/beliefpropagation.jl | 12 +++--- src/utils.jl | 12 ++++++ test/test_belief_propagation.jl | 48 ++++++++++++++++++---- 4 files changed, 60 insertions(+), 14 deletions(-) diff --git a/examples/belief_propagation/bpexample.jl b/examples/belief_propagation/bpexample.jl index 6fc96ed1..a86a0f53 100644 --- a/examples/belief_propagation/bpexample.jl +++ b/examples/belief_propagation/bpexample.jl @@ -4,6 +4,7 @@ using Metis using ITensorNetworks using Random using SplitApplyCombine +using NamedGraphs using ITensorNetworks: belief_propagation, @@ -35,6 +36,7 @@ function main() ) mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact")) + numerator_network = approx_network_region( ψψ, mts, [(v, 1)]; verts_tn=ITensorNetwork([apply(op("Sz", s[v]), ψ[v])]) ) diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index c20402a0..bbd73dd6 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -81,8 +81,8 @@ function belief_propagation_iteration( mts::DataGraph; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), compute_norm=false, - update_sequence::String="parallel", - edges=Graphs.edges(mts), + update_sequence::String="sequential", + edges = edge_update_order(undirected_graph(underlying_graph(mts))), ) new_mts = copy(mts) if update_sequence != "parallel" && update_sequence != "sequential" @@ -97,7 +97,6 @@ function belief_propagation_iteration( incoming_mts[e_in] for e_in in setdiff(boundary_edges(incoming_mts, [src(e)]; dir=:in), [reverse(e)]) ] - new_mts[src(e) => dst(e)] = update_message_tensor( tn, incoming_mts[src(e)], environment_tensornetworks; contract_kwargs ) @@ -120,17 +119,16 @@ function belief_propagation( niters=20, target_precision::Union{Float64,Nothing}=nothing, update_sequence::String="parallel", + edges = edge_update_order(undirected_graph(underlying_graph(mts))) ) compute_norm = target_precision == nothing ? false : true for i in 1:niters mts, c = belief_propagation_iteration( - tn, mts; contract_kwargs, compute_norm, update_sequence + tn, mts; contract_kwargs, compute_norm, update_sequence, edges ) if compute_norm && c <= target_precision println( - "Belief Propagation finished. Reached a canonicalness of " * - string(c) * - " after $i iterations. ", + "BP converged to desired precision after $i iterations.", ) break end diff --git a/src/utils.jl b/src/utils.jl index 9a1be885..ca08d04d 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -25,3 +25,15 @@ function line_to_tree(line::Vector) end return [line_to_tree(line[1:(end - 1)]), line[end]] end + +#Find an optimal ordering of the edges in an undirected graph +function edge_update_order(g::AbstractNamedGraph) + es = [] + for v in vertices(g) + new_es = reverse(reverse.(edges(bfs_tree(g, v)))) + push!(es, setdiff(new_es, es)...) + end + + @assert Set(es) == Set(vcat(edges(g), reverse.(edges(g)))) + return es +end \ No newline at end of file diff --git a/test/test_belief_propagation.jl b/test/test_belief_propagation.jl index a5e54556..4c75bf90 100644 --- a/test/test_belief_propagation.jl +++ b/test/test_belief_propagation.jl @@ -22,7 +22,7 @@ ITensors.disable_warn_order() @testset "belief_propagation" begin - #FIRST TEST SINGLE SITE ON AN MPS, SHOULD BE EXACT + #First test on an MPS, should be exact g_dims = (1, 6) g = named_grid(g_dims) s = siteinds("S=1/2", g) @@ -44,7 +44,7 @@ ITensors.disable_warn_order() ψψ, mts; contract_kwargs=(; alg="exact"), - target_precision=1e-8, + niters = 1, update_sequence="sequential", ) @@ -56,8 +56,40 @@ ITensors.disable_warn_order() @test abs.(bp_sz - exact_sz) <= 1e-14 - #NOW TEST TWO_SITE_EXPEC TAKING ON THE PARTITION FUNCTION OF THE RECTANGULAR ISING. SHOULD BE REASONABLE AND - #INDEPENDENT OF INIT CONDITIONS, FOR SMALL BETA + #Now test on a tree, should also be exact + g = named_comb_tree((6,6)) + s = siteinds("S=1/2", g) + χ = 2 + Random.seed!(1564) + ψ = randomITensorNetwork(s; link_space=χ) + + ψψ = ψ ⊗ prime(dag(ψ); sites=[]) + + v = (1, 3) + + Oψ = copy(ψ) + Oψ[v] = apply(op("Sz", s[v]), ψ[v]) + exact_sz = contract_inner(Oψ, ψ) / contract_inner(ψ, ψ) + + Z = partition(ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ))))) + mts = message_tensors(Z) + mts = belief_propagation( + ψψ, + mts; + contract_kwargs=(; alg="exact"), + niters = 1, + update_sequence="sequential", + ) + + numerator_network = approx_network_region( + ψψ, mts, [(v, 1)]; verts_tn=ITensorNetwork(ITensor[apply(op("Sz", s[v]), ψ[v])]) + ) + denominator_network = approx_network_region(ψψ, mts, [(v, 1)]) + bp_sz = contract(numerator_network)[] / contract(denominator_network)[] + + @test abs.(bp_sz - exact_sz) <= 1e-14 + + #Now test two-site expec taking on the partition function of the Ising model. Not exact, but close g_dims = (3, 4) g = named_grid(g_dims) s = IndsNetwork(g; link_space=2) @@ -84,7 +116,7 @@ ITensors.disable_warn_order() @test abs.(bp_szsz - actual_szsz) <= 0.05 - #TEST FORMING OF A TWO SITE RDM. JUST CHECK THAT IS HAS THE CORRECT SIZE, TRACE AND IS PSD + #Test forming a two-site RDM. Check it has the correct size, trace 1 and is PSD g_dims = (4, 4) g = named_grid(g_dims) s = siteinds("S=1/2", g) @@ -118,7 +150,7 @@ ITensors.disable_warn_order() @test all(>=(0), real(eigs)) && all(==(0), imag(eigs)) @test size(rdm) == (2^length(vs), 2^length(vs)) - #TEST MORE ADVANCED BP WITH ITENSORNETWORK MESSAGE TENSORS. IN THIS CASE IT SHOULD BE LIKE BOUNDARY MPS + #Test more advanced block BP with MPS message tensors on a grid g_dims = (5, 4) g = named_grid(g_dims) s = siteinds("S=1/2", g) @@ -144,6 +176,8 @@ ITensors.disable_warn_order() contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, cutoff=1e-16, maxdim ), + niters = 1, + update_sequence="sequential", ) numerator_network = approx_network_region(ψψ, mts, [v]; verts_tn=ITensorNetwork(ψOψ[v])) @@ -155,4 +189,4 @@ ITensors.disable_warn_order() contract_boundary_mps(ψOψ; cutoff=1e-16) / contract_boundary_mps(ψψ; cutoff=1e-16) @test abs.(bp_sz - exact_sz) <= 1e-5 -end +end \ No newline at end of file From 1cf1daf5a3a0f04d43c4b15477d6395ce18cd231 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Mon, 23 Oct 2023 11:10:52 -0400 Subject: [PATCH 09/37] Added custom edge specification on vidal_itn_isometries --- src/gauging.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/gauging.jl b/src/gauging.jl index 194dce98..0f8d499a 100644 --- a/src/gauging.jl +++ b/src/gauging.jl @@ -173,10 +173,10 @@ function symmetric_to_vidal_gauge( end """Function to measure the 'isometries' of a state in the Vidal Gauge""" -function vidal_itn_isometries(ψ::ITensorNetwork, bond_tensors::DataGraph) +function vidal_itn_isometries(ψ::ITensorNetwork, bond_tensors::DataGraph; edges = vcat(edges(ψ), reverse.(edges(ψ)))) isometries = DataGraph{vertextype(ψ),ITensor,ITensor}(directed_graph(underlying_graph(ψ))) - for e in vcat(edges(ψ), reverse.(edges(ψ))) + for e in edges vsrc, vdst = src(e), dst(e) ψv = copy(ψ[vsrc]) for vn in setdiff(neighbors(ψ, vsrc), [vdst]) From ffa9d5cfa5ca5d4da0ee38976b9959cfc40bb8d2 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Mon, 23 Oct 2023 11:12:34 -0400 Subject: [PATCH 10/37] Added custom edge specification on vidal_gauge --- src/gauging.jl | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/src/gauging.jl b/src/gauging.jl index 0f8d499a..2ef8dd0c 100644 --- a/src/gauging.jl +++ b/src/gauging.jl @@ -17,11 +17,12 @@ function vidal_gauge( bond_tensors::DataGraph; eigen_message_tensor_cutoff=10 * eps(real(scalartype(ψ))), regularization=10 * eps(real(scalartype(ψ))), + edges = NamedGraphs.edges(ψ), svd_kwargs..., ) ψ_vidal = copy(ψ) - for e in edges(ψ_vidal) + for e in edges vsrc, vdst = src(e), dst(e) ψvsrc, ψvdst = ψ_vidal[vsrc], ψ_vidal[vdst] @@ -173,7 +174,7 @@ function symmetric_to_vidal_gauge( end """Function to measure the 'isometries' of a state in the Vidal Gauge""" -function vidal_itn_isometries(ψ::ITensorNetwork, bond_tensors::DataGraph; edges = vcat(edges(ψ), reverse.(edges(ψ)))) +function vidal_itn_isometries(ψ::ITensorNetwork, bond_tensors::DataGraph; edges = vcat(NamedGraphs.edges(ψ), reverse.(NamedGraphs.edges(ψ)))) isometries = DataGraph{vertextype(ψ),ITensor,ITensor}(directed_graph(underlying_graph(ψ))) for e in edges From 6f8d83209ffa7cfa97fcc9b178605d642b524fa5 Mon Sep 17 00:00:00 2001 From: Joey Date: Mon, 23 Oct 2023 11:24:03 -0400 Subject: [PATCH 11/37] New Testing for Sequences --- examples/belief_propagation/bpsequences.jl | 91 ++++++++++++++++++++++ 1 file changed, 91 insertions(+) create mode 100644 examples/belief_propagation/bpsequences.jl diff --git a/examples/belief_propagation/bpsequences.jl b/examples/belief_propagation/bpsequences.jl new file mode 100644 index 00000000..6a9797e0 --- /dev/null +++ b/examples/belief_propagation/bpsequences.jl @@ -0,0 +1,91 @@ +using Compat +using ITensors +using Metis +using ITensorNetworks +using Random +using SplitApplyCombine +using NamedGraphs + +using ITensorNetworks: + belief_propagation, + approx_network_region, + contract_inner, + message_tensors, + nested_graph_leaf_vertices + +function main() + + g = named_comb_tree((6,6)) + s = siteinds("S=1/2", g) + χ = 4 + + Random.seed!(5467) + + ψ = randomITensorNetwork(s; link_space=χ) + ψψ = ψ ⊗ prime(dag(ψ); sites=[]) + + #Initial message tensors for BP + mts_init = message_tensors( + ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ)))) + ) + + println("First testing out a comb tree. Random network with bond dim $χ") + + #Now test out various sequences + print("Parallel updates (sequence is irrelevant): ") + belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "parallel") + print("Sequential updates (sequence is default edge list of the message tensors): ") + belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "sequential", edges = edges(mts_init)) + print("Sequential updates (sequence is our custom sequence finder): ") + belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "sequential") + + g = named_grid((6,6)) + s = siteinds("S=1/2", g) + χ = 2 + + Random.seed!(5467) + + ψ = randomITensorNetwork(s; link_space=χ) + ψψ = ψ ⊗ prime(dag(ψ); sites=[]) + + #Initial message tensors for BP + mts_init = message_tensors( + ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ)))) + ) + + println("Now testing out a 6x6 grid. Random network with bond dim $χ") + + #Now test out various sequences + print("Parallel updates (sequence is irrelevant): ") + belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "parallel") + print("Sequential updates (sequence is default edge list of the message tensors): ") + belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "sequential", edges = edges(mts_init)) + print("Sequential updates (sequence is our custom sequence finder): ") + belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "sequential") + + g = NamedGraphs.hexagonal_lattice_graph(4,4) + s = siteinds("S=1/2", g) + χ = 3 + + Random.seed!(5467) + + ψ = randomITensorNetwork(s; link_space=χ) + ψψ = ψ ⊗ prime(dag(ψ); sites=[]) + + #Initial message tensors for BP + mts_init = message_tensors( + ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ)))) + ) + + println("Now testing out a 4 x 4 hexagonal lattice. Random network with bond dim $χ") + + #Now test out various sequences + print("Parallel updates (sequence is irrelevant): ") + belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "parallel") + print("Sequential updates (sequence is default edge list of the message tensors): ") + belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "sequential", edges = edges(mts_init)) + print("Sequential updates (sequence is our custom sequence finder): ") + belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "sequential") +end + +main() From 0832ad552101736c170a2017243765ff8e80c361 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Mon, 23 Oct 2023 12:56:28 -0400 Subject: [PATCH 12/37] Functions for optimal edge order --- examples/gauging/gauging_itns.jl | 3 +-- src/utils.jl | 24 ++++++++++++++++-------- 2 files changed, 17 insertions(+), 10 deletions(-) diff --git a/examples/gauging/gauging_itns.jl b/examples/gauging/gauging_itns.jl index d561295c..82f6cf8e 100644 --- a/examples/gauging/gauging_itns.jl +++ b/examples/gauging/gauging_itns.jl @@ -4,7 +4,6 @@ using Metis using ITensorNetworks using Random using SplitApplyCombine -using ProfileView using ITensorNetworks: message_tensors, @@ -70,7 +69,7 @@ function benchmark_state_gauging( if mode == "BeliefPropagation" times_iters[i] = @elapsed mts, _ = belief_propagation_iteration( - ψψ, mts; contract_kwargs=(; alg="exact"), update_order=BP_update_order + ψψ, mts; contract_kwargs=(; alg="exact"), update_sequence=BP_update_order ) times_gauging[i] = @elapsed ψ, bond_tensors = vidal_gauge(ψinit, mts) diff --git a/src/utils.jl b/src/utils.jl index ca08d04d..1d821dc8 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -26,14 +26,22 @@ function line_to_tree(line::Vector) return [line_to_tree(line[1:(end - 1)]), line[end]] end -#Find an optimal ordering of the edges in an undirected graph -function edge_update_order(g::AbstractNamedGraph) - es = [] - for v in vertices(g) - new_es = reverse(reverse.(edges(bfs_tree(g, v)))) - push!(es, setdiff(new_es, es)...) +function edge_update_order(g) + forests = NamedGraphs.build_forest_cover(g) + edges = NamedEdge[] + for forest in forests + trees = NamedGraph[forest[vs] for vs in connected_components(forest)] + for tree in trees + push!(edges, tree_edge_update_order(tree)...) + end end - @assert Set(es) == Set(vcat(edges(g), reverse.(edges(g)))) - return es + return edges +end + +#Find an optimal ordering of the edges in a tree +function tree_edge_update_order(g::AbstractNamedGraph; root_vertex = first(keys(sort(eccentricities(g); rev=true)))) + @assert is_tree(g) + es = post_order_dfs_edges(g, root_vertex) + return vcat(es, reverse(reverse.(es))) end \ No newline at end of file From dfec1e271e05a044a89af01de53183e9234048b6 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Thu, 26 Oct 2023 10:57:14 -0400 Subject: [PATCH 13/37] Further changes --- examples/gauging/gauging_itns.jl | 3 ++- src/beliefpropagation/beliefpropagation.jl | 12 +++++++----- src/beliefpropagation/sqrt_beliefpropagation.jl | 6 +++--- src/gauging.jl | 3 ++- test/test_belief_propagation.jl | 10 +++++----- test/test_gauging.jl | 2 +- 6 files changed, 20 insertions(+), 16 deletions(-) diff --git a/examples/gauging/gauging_itns.jl b/examples/gauging/gauging_itns.jl index 82f6cf8e..0f65dc8e 100644 --- a/examples/gauging/gauging_itns.jl +++ b/examples/gauging/gauging_itns.jl @@ -22,6 +22,7 @@ using ITensorNetworks: using NamedGraphs using NamedGraphs: add_edges!, rem_vertex!, hexagonal_lattice_graph using Graphs +using MKL """Eager Gauging""" function eager_gauging(ψ::ITensorNetwork, bond_tensors::DataGraph, mts::DataGraph) @@ -85,7 +86,7 @@ function benchmark_state_gauging( C[i] = vidal_itn_canonicalness(ψ, bond_tensors) end - + @show times_iters, time simulation_times = cumsum(times_iters) + times_gauging return simulation_times, C diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index bbd73dd6..bb7f376f 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -22,12 +22,13 @@ function message_tensors_skeleton(subgraphs::DataGraph) end function message_tensors( - subgraphs::DataGraph; itensor_constructor=inds_e -> dense(delta(inds_e)) + subgraphs::DataGraph; itensor_constructor=x -> ITensor[dense(delta(i)) for i in x] ) mts = message_tensors_skeleton(subgraphs) for e in edges(subgraphs) - inds_e = commoninds(subgraphs[src(e)], subgraphs[dst(e)]) - mts[e] = ITensorNetwork(map(itensor_constructor, inds_e)) + inds_e = [i for i in commoninds(subgraphs[src(e)], subgraphs[dst(e)])] + itensors = itensor_constructor(inds_e) + mts[e] = ITensorNetwork(itensors) mts[reverse(e)] = dag(mts[e]) end return mts @@ -118,7 +119,7 @@ function belief_propagation( contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), niters=20, target_precision::Union{Float64,Nothing}=nothing, - update_sequence::String="parallel", + update_sequence::String="sequential", edges = edge_update_order(undirected_graph(underlying_graph(mts))) ) compute_norm = target_precision == nothing ? false : true @@ -143,10 +144,11 @@ function belief_propagation( npartitions=nothing, subgraph_vertices=nothing, niters=20, + update_sequence::String="sequential", target_precision::Union{Float64,Nothing}=nothing, ) mts = message_tensors(tn; nvertices_per_partition, npartitions, subgraph_vertices) - return belief_propagation(tn, mts; contract_kwargs, niters, target_precision) + return belief_propagation(tn, mts; contract_kwargs, niters, target_precision, update_sequence) end """ diff --git a/src/beliefpropagation/sqrt_beliefpropagation.jl b/src/beliefpropagation/sqrt_beliefpropagation.jl index 07b074da..c9ba053e 100644 --- a/src/beliefpropagation/sqrt_beliefpropagation.jl +++ b/src/beliefpropagation/sqrt_beliefpropagation.jl @@ -5,7 +5,7 @@ function sqrt_belief_propagation( tn::ITensorNetwork, mts::DataGraph; niters=20, - update_sequence::String="parallel", + update_sequence::String="sequential", # target_precision::Union{Float64,Nothing}=nothing, ) # compute_norm = target_precision == nothing ? false : true @@ -27,8 +27,8 @@ end function sqrt_belief_propagation_iteration( tn::ITensorNetwork, sqrt_mts::DataGraph; - update_sequence::String="parallel", - edges=Graphs.edges(sqrt_mts), + update_sequence::String="sequential", + edges=edge_update_order(undirected_graph(underlying_graph(mts))), # compute_norm=false, ) diff --git a/src/gauging.jl b/src/gauging.jl index 2ef8dd0c..79aca813 100644 --- a/src/gauging.jl +++ b/src/gauging.jl @@ -95,6 +95,7 @@ function vidal_gauge( regularization=10 * eps(real(scalartype(ψ))), niters=30, target_canonicalness::Union{Nothing,Float64}=nothing, + update_sequence = "sequential", svd_kwargs..., ) ψψ = norm_network(ψ) @@ -102,7 +103,7 @@ function vidal_gauge( mts = message_tensors(Z) mts = belief_propagation( - ψψ, mts; contract_kwargs=(; alg="exact"), niters, target_precision=target_canonicalness + ψψ, mts; contract_kwargs=(; alg="exact"), niters, target_precision=target_canonicalness, update_sequence ) return vidal_gauge( ψ, mts; eigen_message_tensor_cutoff, regularization, niters, svd_kwargs... diff --git a/test/test_belief_propagation.jl b/test/test_belief_propagation.jl index 4c75bf90..dff03220 100644 --- a/test/test_belief_propagation.jl +++ b/test/test_belief_propagation.jl @@ -52,7 +52,7 @@ ITensors.disable_warn_order() ψψ, mts, [(v, 1)]; verts_tn=ITensorNetwork(ITensor[apply(op("Sz", s[v]), ψ[v])]) ) denominator_network = approx_network_region(ψψ, mts, [(v, 1)]) - bp_sz = contract(numerator_network)[] / contract(denominator_network)[] + bp_sz = ITensors.contract(numerator_network)[] / ITensors.contract(denominator_network)[] @test abs.(bp_sz - exact_sz) <= 1e-14 @@ -85,7 +85,7 @@ ITensors.disable_warn_order() ψψ, mts, [(v, 1)]; verts_tn=ITensorNetwork(ITensor[apply(op("Sz", s[v]), ψ[v])]) ) denominator_network = approx_network_region(ψψ, mts, [(v, 1)]) - bp_sz = contract(numerator_network)[] / contract(denominator_network)[] + bp_sz = ITensors.contract(numerator_network)[] / ITensors.contract(denominator_network)[] @test abs.(bp_sz - exact_sz) <= 1e-14 @@ -112,7 +112,7 @@ ITensors.disable_warn_order() ) denominator_network = approx_network_region(ψψ, mts, vs) - bp_szsz = contract(numerator_network)[] / contract(denominator_network)[] + bp_szsz = ITensors.contract(numerator_network)[] / ITensors.contract(denominator_network)[] @test abs.(bp_szsz - actual_szsz) <= 0.05 @@ -134,7 +134,7 @@ ITensors.disable_warn_order() mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact")) ψψsplit = split_index(ψψ, NamedEdge.([(v, 1) => (v, 2) for v in vs])) - rdm = contract( + rdm = ITensors.contract( approx_network_region( ψψ, mts, @@ -183,7 +183,7 @@ ITensors.disable_warn_order() numerator_network = approx_network_region(ψψ, mts, [v]; verts_tn=ITensorNetwork(ψOψ[v])) denominator_network = approx_network_region(ψψ, mts, [v]) - bp_sz = contract(numerator_network)[] / contract(denominator_network)[] + bp_sz = ITensors.contract(numerator_network)[] / ITensors.contract(denominator_network)[] exact_sz = contract_boundary_mps(ψOψ; cutoff=1e-16) / contract_boundary_mps(ψψ; cutoff=1e-16) diff --git a/test/test_gauging.jl b/test/test_gauging.jl index b8d9fb5f..5d5e6c77 100644 --- a/test/test_gauging.jl +++ b/test/test_gauging.jl @@ -21,7 +21,7 @@ using SplitApplyCombine dims = (n, n) g = named_grid(dims) s = siteinds("S=1/2", g) - χ = 3 + χ = 6 Random.seed!(5467) ψ = randomITensorNetwork(s; link_space=χ) From 1e6ec2b0f3d0d19ed68687ba3831ee8438dd78ef Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Thu, 26 Oct 2023 16:11:01 -0400 Subject: [PATCH 14/37] Better specification of update sequence for BP --- examples/belief_propagation/bpsequences.jl | 79 +++++++++--- src/beliefpropagation/beliefpropagation.jl | 114 ++++++++++++++---- .../sqrt_beliefpropagation.jl | 94 +++++++++------ src/gauging.jl | 14 ++- src/utils.jl | 6 +- test/test_belief_propagation.jl | 26 ++-- 6 files changed, 227 insertions(+), 106 deletions(-) diff --git a/examples/belief_propagation/bpsequences.jl b/examples/belief_propagation/bpsequences.jl index 6a9797e0..5691f6ae 100644 --- a/examples/belief_propagation/bpsequences.jl +++ b/examples/belief_propagation/bpsequences.jl @@ -14,8 +14,7 @@ using ITensorNetworks: nested_graph_leaf_vertices function main() - - g = named_comb_tree((6,6)) + g = named_comb_tree((6, 6)) s = siteinds("S=1/2", g) χ = 4 @@ -29,17 +28,33 @@ function main() ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ)))) ) - println("First testing out a comb tree. Random network with bond dim $χ") + println("\nFirst testing out a comb tree. Random network with bond dim $χ") #Now test out various sequences print("Parallel updates (sequence is irrelevant): ") - belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "parallel") + belief_propagation( + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + edges=[[e] for e in edges(mts_init)], + ) print("Sequential updates (sequence is default edge list of the message tensors): ") - belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "sequential", edges = edges(mts_init)) + belief_propagation( + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + edges=[e for e in edges(mts_init)], + ) print("Sequential updates (sequence is our custom sequence finder): ") - belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "sequential") + belief_propagation( + ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision=1e-10, niters=100 + ) - g = named_grid((6,6)) + g = named_grid((6, 6)) s = siteinds("S=1/2", g) χ = 2 @@ -53,17 +68,33 @@ function main() ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ)))) ) - println("Now testing out a 6x6 grid. Random network with bond dim $χ") + println("\nNow testing out a 6x6 grid. Random network with bond dim $χ") #Now test out various sequences print("Parallel updates (sequence is irrelevant): ") - belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "parallel") + belief_propagation( + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + edges=[[e] for e in edges(mts_init)], + ) print("Sequential updates (sequence is default edge list of the message tensors): ") - belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "sequential", edges = edges(mts_init)) + belief_propagation( + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + edges=[e for e in edges(mts_init)], + ) print("Sequential updates (sequence is our custom sequence finder): ") - belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "sequential") + belief_propagation( + ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision=1e-10, niters=100 + ) - g = NamedGraphs.hexagonal_lattice_graph(4,4) + g = NamedGraphs.hexagonal_lattice_graph(4, 4) s = siteinds("S=1/2", g) χ = 3 @@ -77,15 +108,31 @@ function main() ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ)))) ) - println("Now testing out a 4 x 4 hexagonal lattice. Random network with bond dim $χ") + println("\nNow testing out a 4 x 4 hexagonal lattice. Random network with bond dim $χ") #Now test out various sequences print("Parallel updates (sequence is irrelevant): ") - belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "parallel") + belief_propagation( + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + edges=[[e] for e in edges(mts_init)], + ) print("Sequential updates (sequence is default edge list of the message tensors): ") - belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "sequential", edges = edges(mts_init)) + belief_propagation( + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + edges=[e for e in edges(mts_init)], + ) print("Sequential updates (sequence is our custom sequence finder): ") - belief_propagation(ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision = 1e-10, niters = 100, update_sequence = "sequential") + return belief_propagation( + ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision=1e-10, niters=100 + ) end main() diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index bb7f376f..b4558993 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -75,31 +75,24 @@ function update_message_tensor( end """ -Do an update of all message tensors for a given ITensornetwork and its partition into sub graphs +Do a sequential update of message tensors on `edges` for a given ITensornetwork and its partition into sub graphs """ function belief_propagation_iteration( tn::ITensorNetwork, - mts::DataGraph; + mts::DataGraph, + edges::Vector{E}; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), compute_norm=false, - update_sequence::String="sequential", - edges = edge_update_order(undirected_graph(underlying_graph(mts))), -) +) where {E<:NamedEdge} new_mts = copy(mts) - if update_sequence != "parallel" && update_sequence != "sequential" - error( - "Specified update order is not currently implemented. Choose parallel or sequential." - ) - end - incoming_mts = update_sequence == "parallel" ? mts : new_mts c = 0 for e in edges environment_tensornetworks = ITensorNetwork[ - incoming_mts[e_in] for - e_in in setdiff(boundary_edges(incoming_mts, [src(e)]; dir=:in), [reverse(e)]) + new_mts[e_in] for + e_in in setdiff(boundary_edges(new_mts, [src(e)]; dir=:in), [reverse(e)]) ] new_mts[src(e) => dst(e)] = update_message_tensor( - tn, incoming_mts[src(e)], environment_tensornetworks; contract_kwargs + tn, new_mts[src(e)], environment_tensornetworks; contract_kwargs ) if compute_norm @@ -113,24 +106,96 @@ function belief_propagation_iteration( return new_mts, c / (length(edges)) end +""" +Do parallel updates between groups of edges of all message tensors for a given ITensornetwork and its partition into sub graphs +""" +function belief_propagation_iteration( + tn::ITensorNetwork, + mts::DataGraph, + edge_groups::Vector{Vector{E}}; + contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), + compute_norm=false, +) where {E<:NamedEdge} + new_mts = copy(mts) + c = 0 + for edges in edge_groups + updated_mts, ct = belief_propagation_iteration( + tn, mts, edges; contract_kwargs, compute_norm + ) + for e in edges + new_mts[e] = updated_mts[e] + end + c += ct + end + return new_mts, c / (length(edge_groups)) +end + +function belief_propagation_iteration( + tn::ITensorNetwork, + mts::DataGraph; + contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), + compute_norm=false, + edges::Union{Vector{Vector{E}},Vector{E}}=edge_update_order( + undirected_graph(underlying_graph(mts)) + ), +) where {E<:NamedEdge} + return belief_propagation_iteration(tn, mts, edges; contract_kwargs, compute_norm) +end + +# """ +# Do an update of all message tensors for a given ITensornetwork and its partition into sub graphs +# """ +# function belief_propagation_iteration( +# tn::ITensorNetwork, +# mts::DataGraph; +# contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), +# compute_norm=false, +# update_sequence::String="sequential", +# edges::Vector{Vector{}} = edge_update_order(undirected_graph(underlying_graph(mts))), +# ) +# new_mts = copy(mts) +# if update_sequence != "parallel" && update_sequence != "sequential" +# error( +# "Specified update order is not currently implemented. Choose parallel or sequential." +# ) +# end +# incoming_mts = update_sequence == "parallel" ? mts : new_mts +# c = 0 +# for e in edges +# environment_tensornetworks = ITensorNetwork[ +# incoming_mts[e_in] for +# e_in in setdiff(boundary_edges(incoming_mts, [src(e)]; dir=:in), [reverse(e)]) +# ] +# new_mts[src(e) => dst(e)] = update_message_tensor( +# tn, incoming_mts[src(e)], environment_tensornetworks; contract_kwargs +# ) + +# if compute_norm +# LHS, RHS = ITensors.contract(ITensor(mts[src(e) => dst(e)])), +# ITensors.contract(ITensor(new_mts[src(e) => dst(e)])) +# LHS /= sum(diag(LHS)) +# RHS /= sum(diag(RHS)) +# c += 0.5 * norm(denseblocks(LHS) - denseblocks(RHS)) +# end +# end +# return new_mts, c / (length(edges)) +# end + function belief_propagation( tn::ITensorNetwork, mts::DataGraph; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), niters=20, target_precision::Union{Float64,Nothing}=nothing, - update_sequence::String="sequential", - edges = edge_update_order(undirected_graph(underlying_graph(mts))) -) + edges::Union{Vector{Vector{E}},Vector{E}}=edge_update_order( + undirected_graph(underlying_graph(mts)) + ), +) where {E<:NamedEdge} compute_norm = target_precision == nothing ? false : true for i in 1:niters - mts, c = belief_propagation_iteration( - tn, mts; contract_kwargs, compute_norm, update_sequence, edges - ) + mts, c = belief_propagation_iteration(tn, mts, edges; contract_kwargs, compute_norm) if compute_norm && c <= target_precision - println( - "BP converged to desired precision after $i iterations.", - ) + println("BP converged to desired precision after $i iterations.") break end end @@ -144,11 +209,10 @@ function belief_propagation( npartitions=nothing, subgraph_vertices=nothing, niters=20, - update_sequence::String="sequential", target_precision::Union{Float64,Nothing}=nothing, ) mts = message_tensors(tn; nvertices_per_partition, npartitions, subgraph_vertices) - return belief_propagation(tn, mts; contract_kwargs, niters, target_precision, update_sequence) + return belief_propagation(tn, mts; contract_kwargs, niters, target_precision) end """ diff --git a/src/beliefpropagation/sqrt_beliefpropagation.jl b/src/beliefpropagation/sqrt_beliefpropagation.jl index c9ba053e..c4eff09e 100644 --- a/src/beliefpropagation/sqrt_beliefpropagation.jl +++ b/src/beliefpropagation/sqrt_beliefpropagation.jl @@ -1,53 +1,19 @@ # using ITensors: scalartype # using ITensorNetworks: find_subgraph, map_diag, sqrt_diag, boundary_edges -function sqrt_belief_propagation( - tn::ITensorNetwork, - mts::DataGraph; - niters=20, - update_sequence::String="sequential", - # target_precision::Union{Float64,Nothing}=nothing, -) - # compute_norm = target_precision == nothing ? false : true - sqrt_mts = sqrt_message_tensors(tn, mts) - for i in 1:niters - sqrt_mts, c = sqrt_belief_propagation_iteration(tn, sqrt_mts; update_sequence) #; compute_norm) - # if compute_norm && c <= target_precision - # println( - # "Belief Propagation finished. Reached a canonicalness of " * - # string(c) * - # " after $i iterations. ", - # ) - # break - # end - end - return sqr_message_tensors(sqrt_mts) -end - function sqrt_belief_propagation_iteration( - tn::ITensorNetwork, - sqrt_mts::DataGraph; - update_sequence::String="sequential", - edges=edge_update_order(undirected_graph(underlying_graph(mts))), - - # compute_norm=false, -) + tn::ITensorNetwork, sqrt_mts::DataGraph, edges::Vector{E} +) where {E<:NamedEdge} new_sqrt_mts = copy(sqrt_mts) - if update_sequence != "parallel" && update_sequence != "sequential" - error( - "Specified update order is not currently implemented. Choose parallel or sequential." - ) - end - incoming_sqrt_mts = update_sequence == "parallel" ? sqrt_mts : new_sqrt_mts c = 0.0 for e in edges environment_tensornetworks = ITensorNetwork[ - incoming_sqrt_mts[e_in] for - e_in in setdiff(boundary_edges(incoming_sqrt_mts, [src(e)]; dir=:in), [reverse(e)]) + new_sqrt_mts[e_in] for + e_in in setdiff(boundary_edges(new_sqrt_mts, [src(e)]; dir=:in), [reverse(e)]) ] new_sqrt_mts[src(e) => dst(e)] = update_sqrt_message_tensor( - tn, incoming_sqrt_mts[src(e)], environment_tensornetworks; + tn, new_sqrt_mts[src(e)], environment_tensornetworks; ) # if compute_norm @@ -61,6 +27,56 @@ function sqrt_belief_propagation_iteration( return new_sqrt_mts, c / (length(edges)) end +function sqrt_belief_propagation_iteration( + tn::ITensorNetwork, sqrt_mts::DataGraph, edges::Vector{Vector{E}} +) where {E<:NamedEdge} + new_sqrt_mts = copy(sqrt_mts) + c = 0.0 + for e_group in edges + updated_sqrt_mts, ct = sqrt_belief_propagation_iteration(tn, sqr_mts, e_group) + for e in e_group + new_sqrt_mts[e] = updated_sqrt_mts[e] + end + c += ct + end + return new_sqrt_mts, c / (length(edges)) +end + +function sqrt_belief_propagation_iteration( + tn::ITensorNetwork, + sqrt_mts::DataGraph; + edges::Union{Vector{Vector{E}},Vector{E}}=edge_update_order( + undirected_graph(underlying_graph(mts)) + ), +) where {E<:NamedEdge} + return sqrt_belief_propagation_iteration(tn, sqrt_mts, edges) +end + +function sqrt_belief_propagation( + tn::ITensorNetwork, + mts::DataGraph; + niters=20, + edges::Union{Vector{Vector{E}},Vector{E}}=edge_update_order( + undirected_graph(underlying_graph(mts)) + ), + # target_precision::Union{Float64,Nothing}=nothing, +) where {E<:NamedEdge} + # compute_norm = target_precision == nothing ? false : true + sqrt_mts = sqrt_message_tensors(tn, mts) + for i in 1:niters + sqrt_mts, c = sqrt_belief_propagation_iteration(tn, sqrt_mts, edges) #; compute_norm) + # if compute_norm && c <= target_precision + # println( + # "Belief Propagation finished. Reached a canonicalness of " * + # string(c) * + # " after $i iterations. ", + # ) + # break + # end + end + return sqr_message_tensors(sqrt_mts) +end + function update_sqrt_message_tensor( tn::ITensorNetwork, subgraph_vertices::Vector, sqrt_mts::Vector{ITensorNetwork}; ) diff --git a/src/gauging.jl b/src/gauging.jl index 79aca813..4086198c 100644 --- a/src/gauging.jl +++ b/src/gauging.jl @@ -17,7 +17,7 @@ function vidal_gauge( bond_tensors::DataGraph; eigen_message_tensor_cutoff=10 * eps(real(scalartype(ψ))), regularization=10 * eps(real(scalartype(ψ))), - edges = NamedGraphs.edges(ψ), + edges=NamedGraphs.edges(ψ), svd_kwargs..., ) ψ_vidal = copy(ψ) @@ -80,11 +80,12 @@ function vidal_gauge( mts::DataGraph; eigen_message_tensor_cutoff=10 * eps(real(scalartype(ψ))), regularization=10 * eps(real(scalartype(ψ))), + edges=NamedGraphs.edges(ψ), svd_kwargs..., ) bond_tensors = initialize_bond_tensors(ψ) return vidal_gauge( - ψ, mts, bond_tensors; eigen_message_tensor_cutoff, regularization, svd_kwargs... + ψ, mts, bond_tensors; eigen_message_tensor_cutoff, regularization, edges, svd_kwargs... ) end @@ -95,7 +96,6 @@ function vidal_gauge( regularization=10 * eps(real(scalartype(ψ))), niters=30, target_canonicalness::Union{Nothing,Float64}=nothing, - update_sequence = "sequential", svd_kwargs..., ) ψψ = norm_network(ψ) @@ -103,7 +103,7 @@ function vidal_gauge( mts = message_tensors(Z) mts = belief_propagation( - ψψ, mts; contract_kwargs=(; alg="exact"), niters, target_precision=target_canonicalness, update_sequence + ψψ, mts; contract_kwargs=(; alg="exact"), niters, target_precision=target_canonicalness ) return vidal_gauge( ψ, mts; eigen_message_tensor_cutoff, regularization, niters, svd_kwargs... @@ -175,7 +175,11 @@ function symmetric_to_vidal_gauge( end """Function to measure the 'isometries' of a state in the Vidal Gauge""" -function vidal_itn_isometries(ψ::ITensorNetwork, bond_tensors::DataGraph; edges = vcat(NamedGraphs.edges(ψ), reverse.(NamedGraphs.edges(ψ)))) +function vidal_itn_isometries( + ψ::ITensorNetwork, + bond_tensors::DataGraph; + edges=vcat(NamedGraphs.edges(ψ), reverse.(NamedGraphs.edges(ψ))), +) isometries = DataGraph{vertextype(ψ),ITensor,ITensor}(directed_graph(underlying_graph(ψ))) for e in edges diff --git a/src/utils.jl b/src/utils.jl index 1d821dc8..e3e2bc82 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -40,8 +40,10 @@ function edge_update_order(g) end #Find an optimal ordering of the edges in a tree -function tree_edge_update_order(g::AbstractNamedGraph; root_vertex = first(keys(sort(eccentricities(g); rev=true)))) +function tree_edge_update_order( + g::AbstractNamedGraph; root_vertex=first(keys(sort(eccentricities(g); rev=true))) +) @assert is_tree(g) es = post_order_dfs_edges(g, root_vertex) return vcat(es, reverse(reverse.(es))) -end \ No newline at end of file +end diff --git a/test/test_belief_propagation.jl b/test/test_belief_propagation.jl index dff03220..82a2db66 100644 --- a/test/test_belief_propagation.jl +++ b/test/test_belief_propagation.jl @@ -40,13 +40,7 @@ ITensors.disable_warn_order() Z = partition(ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ))))) mts = message_tensors(Z) - mts = belief_propagation( - ψψ, - mts; - contract_kwargs=(; alg="exact"), - niters = 1, - update_sequence="sequential", - ) + mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact"), niters=1) numerator_network = approx_network_region( ψψ, mts, [(v, 1)]; verts_tn=ITensorNetwork(ITensor[apply(op("Sz", s[v]), ψ[v])]) @@ -57,7 +51,7 @@ ITensors.disable_warn_order() @test abs.(bp_sz - exact_sz) <= 1e-14 #Now test on a tree, should also be exact - g = named_comb_tree((6,6)) + g = named_comb_tree((6, 6)) s = siteinds("S=1/2", g) χ = 2 Random.seed!(1564) @@ -73,13 +67,7 @@ ITensors.disable_warn_order() Z = partition(ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ))))) mts = message_tensors(Z) - mts = belief_propagation( - ψψ, - mts; - contract_kwargs=(; alg="exact"), - niters = 1, - update_sequence="sequential", - ) + mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact"), niters=1) numerator_network = approx_network_region( ψψ, mts, [(v, 1)]; verts_tn=ITensorNetwork(ITensor[apply(op("Sz", s[v]), ψ[v])]) @@ -112,7 +100,8 @@ ITensors.disable_warn_order() ) denominator_network = approx_network_region(ψψ, mts, vs) - bp_szsz = ITensors.contract(numerator_network)[] / ITensors.contract(denominator_network)[] + bp_szsz = + ITensors.contract(numerator_network)[] / ITensors.contract(denominator_network)[] @test abs.(bp_szsz - actual_szsz) <= 0.05 @@ -176,8 +165,7 @@ ITensors.disable_warn_order() contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, cutoff=1e-16, maxdim ), - niters = 1, - update_sequence="sequential", + niters=1, ) numerator_network = approx_network_region(ψψ, mts, [v]; verts_tn=ITensorNetwork(ψOψ[v])) @@ -189,4 +177,4 @@ ITensors.disable_warn_order() contract_boundary_mps(ψOψ; cutoff=1e-16) / contract_boundary_mps(ψψ; cutoff=1e-16) @test abs.(bp_sz - exact_sz) <= 1e-5 -end \ No newline at end of file +end From 4c352ab83bd525426eed239b2e5ee5e49cdb9946 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Thu, 2 Nov 2023 11:16:52 -0400 Subject: [PATCH 15/37] Fixed IBM processor construction to reflect row and column name swapping in NamedGraphs.hexagonal_lattice_graph --- examples/dynamics/heavy_hex_ising_real_tebd.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/examples/dynamics/heavy_hex_ising_real_tebd.jl b/examples/dynamics/heavy_hex_ising_real_tebd.jl index 4ce55b20..e8dcc863 100644 --- a/examples/dynamics/heavy_hex_ising_real_tebd.jl +++ b/examples/dynamics/heavy_hex_ising_real_tebd.jl @@ -24,10 +24,10 @@ end function ibm_processor_graph(n::Int64, m::Int64) g = heavy_hex_lattice_graph(n, m) dims = maximum(vertices(hexagonal_lattice_graph(n, m))) - v1, v2 = (1, dims[2]), (dims[1], 1) + v1, v2 = (dims[1], 1), (1, dims[2]) add_vertices!(g, [v1, v2]) - add_edge!(g, v1 => v1 .- (0, 1)) - add_edge!(g, v2 => v2 .+ (0, 1)) + add_edge!(g, v1 => v1 .- (1, 0)) + add_edge!(g, v2 => v2 .+ (1, 0)) return g end From 3ad3cbc8710eba03a8c265989718f9383dbf96a2 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Thu, 2 Nov 2023 12:26:14 -0400 Subject: [PATCH 16/37] Forest cover for specifying edge update order. Better specification of parallel vs sequential via edge kwarg. Further examples in BPSequences --- examples/belief_propagation/bpsequences.jl | 75 ++++++++++++++++++- examples/gauging/gauging_itns.jl | 22 ++++-- src/beliefpropagation/beliefpropagation.jl | 72 +++++++----------- .../sqrt_beliefpropagation.jl | 4 +- src/gauging.jl | 8 +- src/utils.jl | 18 ++--- 6 files changed, 131 insertions(+), 68 deletions(-) diff --git a/examples/belief_propagation/bpsequences.jl b/examples/belief_propagation/bpsequences.jl index 5691f6ae..fbcf0c14 100644 --- a/examples/belief_propagation/bpsequences.jl +++ b/examples/belief_propagation/bpsequences.jl @@ -4,6 +4,7 @@ using Metis using ITensorNetworks using Random using SplitApplyCombine +using Graphs using NamedGraphs using ITensorNetworks: @@ -39,6 +40,7 @@ function main() target_precision=1e-10, niters=100, edges=[[e] for e in edges(mts_init)], + verbose=true, ) print("Sequential updates (sequence is default edge list of the message tensors): ") belief_propagation( @@ -48,10 +50,63 @@ function main() target_precision=1e-10, niters=100, edges=[e for e in edges(mts_init)], + verbose=true, ) print("Sequential updates (sequence is our custom sequence finder): ") belief_propagation( - ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision=1e-10, niters=100 + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + verbose=true, + ) + + g = NamedGraph(Graphs.random_regular_graph(100, 3)) + s = siteinds("S=1/2", g) + χ = 4 + + Random.seed!(5467) + + ψ = randomITensorNetwork(s; link_space=χ) + ψψ = ψ ⊗ prime(dag(ψ); sites=[]) + + #Initial message tensors for BP + mts_init = message_tensors( + ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ)))) + ) + + println("\nNow testing out a z = 3 random regular graph. Random network with bond dim $χ") + + #Now test out various sequences + print("Parallel updates (sequence is irrelevant): ") + belief_propagation( + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + edges=[[e] for e in edges(mts_init)], + verbose=true, + ) + print("Sequential updates (sequence is default edge list of the message tensors): ") + belief_propagation( + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + edges=[e for e in edges(mts_init)], + verbose=true, + ) + print("Sequential updates (sequence is our custom sequence finder): ") + belief_propagation( + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + verbose=true, ) g = named_grid((6, 6)) @@ -79,6 +134,7 @@ function main() target_precision=1e-10, niters=100, edges=[[e] for e in edges(mts_init)], + verbose=true, ) print("Sequential updates (sequence is default edge list of the message tensors): ") belief_propagation( @@ -88,10 +144,16 @@ function main() target_precision=1e-10, niters=100, edges=[e for e in edges(mts_init)], + verbose=true, ) print("Sequential updates (sequence is our custom sequence finder): ") belief_propagation( - ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision=1e-10, niters=100 + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + verbose=true, ) g = NamedGraphs.hexagonal_lattice_graph(4, 4) @@ -119,6 +181,7 @@ function main() target_precision=1e-10, niters=100, edges=[[e] for e in edges(mts_init)], + verbose=true, ) print("Sequential updates (sequence is default edge list of the message tensors): ") belief_propagation( @@ -128,10 +191,16 @@ function main() target_precision=1e-10, niters=100, edges=[e for e in edges(mts_init)], + verbose=true, ) print("Sequential updates (sequence is our custom sequence finder): ") return belief_propagation( - ψψ, mts_init; contract_kwargs=(; alg="exact"), target_precision=1e-10, niters=100 + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + verbose=true, ) end diff --git a/examples/gauging/gauging_itns.jl b/examples/gauging/gauging_itns.jl index 0f65dc8e..b420d45e 100644 --- a/examples/gauging/gauging_itns.jl +++ b/examples/gauging/gauging_itns.jl @@ -48,7 +48,7 @@ function benchmark_state_gauging( ψ::ITensorNetwork; mode="BeliefPropagation", no_iterations=50, - BP_update_order::String="parallel", + BP_update_order::String="sequential", ) s = siteinds(ψ) @@ -69,9 +69,15 @@ function benchmark_state_gauging( println("On Iteration " * string(i)) if mode == "BeliefPropagation" - times_iters[i] = @elapsed mts, _ = belief_propagation_iteration( - ψψ, mts; contract_kwargs=(; alg="exact"), update_sequence=BP_update_order - ) + if BP_update_order != "parallel" + times_iters[i] = @elapsed mts, _ = belief_propagation_iteration( + ψψ, mts; contract_kwargs=(; alg="exact") + ) + else + times_iters[i] = @elapsed mts, _ = belief_propagation_iteration( + ψψ, mts; contract_kwargs=(; alg="exact"), edges=[[e] for e in edges(mts)] + ) + end times_gauging[i] = @elapsed ψ, bond_tensors = vidal_gauge(ψinit, mts) elseif mode == "Eager" @@ -98,14 +104,16 @@ s = siteinds("S=1/2", g) ψ = randomITensorNetwork(s; link_space=χ) no_iterations = 30 -BPG_simulation_times, BPG_Cs = benchmark_state_gauging(ψ; no_iterations) +BPG_simulation_times, BPG_Cs = benchmark_state_gauging( + ψ; no_iterations, BP_update_order="parallel" +) BPG_sequential_simulation_times, BPG_sequential_Cs = benchmark_state_gauging( - ψ; no_iterations, BP_update_order="sequential" + ψ; no_iterations ) Eager_simulation_times, Eager_Cs = benchmark_state_gauging(ψ; mode="Eager", no_iterations) SU_simulation_times, SU_Cs = benchmark_state_gauging(ψ; mode="SU", no_iterations) -epsilon = 1e-6 +epsilon = 1e-10 println( "Time for BPG (with parallel updates) to reach C < epsilon was " * diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index b4558993..64581e55 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -135,67 +135,31 @@ function belief_propagation_iteration( mts::DataGraph; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), compute_norm=false, - edges::Union{Vector{Vector{E}},Vector{E}}=edge_update_order( + edges::Union{Vector{Vector{E}},Vector{E}}=belief_propagation_edge_sequence( undirected_graph(underlying_graph(mts)) ), ) where {E<:NamedEdge} return belief_propagation_iteration(tn, mts, edges; contract_kwargs, compute_norm) end -# """ -# Do an update of all message tensors for a given ITensornetwork and its partition into sub graphs -# """ -# function belief_propagation_iteration( -# tn::ITensorNetwork, -# mts::DataGraph; -# contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), -# compute_norm=false, -# update_sequence::String="sequential", -# edges::Vector{Vector{}} = edge_update_order(undirected_graph(underlying_graph(mts))), -# ) -# new_mts = copy(mts) -# if update_sequence != "parallel" && update_sequence != "sequential" -# error( -# "Specified update order is not currently implemented. Choose parallel or sequential." -# ) -# end -# incoming_mts = update_sequence == "parallel" ? mts : new_mts -# c = 0 -# for e in edges -# environment_tensornetworks = ITensorNetwork[ -# incoming_mts[e_in] for -# e_in in setdiff(boundary_edges(incoming_mts, [src(e)]; dir=:in), [reverse(e)]) -# ] -# new_mts[src(e) => dst(e)] = update_message_tensor( -# tn, incoming_mts[src(e)], environment_tensornetworks; contract_kwargs -# ) - -# if compute_norm -# LHS, RHS = ITensors.contract(ITensor(mts[src(e) => dst(e)])), -# ITensors.contract(ITensor(new_mts[src(e) => dst(e)])) -# LHS /= sum(diag(LHS)) -# RHS /= sum(diag(RHS)) -# c += 0.5 * norm(denseblocks(LHS) - denseblocks(RHS)) -# end -# end -# return new_mts, c / (length(edges)) -# end - function belief_propagation( tn::ITensorNetwork, mts::DataGraph; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), niters=20, target_precision::Union{Float64,Nothing}=nothing, - edges::Union{Vector{Vector{E}},Vector{E}}=edge_update_order( + edges::Union{Vector{Vector{E}},Vector{E}}=belief_propagation_edge_sequence( undirected_graph(underlying_graph(mts)) ), + verbose=false, ) where {E<:NamedEdge} compute_norm = target_precision == nothing ? false : true for i in 1:niters mts, c = belief_propagation_iteration(tn, mts, edges; contract_kwargs, compute_norm) if compute_norm && c <= target_precision - println("BP converged to desired precision after $i iterations.") + if verbose + println("BP converged to desired precision after $i iterations.") + end break end end @@ -210,9 +174,10 @@ function belief_propagation( subgraph_vertices=nothing, niters=20, target_precision::Union{Float64,Nothing}=nothing, + verbose=false, ) mts = message_tensors(tn; nvertices_per_partition, npartitions, subgraph_vertices) - return belief_propagation(tn, mts; contract_kwargs, niters, target_precision) + return belief_propagation(tn, mts; contract_kwargs, niters, target_precision, verbose) end """ @@ -247,3 +212,24 @@ function approx_network_region( return environment_tn ⊗ verts_tn end + +""" +Return a custom edge order for how how to update all BP message tensors on a general undirected graph. +On a tree this will yield a sequence which only needs to be performed once. Based on forest covers and depth first searches amongst the forests. +""" +function belief_propagation_edge_sequence( + g::NamedGraph; root_vertex=NamedGraphs.default_root_vertex +) + @assert !is_directed(g) + forests = NamedGraphs.forest_cover(g) + edges = NamedEdge[] + for forest in forests + trees = NamedGraph[forest[vs] for vs in connected_components(forest)] + for tree in trees + tree_edges = post_order_dfs_edges(tree, root_vertex(tree)) + push!(edges, vcat(tree_edges, reverse(reverse.(tree_edges)))...) + end + end + + return edges +end diff --git a/src/beliefpropagation/sqrt_beliefpropagation.jl b/src/beliefpropagation/sqrt_beliefpropagation.jl index c4eff09e..f28ab9a0 100644 --- a/src/beliefpropagation/sqrt_beliefpropagation.jl +++ b/src/beliefpropagation/sqrt_beliefpropagation.jl @@ -45,7 +45,7 @@ end function sqrt_belief_propagation_iteration( tn::ITensorNetwork, sqrt_mts::DataGraph; - edges::Union{Vector{Vector{E}},Vector{E}}=edge_update_order( + edges::Union{Vector{Vector{E}},Vector{E}}=belief_propagation_edge_sequence( undirected_graph(underlying_graph(mts)) ), ) where {E<:NamedEdge} @@ -56,7 +56,7 @@ function sqrt_belief_propagation( tn::ITensorNetwork, mts::DataGraph; niters=20, - edges::Union{Vector{Vector{E}},Vector{E}}=edge_update_order( + edges::Union{Vector{Vector{E}},Vector{E}}=belief_propagation_edge_sequence( undirected_graph(underlying_graph(mts)) ), # target_precision::Union{Float64,Nothing}=nothing, diff --git a/src/gauging.jl b/src/gauging.jl index 4086198c..17b8f748 100644 --- a/src/gauging.jl +++ b/src/gauging.jl @@ -96,6 +96,7 @@ function vidal_gauge( regularization=10 * eps(real(scalartype(ψ))), niters=30, target_canonicalness::Union{Nothing,Float64}=nothing, + verbose=false, svd_kwargs..., ) ψψ = norm_network(ψ) @@ -103,7 +104,12 @@ function vidal_gauge( mts = message_tensors(Z) mts = belief_propagation( - ψψ, mts; contract_kwargs=(; alg="exact"), niters, target_precision=target_canonicalness + ψψ, + mts; + contract_kwargs=(; alg="exact"), + niters, + target_precision=target_canonicalness, + verbose, ) return vidal_gauge( ψ, mts; eigen_message_tensor_cutoff, regularization, niters, svd_kwargs... diff --git a/src/utils.jl b/src/utils.jl index e3e2bc82..057822cb 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -26,24 +26,18 @@ function line_to_tree(line::Vector) return [line_to_tree(line[1:(end - 1)]), line[end]] end -function edge_update_order(g) - forests = NamedGraphs.build_forest_cover(g) +#Custom edge order for updating all BP message tensors on a general undirected graph. On a tree this will yield a sequence which only needs to be performed once. +function BP_edge_update_order(g::NamedGraph; root_vertex=NamedGraphs.default_root_vertex) + @assert !is_directed(g) + forests = NamedGraphs.forest_cover(g) edges = NamedEdge[] for forest in forests trees = NamedGraph[forest[vs] for vs in connected_components(forest)] for tree in trees - push!(edges, tree_edge_update_order(tree)...) + tree_edges = post_order_dfs_edges(tree, root_vertex(tree)) + push!(edges, vcat(tree_edges, reverse(reverse.(tree_edges)))...) end end return edges end - -#Find an optimal ordering of the edges in a tree -function tree_edge_update_order( - g::AbstractNamedGraph; root_vertex=first(keys(sort(eccentricities(g); rev=true))) -) - @assert is_tree(g) - es = post_order_dfs_edges(g, root_vertex) - return vcat(es, reverse(reverse.(es))) -end From 6df26bda983f6fa7113899dabbda2acbac9efa62 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Mon, 6 Nov 2023 10:12:38 -0500 Subject: [PATCH 17/37] Improvement --- src/utils.jl | 16 ---------------- 1 file changed, 16 deletions(-) diff --git a/src/utils.jl b/src/utils.jl index 057822cb..9a1be885 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -25,19 +25,3 @@ function line_to_tree(line::Vector) end return [line_to_tree(line[1:(end - 1)]), line[end]] end - -#Custom edge order for updating all BP message tensors on a general undirected graph. On a tree this will yield a sequence which only needs to be performed once. -function BP_edge_update_order(g::NamedGraph; root_vertex=NamedGraphs.default_root_vertex) - @assert !is_directed(g) - forests = NamedGraphs.forest_cover(g) - edges = NamedEdge[] - for forest in forests - trees = NamedGraph[forest[vs] for vs in connected_components(forest)] - for tree in trees - tree_edges = post_order_dfs_edges(tree, root_vertex(tree)) - push!(edges, vcat(tree_edges, reverse(reverse.(tree_edges)))...) - end - end - - return edges -end From 6ba7828c880099e27fec038b9dad14fb83d231f3 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Mon, 6 Nov 2023 14:55:29 -0500 Subject: [PATCH 18/37] Added BP sequences to test examples. Removed Sqrt_BP as already tested --- test/test_examples/test_examples.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/test_examples/test_examples.jl b/test/test_examples/test_examples.jl index d4532451..8afc3be4 100644 --- a/test/test_examples/test_examples.jl +++ b/test/test_examples/test_examples.jl @@ -14,7 +14,7 @@ using Test "peps.jl", "steiner_tree.jl", joinpath("belief_propagation", "bpexample.jl"), - joinpath("belief_propagation", "sqrt_bp.jl"), + joinpath("belief_propagation", "bpsequences.jl"), joinpath("dynamics", "2d_ising_imag_tebd.jl"), joinpath("dynamics", "heavy_hex_ising_real_tebd.jl"), joinpath("treetensornetworks", "comb_tree.jl"), From 67fef77fc13ebd4e9bd0347d9661bf145f088be0 Mon Sep 17 00:00:00 2001 From: Joseph Tindall <51231103+JoeyT1994@users.noreply.github.com> Date: Wed, 8 Nov 2023 15:09:56 -0500 Subject: [PATCH 19/37] Update examples/gauging/gauging_itns.jl --- examples/gauging/gauging_itns.jl | 1 - 1 file changed, 1 deletion(-) diff --git a/examples/gauging/gauging_itns.jl b/examples/gauging/gauging_itns.jl index b420d45e..28fa6f63 100644 --- a/examples/gauging/gauging_itns.jl +++ b/examples/gauging/gauging_itns.jl @@ -22,7 +22,6 @@ using ITensorNetworks: using NamedGraphs using NamedGraphs: add_edges!, rem_vertex!, hexagonal_lattice_graph using Graphs -using MKL """Eager Gauging""" function eager_gauging(ψ::ITensorNetwork, bond_tensors::DataGraph, mts::DataGraph) From 57dd13ff4de050ee7155baf66f508f5c111fa3c8 Mon Sep 17 00:00:00 2001 From: Joseph Tindall <51231103+JoeyT1994@users.noreply.github.com> Date: Wed, 8 Nov 2023 15:14:49 -0500 Subject: [PATCH 20/37] Update src/beliefpropagation/beliefpropagation.jl --- src/beliefpropagation/beliefpropagation.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index 64581e55..79828d25 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -214,7 +214,7 @@ function approx_network_region( end """ -Return a custom edge order for how how to update all BP message tensors on a general undirected graph. +Return a custom edge order for how to update all BP message tensors on a general undirected graph. On a tree this will yield a sequence which only needs to be performed once. Based on forest covers and depth first searches amongst the forests. """ function belief_propagation_edge_sequence( From cedb1d6e9266de88c156b15cb16fa6cfda4d4593 Mon Sep 17 00:00:00 2001 From: Joseph Tindall <51231103+JoeyT1994@users.noreply.github.com> Date: Thu, 9 Nov 2023 11:01:38 -0500 Subject: [PATCH 21/37] Update src/beliefpropagation/beliefpropagation.jl --- src/beliefpropagation/beliefpropagation.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index 79828d25..4e254f16 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -83,7 +83,7 @@ function belief_propagation_iteration( edges::Vector{E}; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), compute_norm=false, -) where {E<:NamedEdge} +) new_mts = copy(mts) c = 0 for e in edges From a8da3dc1cca011e56c691448f8a32dd6795f2142 Mon Sep 17 00:00:00 2001 From: Joseph Tindall <51231103+JoeyT1994@users.noreply.github.com> Date: Thu, 9 Nov 2023 11:01:43 -0500 Subject: [PATCH 22/37] Update src/beliefpropagation/beliefpropagation.jl --- src/beliefpropagation/beliefpropagation.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index 4e254f16..edd80ce9 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -80,7 +80,7 @@ Do a sequential update of message tensors on `edges` for a given ITensornetwork function belief_propagation_iteration( tn::ITensorNetwork, mts::DataGraph, - edges::Vector{E}; + edges::Vector{<:NamedEdge}; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), compute_norm=false, ) From a5008ee5816a5b3e629524418519a8eb75c7a03e Mon Sep 17 00:00:00 2001 From: Joseph Tindall <51231103+JoeyT1994@users.noreply.github.com> Date: Thu, 9 Nov 2023 11:01:49 -0500 Subject: [PATCH 23/37] Update src/beliefpropagation/beliefpropagation.jl --- src/beliefpropagation/beliefpropagation.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index edd80ce9..5be6cfaa 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -26,7 +26,7 @@ function message_tensors( ) mts = message_tensors_skeleton(subgraphs) for e in edges(subgraphs) - inds_e = [i for i in commoninds(subgraphs[src(e)], subgraphs[dst(e)])] + inds_e = commoninds(subgraphs[src(e)], subgraphs[dst(e)]) itensors = itensor_constructor(inds_e) mts[e] = ITensorNetwork(itensors) mts[reverse(e)] = dag(mts[e]) From db9a474b7a8a2b2b81d1c4b0c4cd46d89906e707 Mon Sep 17 00:00:00 2001 From: Joseph Tindall <51231103+JoeyT1994@users.noreply.github.com> Date: Thu, 9 Nov 2023 11:01:56 -0500 Subject: [PATCH 24/37] Update src/beliefpropagation/beliefpropagation.jl --- src/beliefpropagation/beliefpropagation.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index 5be6cfaa..37231096 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -153,7 +153,7 @@ function belief_propagation( ), verbose=false, ) where {E<:NamedEdge} - compute_norm = target_precision == nothing ? false : true + compute_norm = !isnothing(target_precision) for i in 1:niters mts, c = belief_propagation_iteration(tn, mts, edges; contract_kwargs, compute_norm) if compute_norm && c <= target_precision From f183ec778f377ccb978413ead044d26518fa0f3a Mon Sep 17 00:00:00 2001 From: Joseph Tindall <51231103+JoeyT1994@users.noreply.github.com> Date: Thu, 9 Nov 2023 11:02:03 -0500 Subject: [PATCH 25/37] Update src/beliefpropagation/beliefpropagation.jl --- src/beliefpropagation/beliefpropagation.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index 37231096..1c5bb98a 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -135,7 +135,7 @@ function belief_propagation_iteration( mts::DataGraph; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), compute_norm=false, - edges::Union{Vector{Vector{E}},Vector{E}}=belief_propagation_edge_sequence( + edges::Union{Vector{Vector{<:NamedEdge}},Vector{<:NamedEdge}}=belief_propagation_edge_sequence( undirected_graph(underlying_graph(mts)) ), ) where {E<:NamedEdge} From 247f2af8433484d4b9bf5352b7f672de8129cc73 Mon Sep 17 00:00:00 2001 From: Joseph Tindall <51231103+JoeyT1994@users.noreply.github.com> Date: Thu, 9 Nov 2023 11:02:09 -0500 Subject: [PATCH 26/37] Update src/beliefpropagation/beliefpropagation.jl --- src/beliefpropagation/beliefpropagation.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index 1c5bb98a..d1c4aead 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -22,7 +22,7 @@ function message_tensors_skeleton(subgraphs::DataGraph) end function message_tensors( - subgraphs::DataGraph; itensor_constructor=x -> ITensor[dense(delta(i)) for i in x] + subgraphs::DataGraph; itensor_constructor=inds_e -> ITensor[dense(delta(i)) for i in inds_e] ) mts = message_tensors_skeleton(subgraphs) for e in edges(subgraphs) From a196c877d9f0d42596798812967c505d23d75423 Mon Sep 17 00:00:00 2001 From: Joseph Tindall <51231103+JoeyT1994@users.noreply.github.com> Date: Thu, 9 Nov 2023 11:02:15 -0500 Subject: [PATCH 27/37] Update src/beliefpropagation/beliefpropagation.jl --- src/beliefpropagation/beliefpropagation.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index d1c4aead..e168e5e2 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -138,7 +138,7 @@ function belief_propagation_iteration( edges::Union{Vector{Vector{<:NamedEdge}},Vector{<:NamedEdge}}=belief_propagation_edge_sequence( undirected_graph(underlying_graph(mts)) ), -) where {E<:NamedEdge} +) return belief_propagation_iteration(tn, mts, edges; contract_kwargs, compute_norm) end From 045adc3f7ea8305cdd08ea9bb2c777ab88c236b9 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Thu, 9 Nov 2023 13:01:50 -0500 Subject: [PATCH 28/37] New File forScheduling and defaults for edge sequencing --- examples/belief_propagation/bpexample.jl | 4 +- examples/belief_propagation/bpsequences.jl | 250 +++++------------- examples/gauging/gauging_itns.jl | 5 +- src/ITensorNetworks.jl | 1 + src/beliefpropagation/beliefpropagation.jl | 55 ++-- .../beliefpropagation_schedule.jl | 51 ++++ .../sqrt_beliefpropagation.jl | 27 +- test/test_belief_propagation.jl | 9 +- 8 files changed, 154 insertions(+), 248 deletions(-) create mode 100644 src/beliefpropagation/beliefpropagation_schedule.jl diff --git a/examples/belief_propagation/bpexample.jl b/examples/belief_propagation/bpexample.jl index a86a0f53..a0af8818 100644 --- a/examples/belief_propagation/bpexample.jl +++ b/examples/belief_propagation/bpexample.jl @@ -35,7 +35,7 @@ function main() ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ)))) ) - mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact")) + mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact"), niters=20) numerator_network = approx_network_region( ψψ, mts, [(v, 1)]; verts_tn=ITensorNetwork([apply(op("Sz", s[v]), ψ[v])]) @@ -54,7 +54,7 @@ function main() ) Zpp = partition(ψψ; subgraph_vertices=nested_graph_leaf_vertices(Zp)) mts = message_tensors(Zpp) - mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact")) + mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact"), niters=20) numerator_network = approx_network_region( ψψ, mts, [(v, 1)]; verts_tn=ITensorNetwork([apply(op("Sz", s[v]), ψ[v])]) ) diff --git a/examples/belief_propagation/bpsequences.jl b/examples/belief_propagation/bpsequences.jl index fbcf0c14..819e55a0 100644 --- a/examples/belief_propagation/bpsequences.jl +++ b/examples/belief_propagation/bpsequences.jl @@ -12,196 +12,70 @@ using ITensorNetworks: approx_network_region, contract_inner, message_tensors, - nested_graph_leaf_vertices + nested_graph_leaf_vertices, + edge_sequence function main() - g = named_comb_tree((6, 6)) - s = siteinds("S=1/2", g) - χ = 4 - - Random.seed!(5467) - - ψ = randomITensorNetwork(s; link_space=χ) - ψψ = ψ ⊗ prime(dag(ψ); sites=[]) - - #Initial message tensors for BP - mts_init = message_tensors( - ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ)))) - ) - - println("\nFirst testing out a comb tree. Random network with bond dim $χ") - - #Now test out various sequences - print("Parallel updates (sequence is irrelevant): ") - belief_propagation( - ψψ, - mts_init; - contract_kwargs=(; alg="exact"), - target_precision=1e-10, - niters=100, - edges=[[e] for e in edges(mts_init)], - verbose=true, - ) - print("Sequential updates (sequence is default edge list of the message tensors): ") - belief_propagation( - ψψ, - mts_init; - contract_kwargs=(; alg="exact"), - target_precision=1e-10, - niters=100, - edges=[e for e in edges(mts_init)], - verbose=true, - ) - print("Sequential updates (sequence is our custom sequence finder): ") - belief_propagation( - ψψ, - mts_init; - contract_kwargs=(; alg="exact"), - target_precision=1e-10, - niters=100, - verbose=true, - ) - - g = NamedGraph(Graphs.random_regular_graph(100, 3)) - s = siteinds("S=1/2", g) - χ = 4 - - Random.seed!(5467) - - ψ = randomITensorNetwork(s; link_space=χ) - ψψ = ψ ⊗ prime(dag(ψ); sites=[]) - - #Initial message tensors for BP - mts_init = message_tensors( - ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ)))) - ) - - println("\nNow testing out a z = 3 random regular graph. Random network with bond dim $χ") - - #Now test out various sequences - print("Parallel updates (sequence is irrelevant): ") - belief_propagation( - ψψ, - mts_init; - contract_kwargs=(; alg="exact"), - target_precision=1e-10, - niters=100, - edges=[[e] for e in edges(mts_init)], - verbose=true, - ) - print("Sequential updates (sequence is default edge list of the message tensors): ") - belief_propagation( - ψψ, - mts_init; - contract_kwargs=(; alg="exact"), - target_precision=1e-10, - niters=100, - edges=[e for e in edges(mts_init)], - verbose=true, - ) - print("Sequential updates (sequence is our custom sequence finder): ") - belief_propagation( - ψψ, - mts_init; - contract_kwargs=(; alg="exact"), - target_precision=1e-10, - niters=100, - verbose=true, - ) - - g = named_grid((6, 6)) - s = siteinds("S=1/2", g) - χ = 2 - - Random.seed!(5467) - - ψ = randomITensorNetwork(s; link_space=χ) - ψψ = ψ ⊗ prime(dag(ψ); sites=[]) - - #Initial message tensors for BP - mts_init = message_tensors( - ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ)))) - ) - - println("\nNow testing out a 6x6 grid. Random network with bond dim $χ") - - #Now test out various sequences - print("Parallel updates (sequence is irrelevant): ") - belief_propagation( - ψψ, - mts_init; - contract_kwargs=(; alg="exact"), - target_precision=1e-10, - niters=100, - edges=[[e] for e in edges(mts_init)], - verbose=true, - ) - print("Sequential updates (sequence is default edge list of the message tensors): ") - belief_propagation( - ψψ, - mts_init; - contract_kwargs=(; alg="exact"), - target_precision=1e-10, - niters=100, - edges=[e for e in edges(mts_init)], - verbose=true, - ) - print("Sequential updates (sequence is our custom sequence finder): ") - belief_propagation( - ψψ, - mts_init; - contract_kwargs=(; alg="exact"), - target_precision=1e-10, - niters=100, - verbose=true, - ) - - g = NamedGraphs.hexagonal_lattice_graph(4, 4) - s = siteinds("S=1/2", g) - χ = 3 - - Random.seed!(5467) - - ψ = randomITensorNetwork(s; link_space=χ) - ψψ = ψ ⊗ prime(dag(ψ); sites=[]) - - #Initial message tensors for BP - mts_init = message_tensors( - ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ)))) - ) - - println("\nNow testing out a 4 x 4 hexagonal lattice. Random network with bond dim $χ") - - #Now test out various sequences - print("Parallel updates (sequence is irrelevant): ") - belief_propagation( - ψψ, - mts_init; - contract_kwargs=(; alg="exact"), - target_precision=1e-10, - niters=100, - edges=[[e] for e in edges(mts_init)], - verbose=true, - ) - print("Sequential updates (sequence is default edge list of the message tensors): ") - belief_propagation( - ψψ, - mts_init; - contract_kwargs=(; alg="exact"), - target_precision=1e-10, - niters=100, - edges=[e for e in edges(mts_init)], - verbose=true, - ) - print("Sequential updates (sequence is our custom sequence finder): ") - return belief_propagation( - ψψ, - mts_init; - contract_kwargs=(; alg="exact"), - target_precision=1e-10, - niters=100, - verbose=true, - ) + g_labels = [ + "Comb Tree", + "100 Site Random Regular Graph z = 3", + "6x6 Square Grid", + "4x4 Hexagonal Lattice", + ] + gs = [ + named_comb_tree((6, 6)), + NamedGraph(Graphs.random_regular_graph(100, 3)), + named_grid((6, 6)), + NamedGraphs.hexagonal_lattice_graph(4, 4), + ] + χs = [4, 4, 2, 3] + + for (i, g) in enumerate(gs) + Random.seed!(5467) + g_label = g_labels[i] + χ = χs[i] + s = siteinds("S=1/2", g) + ψ = randomITensorNetwork(s; link_space=χ) + ψψ = ψ ⊗ prime(dag(ψ); sites=[]) + + #Initial message tensors for BP + mts_init = message_tensors( + ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ)))) + ) + + println("\nFirst testing out a $g_label. Random network with bond dim $χ") + + #Now test out various sequences + print("Parallel updates (sequence is irrelevant): ") + belief_propagation( + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + edges=edge_sequence(mts_init; alg=ITensorNetworks.Parallel()), + verbose=true, + ) + print("Sequential updates (sequence is default edge list of the message tensors): ") + belief_propagation( + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + edges=[e for e in edges(mts_init)], + verbose=true, + ) + print("Sequential updates (sequence is our custom sequence finder): ") + belief_propagation( + ψψ, + mts_init; + contract_kwargs=(; alg="exact"), + target_precision=1e-10, + niters=100, + verbose=true, + ) + end end main() diff --git a/examples/gauging/gauging_itns.jl b/examples/gauging/gauging_itns.jl index 28fa6f63..3c217452 100644 --- a/examples/gauging/gauging_itns.jl +++ b/examples/gauging/gauging_itns.jl @@ -74,7 +74,10 @@ function benchmark_state_gauging( ) else times_iters[i] = @elapsed mts, _ = belief_propagation_iteration( - ψψ, mts; contract_kwargs=(; alg="exact"), edges=[[e] for e in edges(mts)] + ψψ, + mts; + contract_kwargs=(; alg="exact"), + edges=edge_sequence(mts; alg=ITensorNetworks.Parallel()), ) end diff --git a/src/ITensorNetworks.jl b/src/ITensorNetworks.jl index 55f261ed..5b660bc9 100644 --- a/src/ITensorNetworks.jl +++ b/src/ITensorNetworks.jl @@ -97,6 +97,7 @@ include("specialitensornetworks.jl") include("renameitensornetwork.jl") include("boundarymps.jl") include(joinpath("beliefpropagation", "beliefpropagation.jl")) +include(joinpath("beliefpropagation", "beliefpropagation_schedule.jl")) include(joinpath("beliefpropagation", "sqrt_beliefpropagation.jl")) include("contraction_tree_to_graph.jl") include("gauging.jl") diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index e168e5e2..9e55fcb2 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -22,7 +22,8 @@ function message_tensors_skeleton(subgraphs::DataGraph) end function message_tensors( - subgraphs::DataGraph; itensor_constructor=inds_e -> ITensor[dense(delta(i)) for i in inds_e] + subgraphs::DataGraph; + itensor_constructor=inds_e -> ITensor[dense(delta(i)) for i in inds_e], ) mts = message_tensors_skeleton(subgraphs) for e in edges(subgraphs) @@ -80,7 +81,7 @@ Do a sequential update of message tensors on `edges` for a given ITensornetwork function belief_propagation_iteration( tn::ITensorNetwork, mts::DataGraph, - edges::Vector{<:NamedEdge}; + edges::Vector{<:AbstractEdge}; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), compute_norm=false, ) @@ -107,15 +108,17 @@ function belief_propagation_iteration( end """ -Do parallel updates between groups of edges of all message tensors for a given ITensornetwork and its partition into sub graphs +Do parallel updates between groups of edges of all message tensors for a given ITensornetwork and its partition into sub graphs. +Currently we send the full message tensor data struct to belief_propagation_iteration for each subgraph. But really we only need the +mts relevant to that subgraph. """ function belief_propagation_iteration( tn::ITensorNetwork, mts::DataGraph, - edge_groups::Vector{Vector{E}}; + edge_groups::Vector{<:Vector{<:AbstractEdge}}; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), compute_norm=false, -) where {E<:NamedEdge} +) new_mts = copy(mts) c = 0 for edges in edge_groups @@ -135,9 +138,7 @@ function belief_propagation_iteration( mts::DataGraph; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), compute_norm=false, - edges::Union{Vector{Vector{<:NamedEdge}},Vector{<:NamedEdge}}=belief_propagation_edge_sequence( - undirected_graph(underlying_graph(mts)) - ), + edges=edge_sequence(mts), ) return belief_propagation_iteration(tn, mts, edges; contract_kwargs, compute_norm) end @@ -146,14 +147,15 @@ function belief_propagation( tn::ITensorNetwork, mts::DataGraph; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), - niters=20, - target_precision::Union{Float64,Nothing}=nothing, - edges::Union{Vector{Vector{E}},Vector{E}}=belief_propagation_edge_sequence( - undirected_graph(underlying_graph(mts)) - ), + niters::Union{Int64,Nothing}=default_bp_niters(mts), + target_precision=nothing, + edges=edge_sequence(mts), verbose=false, -) where {E<:NamedEdge} +) compute_norm = !isnothing(target_precision) + if isnothing(niters) + error("You need to specify a number of iterations for BP!") + end for i in 1:niters mts, c = belief_propagation_iteration(tn, mts, edges; contract_kwargs, compute_norm) if compute_norm && c <= target_precision @@ -172,8 +174,8 @@ function belief_propagation( nvertices_per_partition=nothing, npartitions=nothing, subgraph_vertices=nothing, - niters=20, - target_precision::Union{Float64,Nothing}=nothing, + niters=default_bp_niters(mts), + target_precision=nothing, verbose=false, ) mts = message_tensors(tn; nvertices_per_partition, npartitions, subgraph_vertices) @@ -212,24 +214,3 @@ function approx_network_region( return environment_tn ⊗ verts_tn end - -""" -Return a custom edge order for how to update all BP message tensors on a general undirected graph. -On a tree this will yield a sequence which only needs to be performed once. Based on forest covers and depth first searches amongst the forests. -""" -function belief_propagation_edge_sequence( - g::NamedGraph; root_vertex=NamedGraphs.default_root_vertex -) - @assert !is_directed(g) - forests = NamedGraphs.forest_cover(g) - edges = NamedEdge[] - for forest in forests - trees = NamedGraph[forest[vs] for vs in connected_components(forest)] - for tree in trees - tree_edges = post_order_dfs_edges(tree, root_vertex(tree)) - push!(edges, vcat(tree_edges, reverse(reverse.(tree_edges)))...) - end - end - - return edges -end diff --git a/src/beliefpropagation/beliefpropagation_schedule.jl b/src/beliefpropagation/beliefpropagation_schedule.jl new file mode 100644 index 00000000..a8e3bf66 --- /dev/null +++ b/src/beliefpropagation/beliefpropagation_schedule.jl @@ -0,0 +1,51 @@ +abstract type EdgeSequenceAlgorithm end + +struct ForestCover <: EdgeSequenceAlgorithm end +struct Parallel <: EdgeSequenceAlgorithm end + +default_edge_sequence_alg() = ForestCover() + +default_bp_niters(g::NamedGraph) = is_tree(g) ? 1 : nothing +function default_bp_niters(g::Union{DataGraph,AbstractITensorNetwork}) + return default_bp_niters(undirected_graph(underlying_graph(g))) +end + +function edge_sequence(g::NamedGraph; alg=default_edge_sequence_alg()) + return edge_sequence(alg, g) +end + +function edge_sequence(mts::DataGraph; alg=default_edge_sequence_alg()) + return edge_sequence(alg, undirected_graph(underlying_graph(mts))) +end + +function edge_sequence( + ::ForestCover, g::NamedGraph; root_vertex=NamedGraphs.default_root_vertex +) + @assert !is_directed(g) + forests = NamedGraphs.forest_cover(g) + edges = edgetype(g)[] + for forest in forests + trees = [forest[vs] for vs in connected_components(forest)] + for tree in trees + tree_edges = post_order_dfs_edges(tree, root_vertex(tree)) + push!(edges, vcat(tree_edges, reverse(reverse.(tree_edges)))...) + end + end + + return edges +end + +function edge_sequence( + ::ForestCover, mts::DataGraph; root_vertex=NamedGraphs.default_root_vertex +) + return edge_sequence(ForestCover, undirected_graph(underlying_graph(mts)); root_vertex) +end + +function edge_sequence(::Parallel, g::NamedGraph) + @assert !is_directed(g) + return [[e] for e in vcat(edges(g), reverse.(edges(g)))] +end + +function edge_sequence(::Parallel, mts::DataGraph) + return edge_sequence(Parallel, undirected_graph(underlying_graph(mts))) +end diff --git a/src/beliefpropagation/sqrt_beliefpropagation.jl b/src/beliefpropagation/sqrt_beliefpropagation.jl index f28ab9a0..6459cfc3 100644 --- a/src/beliefpropagation/sqrt_beliefpropagation.jl +++ b/src/beliefpropagation/sqrt_beliefpropagation.jl @@ -2,8 +2,8 @@ # using ITensorNetworks: find_subgraph, map_diag, sqrt_diag, boundary_edges function sqrt_belief_propagation_iteration( - tn::ITensorNetwork, sqrt_mts::DataGraph, edges::Vector{E} -) where {E<:NamedEdge} + tn::ITensorNetwork, sqrt_mts::DataGraph, edges::Vector{<:AbstractEdge} +) new_sqrt_mts = copy(sqrt_mts) c = 0.0 for e in edges @@ -28,8 +28,8 @@ function sqrt_belief_propagation_iteration( end function sqrt_belief_propagation_iteration( - tn::ITensorNetwork, sqrt_mts::DataGraph, edges::Vector{Vector{E}} -) where {E<:NamedEdge} + tn::ITensorNetwork, sqrt_mts::DataGraph, edges::Vector{<:Vector{<:AbstractEdge}} +) new_sqrt_mts = copy(sqrt_mts) c = 0.0 for e_group in edges @@ -43,26 +43,23 @@ function sqrt_belief_propagation_iteration( end function sqrt_belief_propagation_iteration( - tn::ITensorNetwork, - sqrt_mts::DataGraph; - edges::Union{Vector{Vector{E}},Vector{E}}=belief_propagation_edge_sequence( - undirected_graph(underlying_graph(mts)) - ), -) where {E<:NamedEdge} + tn::ITensorNetwork, sqrt_mts::DataGraph; edges=edge_sequence(mts) +) return sqrt_belief_propagation_iteration(tn, sqrt_mts, edges) end function sqrt_belief_propagation( tn::ITensorNetwork, mts::DataGraph; - niters=20, - edges::Union{Vector{Vector{E}},Vector{E}}=belief_propagation_edge_sequence( - undirected_graph(underlying_graph(mts)) - ), + niters::Union{Int64,Nothing}=default_bp_niters(mts), + edges=edge_sequence(mts), # target_precision::Union{Float64,Nothing}=nothing, -) where {E<:NamedEdge} +) # compute_norm = target_precision == nothing ? false : true sqrt_mts = sqrt_message_tensors(tn, mts) + if isnothing(niters) + error("You need to specify a number of iterations for BP!") + end for i in 1:niters sqrt_mts, c = sqrt_belief_propagation_iteration(tn, sqrt_mts, edges) #; compute_norm) # if compute_norm && c <= target_precision diff --git a/test/test_belief_propagation.jl b/test/test_belief_propagation.jl index 82a2db66..d3747598 100644 --- a/test/test_belief_propagation.jl +++ b/test/test_belief_propagation.jl @@ -40,7 +40,7 @@ ITensors.disable_warn_order() Z = partition(ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ))))) mts = message_tensors(Z) - mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact"), niters=1) + mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact")) numerator_network = approx_network_region( ψψ, mts, [(v, 1)]; verts_tn=ITensorNetwork(ITensor[apply(op("Sz", s[v]), ψ[v])]) @@ -67,7 +67,7 @@ ITensors.disable_warn_order() Z = partition(ψψ; subgraph_vertices=collect(values(group(v -> v[1], vertices(ψψ))))) mts = message_tensors(Z) - mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact"), niters=1) + mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact")) numerator_network = approx_network_region( ψψ, mts, [(v, 1)]; verts_tn=ITensorNetwork(ITensor[apply(op("Sz", s[v]), ψ[v])]) @@ -94,7 +94,7 @@ ITensors.disable_warn_order() nsites = 2 Z = partition(ψψ; nvertices_per_partition=nsites) mts = message_tensors(Z) - mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact")) + mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact"), niters=20) numerator_network = approx_network_region( ψψ, mts, vs; verts_tn=ITensorNetwork(ITensor[ψOψ[v] for v in vs]) ) @@ -120,7 +120,7 @@ ITensors.disable_warn_order() ) Zpp = partition(ψψ; subgraph_vertices=nested_graph_leaf_vertices(Zp)) mts = message_tensors(Zpp) - mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact")) + mts = belief_propagation(ψψ, mts; contract_kwargs=(; alg="exact"), niters=20) ψψsplit = split_index(ψψ, NamedEdge.([(v, 1) => (v, 2) for v in vs])) rdm = ITensors.contract( @@ -165,7 +165,6 @@ ITensors.disable_warn_order() contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, cutoff=1e-16, maxdim ), - niters=1, ) numerator_network = approx_network_region(ψψ, mts, [v]; verts_tn=ITensorNetwork(ψOψ[v])) From e7ccf43cc359a3099c302ddf02a78aadae10fa27 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Thu, 9 Nov 2023 15:49:58 -0500 Subject: [PATCH 29/37] Improved Schedule Code --- examples/belief_propagation/bpsequences.jl | 2 +- examples/gauging/gauging_itns.jl | 5 +-- .../beliefpropagation_schedule.jl | 31 ++++++------------- .../sqrt_beliefpropagation.jl | 2 +- .../test_solvers/test_dmrg.jl | 4 +-- 5 files changed, 17 insertions(+), 27 deletions(-) diff --git a/examples/belief_propagation/bpsequences.jl b/examples/belief_propagation/bpsequences.jl index 819e55a0..cbdf7148 100644 --- a/examples/belief_propagation/bpsequences.jl +++ b/examples/belief_propagation/bpsequences.jl @@ -53,7 +53,7 @@ function main() contract_kwargs=(; alg="exact"), target_precision=1e-10, niters=100, - edges=edge_sequence(mts_init; alg=ITensorNetworks.Parallel()), + edges=edge_sequence(mts_init; alg=ITensors.Algorithm("Parallel")), verbose=true, ) print("Sequential updates (sequence is default edge list of the message tensors): ") diff --git a/examples/gauging/gauging_itns.jl b/examples/gauging/gauging_itns.jl index 3c217452..0cf10931 100644 --- a/examples/gauging/gauging_itns.jl +++ b/examples/gauging/gauging_itns.jl @@ -17,7 +17,8 @@ using ITensorNetworks: vidal_to_symmetric_gauge, initialize_bond_tensors, vidal_itn_isometries, - norm_network + norm_network, + edge_sequence using NamedGraphs using NamedGraphs: add_edges!, rem_vertex!, hexagonal_lattice_graph @@ -77,7 +78,7 @@ function benchmark_state_gauging( ψψ, mts; contract_kwargs=(; alg="exact"), - edges=edge_sequence(mts; alg=ITensorNetworks.Parallel()), + edges=edge_sequence(mts; alg=ITensors.Algorithm("Parallel")), ) end diff --git a/src/beliefpropagation/beliefpropagation_schedule.jl b/src/beliefpropagation/beliefpropagation_schedule.jl index a8e3bf66..6c5c29fc 100644 --- a/src/beliefpropagation/beliefpropagation_schedule.jl +++ b/src/beliefpropagation/beliefpropagation_schedule.jl @@ -1,12 +1,7 @@ -abstract type EdgeSequenceAlgorithm end - -struct ForestCover <: EdgeSequenceAlgorithm end -struct Parallel <: EdgeSequenceAlgorithm end - -default_edge_sequence_alg() = ForestCover() +default_edge_sequence_alg() = Algorithm("ForestCover") default_bp_niters(g::NamedGraph) = is_tree(g) ? 1 : nothing -function default_bp_niters(g::Union{DataGraph,AbstractITensorNetwork}) +function default_bp_niters(g::AbstractGraph) return default_bp_niters(undirected_graph(underlying_graph(g))) end @@ -14,12 +9,16 @@ function edge_sequence(g::NamedGraph; alg=default_edge_sequence_alg()) return edge_sequence(alg, g) end -function edge_sequence(mts::DataGraph; alg=default_edge_sequence_alg()) - return edge_sequence(alg, undirected_graph(underlying_graph(mts))) +function edge_sequence(g::AbstractGraph; alg=default_edge_sequence_alg()) + return edge_sequence(alg, undirected_graph(underlying_graph(g))) +end + +function edge_sequence(alg::Algorithm, g::AbstractGraph; kwargs...) + return edge_sequence(alg, g, kwargs...) end function edge_sequence( - ::ForestCover, g::NamedGraph; root_vertex=NamedGraphs.default_root_vertex + ::Algorithm"ForestCover", g::NamedGraph; root_vertex=NamedGraphs.default_root_vertex ) @assert !is_directed(g) forests = NamedGraphs.forest_cover(g) @@ -35,17 +34,7 @@ function edge_sequence( return edges end -function edge_sequence( - ::ForestCover, mts::DataGraph; root_vertex=NamedGraphs.default_root_vertex -) - return edge_sequence(ForestCover, undirected_graph(underlying_graph(mts)); root_vertex) -end - -function edge_sequence(::Parallel, g::NamedGraph) +function edge_sequence(::Algorithm"Parallel", g::NamedGraph) @assert !is_directed(g) return [[e] for e in vcat(edges(g), reverse.(edges(g)))] end - -function edge_sequence(::Parallel, mts::DataGraph) - return edge_sequence(Parallel, undirected_graph(underlying_graph(mts))) -end diff --git a/src/beliefpropagation/sqrt_beliefpropagation.jl b/src/beliefpropagation/sqrt_beliefpropagation.jl index 6459cfc3..20ff3d8e 100644 --- a/src/beliefpropagation/sqrt_beliefpropagation.jl +++ b/src/beliefpropagation/sqrt_beliefpropagation.jl @@ -51,7 +51,7 @@ end function sqrt_belief_propagation( tn::ITensorNetwork, mts::DataGraph; - niters::Union{Int64,Nothing}=default_bp_niters(mts), + niters=default_bp_niters(mts), edges=edge_sequence(mts), # target_precision::Union{Float64,Nothing}=nothing, ) diff --git a/test/test_treetensornetworks/test_solvers/test_dmrg.jl b/test/test_treetensornetworks/test_solvers/test_dmrg.jl index 61c2febc..939b4fee 100644 --- a/test/test_treetensornetworks/test_solvers/test_dmrg.jl +++ b/test/test_treetensornetworks/test_solvers/test_dmrg.jl @@ -28,7 +28,7 @@ using Observers # Compare to `ITensors.MPO` version of `dmrg` H_mpo = MPO([H[v] for v in 1:nv(H)]) psi_mps = MPS([psi[v] for v in 1:nv(psi)]) - e2, psi2 = dmrg(H_mpo, psi_mps; nsweeps, maxdim, normalize=false, outputlevel=0) + e2, psi2 = dmrg(H_mpo, psi_mps; nsweeps, maxdim, outputlevel=0) psi = dmrg(H, psi; nsweeps, maxdim, cutoff, nsite, solver_krylovdim=3, solver_maxiter=1) @test inner(psi', H, psi) ≈ inner(psi2', H_mpo, psi2) @@ -134,7 +134,7 @@ end sline = only.(collect(vertex_data(s)))[linear_order] Hline = MPO(relabel_sites(os, vmap), sline) psiline = randomMPS(sline; linkdims=20) - e2, psi2 = dmrg(Hline, psiline, sweeps; normalize=false, outputlevel=0) + e2, psi2 = dmrg(Hline, psiline, sweeps; outputlevel=0) @test inner(psi', H, psi) ≈ inner(psi2', Hline, psi2) atol = 1e-5 end From 2b8513a7f663533fbbfd32c4d28bdb314e0043f1 Mon Sep 17 00:00:00 2001 From: Joseph Tindall <51231103+JoeyT1994@users.noreply.github.com> Date: Thu, 9 Nov 2023 16:28:21 -0500 Subject: [PATCH 30/37] Update src/beliefpropagation/beliefpropagation_schedule.jl Co-authored-by: Matt Fishman --- src/beliefpropagation/beliefpropagation_schedule.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/beliefpropagation/beliefpropagation_schedule.jl b/src/beliefpropagation/beliefpropagation_schedule.jl index 6c5c29fc..37417f2c 100644 --- a/src/beliefpropagation/beliefpropagation_schedule.jl +++ b/src/beliefpropagation/beliefpropagation_schedule.jl @@ -10,7 +10,7 @@ function edge_sequence(g::NamedGraph; alg=default_edge_sequence_alg()) end function edge_sequence(g::AbstractGraph; alg=default_edge_sequence_alg()) - return edge_sequence(alg, undirected_graph(underlying_graph(g))) + return edge_sequence(Algorithm(alg), undirected_graph(underlying_graph(g))) end function edge_sequence(alg::Algorithm, g::AbstractGraph; kwargs...) From 1b7189481fcbafd05798ca258a900781e7df4a32 Mon Sep 17 00:00:00 2001 From: Joseph Tindall <51231103+JoeyT1994@users.noreply.github.com> Date: Thu, 9 Nov 2023 16:29:29 -0500 Subject: [PATCH 31/37] Update src/beliefpropagation/beliefpropagation_schedule.jl Co-authored-by: Matt Fishman --- src/beliefpropagation/beliefpropagation_schedule.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/beliefpropagation/beliefpropagation_schedule.jl b/src/beliefpropagation/beliefpropagation_schedule.jl index 37417f2c..5dc1f193 100644 --- a/src/beliefpropagation/beliefpropagation_schedule.jl +++ b/src/beliefpropagation/beliefpropagation_schedule.jl @@ -34,7 +34,7 @@ function edge_sequence( return edges end -function edge_sequence(::Algorithm"Parallel", g::NamedGraph) +function edge_sequence(::Algorithm"parallel", g::NamedGraph) @assert !is_directed(g) return [[e] for e in vcat(edges(g), reverse.(edges(g)))] end From fe06b397dd12d17caf111c45485f297cf683b9c9 Mon Sep 17 00:00:00 2001 From: Joseph Tindall <51231103+JoeyT1994@users.noreply.github.com> Date: Thu, 9 Nov 2023 16:29:37 -0500 Subject: [PATCH 32/37] Update src/beliefpropagation/beliefpropagation_schedule.jl Co-authored-by: Matt Fishman --- src/beliefpropagation/beliefpropagation_schedule.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/beliefpropagation/beliefpropagation_schedule.jl b/src/beliefpropagation/beliefpropagation_schedule.jl index 5dc1f193..53113d44 100644 --- a/src/beliefpropagation/beliefpropagation_schedule.jl +++ b/src/beliefpropagation/beliefpropagation_schedule.jl @@ -18,7 +18,7 @@ function edge_sequence(alg::Algorithm, g::AbstractGraph; kwargs...) end function edge_sequence( - ::Algorithm"ForestCover", g::NamedGraph; root_vertex=NamedGraphs.default_root_vertex + ::Algorithm"forest_cover", g::NamedGraph; root_vertex=NamedGraphs.default_root_vertex ) @assert !is_directed(g) forests = NamedGraphs.forest_cover(g) From b21f629c6052e45be4eaf88cf4a100563bfee9b3 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Thu, 9 Nov 2023 16:41:09 -0500 Subject: [PATCH 33/37] Imported Algorithmn from ITensors --- examples/belief_propagation/bpsequences.jl | 2 +- examples/gauging/gauging_itns.jl | 5 +---- src/beliefpropagation/beliefpropagation_schedule.jl | 2 +- src/imports.jl | 2 ++ 4 files changed, 5 insertions(+), 6 deletions(-) diff --git a/examples/belief_propagation/bpsequences.jl b/examples/belief_propagation/bpsequences.jl index cbdf7148..456e96c7 100644 --- a/examples/belief_propagation/bpsequences.jl +++ b/examples/belief_propagation/bpsequences.jl @@ -53,7 +53,7 @@ function main() contract_kwargs=(; alg="exact"), target_precision=1e-10, niters=100, - edges=edge_sequence(mts_init; alg=ITensors.Algorithm("Parallel")), + edges=edge_sequence(mts_init; alg="parallel"), verbose=true, ) print("Sequential updates (sequence is default edge list of the message tensors): ") diff --git a/examples/gauging/gauging_itns.jl b/examples/gauging/gauging_itns.jl index 0cf10931..0246e443 100644 --- a/examples/gauging/gauging_itns.jl +++ b/examples/gauging/gauging_itns.jl @@ -75,10 +75,7 @@ function benchmark_state_gauging( ) else times_iters[i] = @elapsed mts, _ = belief_propagation_iteration( - ψψ, - mts; - contract_kwargs=(; alg="exact"), - edges=edge_sequence(mts; alg=ITensors.Algorithm("Parallel")), + ψψ, mts; contract_kwargs=(; alg="exact"), edges=edge_sequence(mts; alg="parallel") ) end diff --git a/src/beliefpropagation/beliefpropagation_schedule.jl b/src/beliefpropagation/beliefpropagation_schedule.jl index 53113d44..a4f5eb84 100644 --- a/src/beliefpropagation/beliefpropagation_schedule.jl +++ b/src/beliefpropagation/beliefpropagation_schedule.jl @@ -1,4 +1,4 @@ -default_edge_sequence_alg() = Algorithm("ForestCover") +default_edge_sequence_alg() = "forest_cover" default_bp_niters(g::NamedGraph) = is_tree(g) ? 1 : nothing function default_bp_niters(g::AbstractGraph) diff --git a/src/imports.jl b/src/imports.jl index 3186de4d..71406878 100644 --- a/src/imports.jl +++ b/src/imports.jl @@ -96,6 +96,8 @@ import ITensors: scalartype, #adding add +#Algorithm +Algorithm using ITensors.ContractionSequenceOptimization: deepmap From c82f92fe55eae4da4bdd4989a610adaf80c17420 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Thu, 9 Nov 2023 16:54:29 -0500 Subject: [PATCH 34/37] Better dispatching on graph types --- src/beliefpropagation/beliefpropagation_schedule.jl | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/src/beliefpropagation/beliefpropagation_schedule.jl b/src/beliefpropagation/beliefpropagation_schedule.jl index a4f5eb84..fb79d011 100644 --- a/src/beliefpropagation/beliefpropagation_schedule.jl +++ b/src/beliefpropagation/beliefpropagation_schedule.jl @@ -1,6 +1,6 @@ default_edge_sequence_alg() = "forest_cover" -default_bp_niters(g::NamedGraph) = is_tree(g) ? 1 : nothing +@traitfn default_bp_niters(g::NamedGraph::(!IsDirected)) = is_tree(g) ? 1 : nothing function default_bp_niters(g::AbstractGraph) return default_bp_niters(undirected_graph(underlying_graph(g))) end @@ -17,10 +17,11 @@ function edge_sequence(alg::Algorithm, g::AbstractGraph; kwargs...) return edge_sequence(alg, g, kwargs...) end -function edge_sequence( - ::Algorithm"forest_cover", g::NamedGraph; root_vertex=NamedGraphs.default_root_vertex +@traitfn function edge_sequence( + ::Algorithm"forest_cover", + g::NamedGraph::(!IsDirected); + root_vertex=NamedGraphs.default_root_vertex, ) - @assert !is_directed(g) forests = NamedGraphs.forest_cover(g) edges = edgetype(g)[] for forest in forests @@ -34,7 +35,6 @@ function edge_sequence( return edges end -function edge_sequence(::Algorithm"parallel", g::NamedGraph) - @assert !is_directed(g) +@traitfn function edge_sequence(::Algorithm"parallel", g::NamedGraph::(!IsDirected)) return [[e] for e in vcat(edges(g), reverse.(edges(g)))] end From 18321b681217f86f45ac65708d9fac0dd3a75a06 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Mon, 20 Nov 2023 17:05:30 -0500 Subject: [PATCH 35/37] Fixed NamedGraph type for certain operations --- examples/gauging/gauging_itns.jl | 8 ++++---- .../beliefpropagation_schedule.jl | 17 ++++++++++------- 2 files changed, 14 insertions(+), 11 deletions(-) diff --git a/examples/gauging/gauging_itns.jl b/examples/gauging/gauging_itns.jl index 0246e443..7463d54a 100644 --- a/examples/gauging/gauging_itns.jl +++ b/examples/gauging/gauging_itns.jl @@ -46,7 +46,7 @@ end """Bring an ITN into the Vidal gauge, various methods possible. Result is timed""" function benchmark_state_gauging( ψ::ITensorNetwork; - mode="BeliefPropagation", + mode="belief_propagation", no_iterations=50, BP_update_order::String="sequential", ) @@ -68,7 +68,7 @@ function benchmark_state_gauging( for i in 1:no_iterations println("On Iteration " * string(i)) - if mode == "BeliefPropagation" + if mode == "belief_propagation" if BP_update_order != "parallel" times_iters[i] = @elapsed mts, _ = belief_propagation_iteration( ψψ, mts; contract_kwargs=(; alg="exact") @@ -80,7 +80,7 @@ function benchmark_state_gauging( end times_gauging[i] = @elapsed ψ, bond_tensors = vidal_gauge(ψinit, mts) - elseif mode == "Eager" + elseif mode == "eager" times_iters[i] = @elapsed ψ, bond_tensors, mts = eager_gauging(ψ, bond_tensors, mts) else times_iters[i] = @elapsed begin @@ -110,7 +110,7 @@ BPG_simulation_times, BPG_Cs = benchmark_state_gauging( BPG_sequential_simulation_times, BPG_sequential_Cs = benchmark_state_gauging( ψ; no_iterations ) -Eager_simulation_times, Eager_Cs = benchmark_state_gauging(ψ; mode="Eager", no_iterations) +Eager_simulation_times, Eager_Cs = benchmark_state_gauging(ψ; mode="eager", no_iterations) SU_simulation_times, SU_Cs = benchmark_state_gauging(ψ; mode="SU", no_iterations) epsilon = 1e-10 diff --git a/src/beliefpropagation/beliefpropagation_schedule.jl b/src/beliefpropagation/beliefpropagation_schedule.jl index fb79d011..7cda8325 100644 --- a/src/beliefpropagation/beliefpropagation_schedule.jl +++ b/src/beliefpropagation/beliefpropagation_schedule.jl @@ -1,20 +1,23 @@ default_edge_sequence_alg() = "forest_cover" -@traitfn default_bp_niters(g::NamedGraph::(!IsDirected)) = is_tree(g) ? 1 : nothing +@traitfn undirected_default_bp_niters(g::AbstractGraph::(!IsDirected)) = + is_tree(g) ? 1 : nothing function default_bp_niters(g::AbstractGraph) - return default_bp_niters(undirected_graph(underlying_graph(g))) + return undirected_default_bp_niters(undirected_graph(underlying_graph(g))) end -function edge_sequence(g::NamedGraph; alg=default_edge_sequence_alg()) - return edge_sequence(alg, g) +@traitfn function edge_sequence( + g::NamedGraph::(!IsDirected); alg=default_edge_sequence_alg(), kwargs... +) + return edge_sequence(Algorithm(alg), g; kwargs...) end -function edge_sequence(g::AbstractGraph; alg=default_edge_sequence_alg()) - return edge_sequence(Algorithm(alg), undirected_graph(underlying_graph(g))) +function edge_sequence(g::AbstractGraph; alg=default_edge_sequence_alg(), kwargs...) + return edge_sequence(Algorithm(alg), undirected_graph(underlying_graph(g)); kwargs...) end function edge_sequence(alg::Algorithm, g::AbstractGraph; kwargs...) - return edge_sequence(alg, g, kwargs...) + return edge_sequence(alg, undirected_graph(underlying_graph(g)); kwargs...) end @traitfn function edge_sequence( From a3db9c01999d45febd41684cf96a64e8f9809b1e Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Tue, 21 Nov 2023 11:26:08 -0500 Subject: [PATCH 36/37] Trait fns now come in pairs --- .../beliefpropagation_schedule.jl | 19 ++++++++----------- 1 file changed, 8 insertions(+), 11 deletions(-) diff --git a/src/beliefpropagation/beliefpropagation_schedule.jl b/src/beliefpropagation/beliefpropagation_schedule.jl index 7cda8325..cbc5f206 100644 --- a/src/beliefpropagation/beliefpropagation_schedule.jl +++ b/src/beliefpropagation/beliefpropagation_schedule.jl @@ -1,29 +1,26 @@ default_edge_sequence_alg() = "forest_cover" -@traitfn undirected_default_bp_niters(g::AbstractGraph::(!IsDirected)) = - is_tree(g) ? 1 : nothing -function default_bp_niters(g::AbstractGraph) - return undirected_default_bp_niters(undirected_graph(underlying_graph(g))) +@traitfn default_bp_niters(g::::(!IsDirected)) = is_tree(g) ? 1 : nothing +@traitfn function default_bp_niters(g::::IsDirected) + return default_bp_niters(undirected_graph(underlying_graph(g))) end @traitfn function edge_sequence( - g::NamedGraph::(!IsDirected); alg=default_edge_sequence_alg(), kwargs... + g::::(!IsDirected); alg=default_edge_sequence_alg(), kwargs... ) return edge_sequence(Algorithm(alg), g; kwargs...) end -function edge_sequence(g::AbstractGraph; alg=default_edge_sequence_alg(), kwargs...) +@traitfn function edge_sequence(g::::IsDirected; alg=default_edge_sequence_alg(), kwargs...) return edge_sequence(Algorithm(alg), undirected_graph(underlying_graph(g)); kwargs...) end -function edge_sequence(alg::Algorithm, g::AbstractGraph; kwargs...) +@traitfn function edge_sequence(alg::Algorithm, g::::IsDirected; kwargs...) return edge_sequence(alg, undirected_graph(underlying_graph(g)); kwargs...) end @traitfn function edge_sequence( - ::Algorithm"forest_cover", - g::NamedGraph::(!IsDirected); - root_vertex=NamedGraphs.default_root_vertex, + ::Algorithm"forest_cover", g::::(!IsDirected); root_vertex=NamedGraphs.default_root_vertex ) forests = NamedGraphs.forest_cover(g) edges = edgetype(g)[] @@ -38,6 +35,6 @@ end return edges end -@traitfn function edge_sequence(::Algorithm"parallel", g::NamedGraph::(!IsDirected)) +@traitfn function edge_sequence(::Algorithm"parallel", g::::(!IsDirected)) return [[e] for e in vcat(edges(g), reverse.(edges(g)))] end From 487c43ffd5c4bb38660bcb658cb899cf48eb6de5 Mon Sep 17 00:00:00 2001 From: Joseph Tindall Date: Tue, 28 Nov 2023 09:57:43 -0500 Subject: [PATCH 37/37] Type assertion removed --- src/beliefpropagation/beliefpropagation.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/beliefpropagation/beliefpropagation.jl b/src/beliefpropagation/beliefpropagation.jl index 9e55fcb2..8a1f17b4 100644 --- a/src/beliefpropagation/beliefpropagation.jl +++ b/src/beliefpropagation/beliefpropagation.jl @@ -147,7 +147,7 @@ function belief_propagation( tn::ITensorNetwork, mts::DataGraph; contract_kwargs=(; alg="density_matrix", output_structure=path_graph_structure, maxdim=1), - niters::Union{Int64,Nothing}=default_bp_niters(mts), + niters=default_bp_niters(mts), target_precision=nothing, edges=edge_sequence(mts), verbose=false,