TensorNetworkQuantumSimulator.jl/test/test_constructors.jl at main · Algorithmiq/TensorNetworkQuantumSimulator.jl · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
@eval module $(gensym())
using Dictionaries: Dictionary
using ITensors: ITensors, Index, dag, inds, prime
using Random
using TensorNetworkQuantumSimulator
using Test: @testset, @test, @test_throws


@testset "Test Constructors" begin
    Random.seed!(123)

    #TensorNetwork construction from tensors
    i, j, k, l = Index(2), Index(2), Index(2), Index(2)
    A, B, C, D = ITensors.random_itensor(i, j), ITensors.random_itensor(j, k), ITensors.random_itensor(k, l), ITensors.random_itensor(l, i)
    t = TensorNetwork([A, B, C, D])
    @test t isa TensorNetwork
    @test scalartype(t) == eltype(A)
    @test maxvirtualdim(t) == 2
    @test graph(t) isa NamedGraph
    @test graph(t) == add_edge(named_path_graph(4), 1 => 4)

    #TensorNetwork pre-defined constructor
    g = named_hexagonal_lattice_graph(3, 3)
    χ = 3
    for eltype in [Float32, Float64, ComplexF32, ComplexF64]
        ψ = random_tensornetwork(eltype, g; bond_dimension = χ)
        @test ψ isa TensorNetwork
        @test scalartype(ψ) == eltype
        @test graph(ψ) == g
        @test maxvirtualdim(ψ) == 3
        @test all([length(inds(ψ[v])) == degree(g, v) for v in vertices(ψ)])

        ψdag = map_virtualinds(prime, map_tensors(dag, ψ))
        @test ψdag isa TensorNetwork
        @test ITensors.contract(ψdag; alg = "exact") ≈ conj(ITensors.contract(ψ; alg = "exact"))

        v = first(vertices(g))
        rem_vertex!(ψ, v)
        @test graph(ψ) == rem_vertex(g, v)
        @test !all([length(inds(ψ[v])) == degree(g, v) for v in vertices(ψ)])
    end

    #SiteInds
    s = siteinds("S=1/2", g)
    @test s isa Dictionary
    @test keys(s) == vertices(g)
    @test all([s[v] isa Vector{<:Index} for v in vertices(g)])
    @test all([length(s[v]) == 1 for v in vertices(g)])

    #TensorNetworkState
    χ = 3
    for eltype in [Float32, Float64, ComplexF32, ComplexF64]
        ψ = random_tensornetworkstate(eltype, g, s; bond_dimension = χ)
        @test ψ isa TensorNetworkState
        @test scalartype(ψ) == eltype
        @test siteinds(ψ) == s
        @test graph(ψ) == g
        @test maxvirtualdim(ψ) == 3
        @test all([length(inds(ψ[v])) == degree(g, v) + 1 for v in vertices(ψ)])
        @test all([siteinds(ψ, v) == s[v] for v in vertices(ψ)])

        ψ = tensornetworkstate(eltype, v -> "X+", g, "S=1/2")
        @test maxvirtualdim(ψ) == 1
        @test scalartype(ψ) == eltype
        @test all([length(inds(ψ[v])) == degree(g, v) + 1 for v in vertices(ψ)])
    end

    #Test GHZ state constructor
    ψ1, ψ2 = tensornetworkstate(Float64, v -> "↑", g, s), tensornetworkstate(Float64, v -> "↓", g, s)
    ψGHZ = ψ1 + ψ2
    @test ψGHZ isa TensorNetworkState
    @test maxvirtualdim(ψGHZ) == 2
    v, vn = first(vertices(g)), first(neighbors(g, first(vertices(g))))
    @test von_neumann_entanglement_entropy(ψGHZ, first(edges(ψGHZ)); alg = "bp") ≈ log(2)

    #Test identity state constructor
    s = siteinds("S=1/2", g; inds_per_site = 2)
    I = identity_tensornetworkstate(Float64, g, s)
    @test maxvirtualdim(I) == 1
    @test all([length(siteinds(I, v)) == 2 for v in vertices(g)])

    @test_throws ErrorException identity_tensornetworkstate(Float64, g, siteinds("S=1/2", g; inds_per_site = 3))

end

end