Add support for specifying RNGs in tensor network constructors, use StableRNGs.jl in tests

Project.toml

@@ -1,7 +1,7 @@
 name = "ITensorNetworks"
 uuid = "2919e153-833c-4bdc-8836-1ea460a35fc7"
 authors = ["Matthew Fishman <[email protected]>, Joseph Tindall <[email protected]> and contributors"]
-version = "0.11.9"
+version = "0.11.10"
 
 [deps]
 AbstractTrees = "1520ce14-60c1-5f80-bbc7-55ef81b5835c"

src/specialitensornetworks.jl

@@ -2,6 +2,7 @@ using ITensors: delta
 using ITensors.NDTensors: dim
 using DataGraphs: IsUnderlyingGraph
 using Distributions: Distribution
+using Random: Random, AbstractRNG
 
 """
 RETURN A TENSOR NETWORK WITH COPY TENSORS ON EACH VERTEX. 
@@ -28,37 +29,67 @@ end
 """
 Build an ITensor network on a graph specified by the inds network s. Bond_dim is given by link_space and entries are randomised (normal distribution, mean 0 std 1)
 """
-function random_tensornetwork(eltype::Type, s::IndsNetwork; kwargs...)
+function random_tensornetwork(rng::AbstractRNG, eltype::Type, s::IndsNetwork; kwargs...)
   return ITensorNetwork(s; kwargs...) do v
-    return inds -> itensor(randn(eltype, dim.(inds)...), inds)
+    return inds -> itensor(randn(rng, eltype, dim.(inds)...), inds)
   end
 end
 
+function random_tensornetwork(eltype::Type, s::IndsNetwork; kwargs...)
+  return random_tensornetwork(Random.default_rng(), eltype, s; kwargs...)
+end
+
+function random_tensornetwork(rng::AbstractRNG, s::IndsNetwork; kwargs...)
+  return random_tensornetwork(rng, Float64, s; kwargs...)
+end
+
 function random_tensornetwork(s::IndsNetwork; kwargs...)
-  return random_tensornetwork(Float64, s; kwargs...)
+  return random_tensornetwork(Random.default_rng(), s; kwargs...)
+end
+
+@traitfn function random_tensornetwork(
+  rng::AbstractRNG, eltype::Type, g::::IsUnderlyingGraph; kwargs...
+)
+  return random_tensornetwork(rng, eltype, IndsNetwork(g); kwargs...)
 end
 
 @traitfn function random_tensornetwork(eltype::Type, g::::IsUnderlyingGraph; kwargs...)
-  return random_tensornetwork(eltype, IndsNetwork(g); kwargs...)
+  return random_tensornetwork(Random.default_rng(), eltype, g; kwargs...)
+end
+
+@traitfn function random_tensornetwork(rng::AbstractRNG, g::::IsUnderlyingGraph; kwargs...)
+  return random_tensornetwork(rng, Float64, g; kwargs...)
 end
 
 @traitfn function random_tensornetwork(g::::IsUnderlyingGraph; kwargs...)
-  return random_tensornetwork(Float64, IndsNetwork(g); kwargs...)
+  return random_tensornetwork(Random.default_rng(), g; kwargs...)
 end
 
 """
 Build an ITensor network on a graph specified by the inds network s.
 Bond_dim is given by link_space and entries are randomized.
 The random distribution is based on the input argument `distribution`.
 """
-function random_tensornetwork(distribution::Distribution, s::IndsNetwork; kwargs...)
+function random_tensornetwork(
+  rng::AbstractRNG, distribution::Distribution, s::IndsNetwork; kwargs...
+)
   return ITensorNetwork(s; kwargs...) do v
-    return inds -> itensor(rand(distribution, dim.(inds)...), inds)
+    return inds -> itensor(rand(rng, distribution, dim.(inds)...), inds)
   end
 end
 
+function random_tensornetwork(distribution::Distribution, s::IndsNetwork; kwargs...)
+  return random_tensornetwork(Random.default_rng(), distribution, s; kwargs...)
+end
+
+@traitfn function random_tensornetwork(
+  rng::AbstractRNG, distribution::Distribution, g::::IsUnderlyingGraph; kwargs...
+)
+  return random_tensornetwork(rng, distribution, IndsNetwork(g); kwargs...)
+end
+
 @traitfn function random_tensornetwork(
   distribution::Distribution, g::::IsUnderlyingGraph; kwargs...
 )
-  return random_tensornetwork(distribution, IndsNetwork(g); kwargs...)
+  return random_tensornetwork(Random.default_rng(), distribution, g; kwargs...)
 end

test/Project.toml

@@ -24,6 +24,7 @@ OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
 Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 SplitApplyCombine = "03a91e81-4c3e-53e1-a0a4-9c0c8f19dd66"
+StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3"
 Suppressor = "fd094767-a336-5f1f-9728-57cf17d0bbfb"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 UnicodePlots = "b8865327-cd53-5732-bb35-84acbb429228"

test/test_additensornetworks.jl

@@ -5,11 +5,9 @@ using NamedGraphs.NamedGraphGenerators: named_grid
 using ITensorNetworks: ITensorNetwork, inner_network, random_tensornetwork, siteinds
 using ITensors: ITensors, apply, op, scalar, inner
 using LinearAlgebra: norm_sqr
-using Random: Random
+using StableRNGs: StableRNG
 using Test: @test, @testset
-
 @testset "add_itensornetworks" begin
-  Random.seed!(5623)
   g = named_grid((2, 2))
   s = siteinds("S=1/2", g)
   ψ1 = ITensorNetwork(v -> "↑", s)
@@ -32,8 +30,9 @@ using Test: @test, @testset
   rem_edge!(s2, NamedEdge((1, 1) => (1, 2)))
 
   v = rand(vertices(g))
-  ψ1 = random_tensornetwork(s1; link_space=χ)
-  ψ2 = random_tensornetwork(s2; link_space=χ)
+  rng = StableRNG(1234)
+  ψ1 = random_tensornetwork(rng, s1; link_space=χ)
+  ψ2 = random_tensornetwork(rng, s2; link_space=χ)
 
   ψ12 = ψ1 + ψ2
 

test/test_apply.jl

@@ -14,39 +14,32 @@ using ITensorNetworks:
 using ITensors: ITensors, inner, op
 using NamedGraphs.NamedGraphGenerators: named_grid
 using NamedGraphs.PartitionedGraphs: PartitionVertex
-using Random: Random
 using SplitApplyCombine: group
+using StableRNGs: StableRNG
 using Test: @test, @testset
-
 @testset "apply" begin
-  Random.seed!(5623)
   g_dims = (2, 2)
   n = prod(g_dims)
   g = named_grid(g_dims)
   s = siteinds("S=1/2", g)
   χ = 2
-  ψ = random_tensornetwork(s; link_space=χ)
+  rng = StableRNG(1234)
+  ψ = random_tensornetwork(rng, s; link_space=χ)
   v1, v2 = (2, 2), (1, 2)
   ψψ = norm_sqr_network(ψ)
   #Simple Belief Propagation Grouping
   bp_cache = BeliefPropagationCache(ψψ, group(v -> v[1], vertices(ψψ)))
   bp_cache = update(bp_cache; maxiter=20)
   envsSBP = environment(bp_cache, PartitionVertex.([v1, v2]))
-
   ψv = VidalITensorNetwork(ψ)
-
   #This grouping will correspond to calculating the environments exactly (each column of the grid is a partition)
   bp_cache = BeliefPropagationCache(ψψ, group(v -> v[1][1], vertices(ψψ)))
   bp_cache = update(bp_cache; maxiter=20)
   envsGBP = environment(bp_cache, [(v1, "bra"), (v1, "ket"), (v2, "bra"), (v2, "ket")])
-
   inner_alg = "exact"
-
   ngates = 5
-
   for i in 1:ngates
     o = op("RandomUnitary", s[v1]..., s[v2]...)
-
     ψOexact = apply(o, ψ; cutoff=1e-16)
     ψOSBP = apply(
       o,
@@ -79,9 +72,7 @@ using Test: @test, @testset
     fGBP =
       inner(ψOGBP, ψOexact; alg=inner_alg) /
       sqrt(inner(ψOexact, ψOexact; alg=inner_alg) * inner(ψOGBP, ψOGBP; alg=inner_alg))
-
     @test real(fGBP * conj(fGBP)) >= real(fSBP * conj(fSBP))
-
     @test isapprox(real(fSBP * conj(fSBP)), real(fVidal * conj(fVidal)); atol=1e-3)
   end
 end

test/test_belief_propagation.jl

@@ -31,31 +31,27 @@ using LinearAlgebra: eigvals, tr
 using NamedGraphs: NamedEdge, NamedGraph
 using NamedGraphs.NamedGraphGenerators: named_comb_tree, named_grid
 using NamedGraphs.PartitionedGraphs: PartitionVertex, partitionedges
-using Random: Random
 using SplitApplyCombine: group
+using StableRNGs: StableRNG
 using Test: @test, @testset
-
 @testset "belief_propagation" begin
   ITensors.disable_warn_order()
-
   g = named_grid((3, 3))
   s = siteinds("S=1/2", g)
   χ = 2
-  Random.seed!(1234)
-  ψ = random_tensornetwork(s; link_space=χ)
+  rng = StableRNG(1234)
+  ψ = random_tensornetwork(rng, s; link_space=χ)
   ψψ = ψ ⊗ prime(dag(ψ); sites=[])
-
   bpc = BeliefPropagationCache(ψψ)
   bpc = update(bpc; maxiter=50, tol=1e-10)
-
   #Test messages are converged
   for pe in partitionedges(partitioned_tensornetwork(bpc))
     @test update_message(bpc, pe) ≈ message(bpc, pe) atol = 1e-8
   end
-
   #Test updating the underlying tensornetwork in the cache
   v = first(vertices(ψψ))
-  new_tensor = random_itensor(inds(ψψ[v]))
+  rng = StableRNG(1234)
+  new_tensor = random_itensor(rng, inds(ψψ[v]))
   bpc_updated = update_factor(bpc, v, new_tensor)
   ψψ_updated = tensornetwork(bpc_updated)
   @test ψψ_updated[v] == new_tensor
@@ -79,7 +75,8 @@ using Test: @test, @testset
   #Test edge case of network which evalutes to 0
   χ = 2
   g = named_grid((3, 1))
-  ψ = random_tensornetwork(ComplexF64, g; link_space=χ)
+  rng = StableRNG(1234)
+  ψ = random_tensornetwork(rng, ComplexF64, g; link_space=χ)
   ψ[(1, 1)] = 0.0 * ψ[(1, 1)]
   @test iszero(scalar(ψ; alg="bp"))
 end

test/test_binary_tree_partition.jl

@@ -23,6 +23,7 @@ using NamedGraphs.NamedGraphGenerators: named_grid
 using NamedGraphs.GraphsExtensions:
   is_binary_arborescence, post_order_dfs_vertices, root_vertex
 using OMEinsumContractionOrders: OMEinsumContractionOrders
+using StableRNGs: StableRNG
 using Test: @test, @testset
 
 @testset "test mincut functions on top of MPS" begin
@@ -60,7 +61,8 @@ end
 @testset "test _binary_tree_partition_inds of a 2D network" begin
   N = (3, 3, 3)
   linkdim = 2
-  network = random_tensornetwork(IndsNetwork(named_grid(N)); link_space=linkdim)
+  rng = StableRNG(1234)
+  network = random_tensornetwork(rng, IndsNetwork(named_grid(N)); link_space=linkdim)
   tn = Array{ITensor,length(N)}(undef, N...)
   for v in vertices(network)
     tn[v...] = network[v...]

test/test_contract_deltas.jl

@@ -17,6 +17,7 @@ using ITensorNetworks:
   random_tensornetwork
 using NamedGraphs.GraphsExtensions: leaf_vertices, root_vertex
 using NamedGraphs.NamedGraphGenerators: named_grid
+using StableRNGs: StableRNG
 using Test: @test, @testset
 
 @testset "test _contract_deltas with no deltas" begin
@@ -41,7 +42,8 @@ end
 @testset "test _contract_deltas over partition" begin
   N = (3, 3, 3)
   linkdim = 2
-  network = random_tensornetwork(IndsNetwork(named_grid(N)); link_space=linkdim)
+  rng = StableRNG(1234)
+  network = random_tensornetwork(rng, IndsNetwork(named_grid(N)); link_space=linkdim)
   tn = Array{ITensor,length(N)}(undef, N...)
   for v in vertices(network)
     tn[v...] = network[v...]

test/test_contraction_sequence.jl

@@ -5,17 +5,16 @@ using ITensorNetworks:
 using ITensors: ITensors, contract
 using NamedGraphs.NamedGraphGenerators: named_grid
 using OMEinsumContractionOrders: OMEinsumContractionOrders
-using Random: Random
+using StableRNGs: StableRNG
 using Test: @test, @testset
-
-Random.seed!(1234)
 @testset "contraction_sequence" begin
   ITensors.@disable_warn_order begin
     dims = (2, 3)
     g = named_grid(dims)
     s = siteinds("S=1/2", g)
     χ = 10
-    ψ = random_tensornetwork(s; link_space=χ)
+    rng = StableRNG(1234)
+    ψ = random_tensornetwork(rng, s; link_space=χ)
     tn = norm_sqr_network(ψ)
     seq_optimal = contraction_sequence(tn; alg="optimal")
     res_optimal = contract(tn; sequence=seq_optimal)[]

test/test_contraction_sequence_to_graph.jl

@@ -6,20 +6,18 @@ using ITensorNetworks:
   contraction_sequence_to_graph,
   contraction_tree_leaf_bipartition,
   random_tensornetwork
-using Test: @test, @testset
 using NamedGraphs.GraphsExtensions:
   is_leaf_vertex, leaf_vertices, non_leaf_edges, root_vertex
 using NamedGraphs.NamedGraphGenerators: named_grid
-
+using StableRNGs: StableRNG
+using Test: @test, @testset
 @testset "contraction_sequence_to_graph" begin
   n = 3
   dims = (n, n)
   g = named_grid(dims)
-
-  tn = random_tensornetwork(g; link_space=2)
-
+  rng = StableRNG(1234)
+  tn = random_tensornetwork(rng, g; link_space=2)
   seq = contraction_sequence(tn)
-
   g_directed_seq = contraction_sequence_to_digraph(seq)
   g_seq_leaves = leaf_vertices(g_directed_seq)
   @test length(g_seq_leaves) == n * n

test/test_expect.jl

@@ -1,5 +1,4 @@
 @eval module $(gensym())
-
 using Graphs: SimpleGraph, uniform_tree
 using NamedGraphs: NamedGraph, vertices
 using NamedGraphs.NamedGraphGenerators: named_grid
@@ -10,18 +9,16 @@ using ITensorNetworks:
   expect,
   random_tensornetwork,
   original_state_vertex
-using Random: Random
 using SplitApplyCombine: group
+using StableRNGs: StableRNG
 using Test: @test, @testset
-
 @testset "Test Expect" begin
-  Random.seed!(1234)
-
   #Test on a tree
   L, χ = 4, 2
   g = NamedGraph(SimpleGraph(uniform_tree(L)))
   s = siteinds("S=1/2", g)
-  ψ = random_tensornetwork(s; link_space=χ)
+  rng = StableRNG(1234)
+  ψ = random_tensornetwork(rng, s; link_space=χ)
   sz_bp = expect(ψ, "Sz"; alg="bp")
   sz_exact = expect(ψ, "Sz"; alg="exact")
   @test sz_bp ≈ sz_exact
@@ -30,7 +27,8 @@ using Test: @test, @testset
   L, χ = 2, 2
   g = named_grid((L, L))
   s = siteinds("S=1/2", g)
-  ψ = random_tensornetwork(s; link_space=χ)
+  rng = StableRNG(1234)
+  ψ = random_tensornetwork(rng, s; link_space=χ)
   cache_construction_function =
     f -> BeliefPropagationCache(
       f; partitioned_vertices=group(v -> (original_state_vertex(f, v)[1]), vertices(f))