correct droplet

benchmarks/tests_grid_from_file.jl

@@ -31,7 +31,7 @@ println("examples = ", examples)
 
 
 β = 3.
-δ = 0.
+δ = 1e-8
 
 number_of_states = 10
 

src/peps_no_types.jl

@@ -173,15 +173,12 @@ structure of the partial solution
 mutable struct Partial_sol{T <: Real}
     spins::Vector{Int}
     objective::T
-    boundary::Vector{Int}
+
     function(::Type{Partial_sol{T}})(spins::Vector{Int}, objective::T) where T <:Real
-        new{T}(spins, objective, Int[])
-    end
-    function(::Type{Partial_sol{T}})(spins::Vector{Int}, objective::T, boundary::Vector{Int}) where T <:Real
-        new{T}(spins, objective, boundary)
+        new{T}(spins, objective)
     end
     function(::Type{Partial_sol{T}})() where T <:Real
-        new{T}(Int[], 1., Int[])
+        new{T}(Int[], 1.)
     end
 end
 
@@ -196,10 +193,6 @@ function update_partial_solution(ps::Partial_sol{T}, s::Int, objective::T) where
     Partial_sol{T}(vcat(ps.spins, [s]), objective)
 end
 
-function update_partial_solution(ps::Partial_sol{T}, s::Int, objective::T, boundary::Vector{Int}) where T <: Real
-    Partial_sol{T}(vcat(ps.spins, [s]), objective, ps.spins[boundary])
-end
-
 """
     spin_indices_from_above(gg::MetaGraph, ps::Partial_sol, j::Int)
 
@@ -326,18 +319,18 @@ function solve(g::MetaGraph, no_sols::Int = 2; node_size::Tuple{Int, Int} = (1,1
 
         for j in grid[row,:]
 
-            ie = indices_on_boundary(grid, j)
+            dX = indices_on_boundary(grid, j)
 
             partial_s_temp = Partial_sol{T}[]
 
-            partial_s = merge_boundaries(partial_s, δ)
+            partial_s = merge_boundaries(partial_s, dX, δ)
             for ps in partial_s
 
                 objectives = conditional_probabs(gg, ps, j, lower_mps, vec_of_T)
 
                 for l in 1:length(objectives)
                     new_objectives = ps.objective*objectives[l]
-                    ps1 = update_partial_solution(ps, l, new_objectives, ie)
+                    ps1 = update_partial_solution(ps, l, new_objectives)
                     push!(partial_s_temp, ps1)
                 end
 
@@ -385,10 +378,11 @@ function return_solutions(partial_s::Vector{Partial_sol{T}}, ns:: MetaGraph)  wh
     return spins, objective
 end
 
-function merge_boundaries(partial_s::Vector{Partial_sol{T}}, δ) where T <:Real
-    if (length(partial_s) > 1) || δ == .0
+function merge_boundaries(partial_s::Vector{Partial_sol{T}}, dX::Vector{Int}, δ) where T <:Real
+    if (length(partial_s) > 1) & (δ != .0)
         leave = [true for _ in partial_s]
-        boundaries = [ps.spins[ps.boundary] for ps in partial_s]
+
+        boundaries = [ps.spins[dX] for ps in partial_s]
 
         unique_bondaries = unique(boundaries)
         if boundaries != unique_bondaries
@@ -401,6 +395,7 @@ function merge_boundaries(partial_s::Vector{Partial_sol{T}}, δ) where T <:Real
                     objectives = objectives./maximum(objectives)
                     for ind in i[objectives .< δ]
                         leave[ind] = false
+                        println(leave)
                     end
                 end
             end

test/peps_tests.jl

@@ -15,16 +15,15 @@ if true
         @test ps.spins == []
         @test ps.objective == 1.
 
-        ps1 = Partial_sol{Float64}([1,1], 1., [1])
+        ps1 = Partial_sol{Float64}([1,1], 1.)
         @test ps1.spins == [1,1]
         @test ps1.objective == 1.
 
-        ps2 = update_partial_solution(ps1, 2, 1., [2])
+        ps2 = update_partial_solution(ps1, 2, 1.)
         @test ps2.spins == [1,1,2]
         @test ps2.objective == 1.
-        @test ps2.boundary == [1.]
 
-        ps3 = Partial_sol{Float64}([1,1,1], .2, [1,2])
+        ps3 = Partial_sol{Float64}([1,1,1], .2)
 
         b = select_best_solutions([ps3, ps2], 1)
         @test b[1].spins == [1, 1, 2]
@@ -33,8 +32,8 @@ if true
 
     @testset "functions of graph" begin
 
-        a = Partial_sol{Float64}([1,1,1,2], 0.2, [1,2,3])
-        b = Partial_sol{Float64}([1,1,2,2], 1., [1,2,3])
+        a = Partial_sol{Float64}([1,1,1,2], 0.2)
+        b = Partial_sol{Float64}([1,1,2,2], 1.)
 
         M = [1. 1. 1. 0.; 1. 1. 0. 1.; 1. 0. 1. 1.; 0. 1. 1. 1.]
 
@@ -50,14 +49,14 @@ if true
         g = M2graph(M)
         gg = graph4peps(g, (2,2))
 
-        ps = Partial_sol{Float64}([6], .2, Int[])
+        ps = Partial_sol{Float64}([6], .2)
         ul,ur = spin_indices_from_above(gg, ps, 2)
         l = spin_index_from_left(gg, ps, 2)
         @test ul == [2]
         @test ur == [1]
         @test l == 2
 
-        ps = Partial_sol{Float64}([4,6], .2, [1])
+        ps = Partial_sol{Float64}([4,6], .2)
         ul,ur = spin_indices_from_above(gg, ps, 3)
         l = spin_index_from_left(gg, ps, 3)
         @test ul == Int[]
@@ -77,29 +76,29 @@ if true
         @test i == [3, 4, 5, 6]
 
         boundary = [2,3]
-        a = Partial_sol{Float64}([1,1,1], 0.2, boundary)
-        b = Partial_sol{Float64}([2,1,1], 0.18, boundary)
-        c = Partial_sol{Float64}([1,1,2], 1., boundary)
-        d = Partial_sol{Float64}([2,1,2], .1, boundary)
+        a = Partial_sol{Float64}([1,1,1], 0.2)
+        b = Partial_sol{Float64}([2,1,1], 0.18)
+        c = Partial_sol{Float64}([1,1,2], 1.)
+        d = Partial_sol{Float64}([2,1,2], .1)
 
         vps = [a,b,c,d]
 
         thershold = 0.15
         # 0.18/0.2 = 0.9
         # 0.1/1. = 0.1
 
-        ps1 = merge_boundaries(vps, thershold)
+        ps1 = merge_boundaries(vps, boundary, thershold)
         println(ps1 == [a,b,c])
 
         thershold = 0.95
 
-        ps1 = merge_boundaries(vps, thershold)
-        println(ps1)
+        ps1 = merge_boundaries(vps, boundary, thershold)
+
         println(ps1 == [a,c])
 
         thershold = 0.
 
-        ps1 = merge_boundaries(vps, thershold)
+        ps1 = merge_boundaries(vps, boundary, thershold)
         println(ps1 == [a,b,c,d])
     end
 end
@@ -279,7 +278,7 @@ Mq[8,9] = Mq[9,8] = -0.05
     println(size(mpo1[3]))
 
     pp = PEPSRow(peps, 2)
-    println(pp)
+    #println(pp)
 
     mpo2 = MPO(PEPSRow(peps, 2))
 
@@ -312,7 +311,7 @@ Mq[8,9] = Mq[9,8] = -0.05
     B = generate_tensor(peps, (1,1))
     println(size(B))
     println(size(vec(cc)))
-    println(sum(cc) ≈ sum(B))
+    @test sum(cc) ≈ sum(B)
 
 
 end
@@ -367,15 +366,15 @@ end
     #println(size(AA[1]))
 
     # marginal prob
-    sol = Partial_sol{Float64}(Int[], 0., Int[])
+    sol = Partial_sol{Float64}(Int[], 0.)
     j = 1
     objective = conditional_probabs(gg, sol, j, lower_mps, A)
     println(objective)
 
     #objective1 = conditional_probabs(gg, sol, j, l_mps, AA)
     #println(objective1)
 
-    sol = Partial_sol{Float64}([1], objective[1], Int[])
+    sol = Partial_sol{Float64}([1], objective[1])
 
     p1 = sum(cc[1,:,:,:,:,:,:,:,:])/su
     p2 = sum(cc[2,:,:,:,:,:,:,:,:])/su
@@ -396,7 +395,7 @@ end
     A =  M[2,:]
 
     # objective value from the previous step is set artificially
-    sol = Partial_sol{Float64}(Int[1,1,1,1], 1., [1,2,3])
+    sol = Partial_sol{Float64}(Int[1,1,1,1], 1.)
     objective = conditional_probabs(gg, sol, j, lower_mps, A)
     #conditional prob
     p1 = sum(cc[1,1,1,1,1,:,:,:,:])/sum(cc[1,1,1,1,:,:,:,:,:])
@@ -408,7 +407,7 @@ end
 end
 
 @testset "test an exemple instance" begin
-    δ = 0.
+    δ = 1e-6
     g = make_interactions_case2()
     spins, objective = solve(g, 10; β = 3., χ = 2, threshold = 1e-11, δ = δ)
     @test spins[1] == [1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1]
@@ -428,7 +427,7 @@ end
 
 
 @testset "test an exemple instance on Float32" begin
-    δ = 0.
+    δ = 1e-6
     g = make_interactions_case2()
     T = Float32
     spins, objective = solve(g, 10; β = T(3.), χ = 2, threshold = 1e-11, δ = δ)