From e01e9a4e0c0198548b1d960c4b595c4272ec1cdd Mon Sep 17 00:00:00 2001 From: Iago Leal Date: Wed, 18 Dec 2024 18:02:27 -0500 Subject: [PATCH 1/4] Allow atol / rtol as stopping criteria --- src/solver.jl | 18 ++++++++++++++++++ 1 file changed, 18 insertions(+) diff --git a/src/solver.jl b/src/solver.jl index c0dd682..31295ec 100644 --- a/src/solver.jl +++ b/src/solver.jl @@ -3,6 +3,7 @@ import Combinatorics: multiset_permutations using ITensorMPS, ITensors import ITensorMPS: MPS, MPO, dmrg, OpSum, OpName, SiteType, StateName +import ITensorMPS: AbstractObserver function issquare(A :: AbstractMatrix) @@ -14,6 +15,21 @@ function maybe(f::Function, mx::Union{T, Nothing}; default=nothing) where T return isnothing(mx) ? default : f(mx) end +mutable struct ConvergenceObserver <: AbstractObserver + atol :: Float64 + rtol :: Float64 + previous_energy :: Float64 +end + +function ITensorMPS.checkdone!(o::ConvergenceObserver; energy, sweep, kwargs...) + abs_err = abs(energy - o.previous_energy) + rel_err = abs_err / abs(energy) + o.previous_energy = energy + + return abs_err < o.atol || rel_err < o.rtol +end + + # Diagonal matrix whose eigenvalues are the ordered feasible values for an integer variable. # For qubits, this is a projection on |1>. Or equivalently, (I - σ_z) / 2. # This looks like type piracy, @@ -145,9 +161,11 @@ function minimize( Q :: AbstractMatrix{T} # Initial product state # Slight entanglement to help DMRG avoid local minima psi0 = random_mps(T, sites; linkdims=2) + observer = ConvergenceObserver(atol, rtol, Inf) energy, tn = dmrg(device(H), device(psi0) ; nsweeps = iterations + , observer = observer , maxdim, cutoff, kwargs...) # The calculated energy has approximation errors compared to the true solution. From e385792b7a2949dfd60c5b2908e57b44e7810f30 Mon Sep 17 00:00:00 2001 From: Iago Leal Date: Thu, 19 Dec 2024 18:05:08 -0500 Subject: [PATCH 2/4] Use isapprox for calculating convergence --- src/solver.jl | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/src/solver.jl b/src/solver.jl index 31295ec..4734dd4 100644 --- a/src/solver.jl +++ b/src/solver.jl @@ -22,11 +22,10 @@ mutable struct ConvergenceObserver <: AbstractObserver end function ITensorMPS.checkdone!(o::ConvergenceObserver; energy, sweep, kwargs...) - abs_err = abs(energy - o.previous_energy) - rel_err = abs_err / abs(energy) + converged = isapprox(energy, o.previous_energy; atol = o.atol, rtol = o.rtol) o.previous_energy = energy - return abs_err < o.atol || rel_err < o.rtol + return converged end From 508e65d5c2596a75f8761c78ee284d2b01442e53 Mon Sep 17 00:00:00 2001 From: Iago Leal Date: Wed, 25 Dec 2024 19:26:00 -0500 Subject: [PATCH 3/4] Allow time limit stopping criterion --- src/solver.jl | 15 +++++++++++---- 1 file changed, 11 insertions(+), 4 deletions(-) diff --git a/src/solver.jl b/src/solver.jl index 4734dd4..f00c50f 100644 --- a/src/solver.jl +++ b/src/solver.jl @@ -18,14 +18,18 @@ end mutable struct ConvergenceObserver <: AbstractObserver atol :: Float64 rtol :: Float64 + time_limit :: Float64 + init_time :: Float64 previous_energy :: Float64 + + ConvergenceObserver(atol, rtol, time_limit = +Inf) = new(atol, rtol, time_limit, time() , +Inf) end function ITensorMPS.checkdone!(o::ConvergenceObserver; energy, sweep, kwargs...) converged = isapprox(energy, o.previous_energy; atol = o.atol, rtol = o.rtol) o.previous_energy = energy - return converged + return converged || time() - o.init_time > o.time_limit end @@ -145,8 +149,11 @@ minimize(Q; device = Metal.mtl) function minimize( Q :: AbstractMatrix{T} , l :: Union{AbstractVector{T}, Nothing} = nothing , c :: T = zero(T) - ; cutoff :: Float64 = 1e-8 - , iterations :: Int = 10 + ; cutoff :: Float64 = 1e-8 + , atol :: Float64 = cutoff + , rtol :: Float64 = atol + , iterations :: Int = 10 + , time_limit :: Float64 = +Inf , maxdim = [10, 20, 100, 100, 200] , device :: Function = identity , kwargs... @@ -160,7 +167,7 @@ function minimize( Q :: AbstractMatrix{T} # Initial product state # Slight entanglement to help DMRG avoid local minima psi0 = random_mps(T, sites; linkdims=2) - observer = ConvergenceObserver(atol, rtol, Inf) + observer = ConvergenceObserver(atol, rtol, time_limit) energy, tn = dmrg(device(H), device(psi0) ; nsweeps = iterations From 3532c5b61089ce01adac50013eda15069b2fa5c8 Mon Sep 17 00:00:00 2001 From: Iago Leal Date: Tue, 7 Jan 2025 16:17:14 -0700 Subject: [PATCH 4/4] Variance stopping criterion --- src/solver.jl | 41 +++++++++++++++++++++++------------------ 1 file changed, 23 insertions(+), 18 deletions(-) diff --git a/src/solver.jl b/src/solver.jl index f00c50f..63bf138 100644 --- a/src/solver.jl +++ b/src/solver.jl @@ -15,21 +15,25 @@ function maybe(f::Function, mx::Union{T, Nothing}; default=nothing) where T return isnothing(mx) ? default : f(mx) end +variance(H::MPO, x::MPS) = inner(H, x, H, x) - inner(x, H, x)^2 + mutable struct ConvergenceObserver <: AbstractObserver atol :: Float64 rtol :: Float64 + variance_tol :: Float64 time_limit :: Float64 init_time :: Float64 previous_energy :: Float64 - ConvergenceObserver(atol, rtol, time_limit = +Inf) = new(atol, rtol, time_limit, time() , +Inf) + ConvergenceObserver(atol, rtol, vtol=0.0, time_limit = +Inf) = new(atol, rtol, vtol, time_limit, time() , +Inf) end -function ITensorMPS.checkdone!(o::ConvergenceObserver; energy, sweep, kwargs...) - converged = isapprox(energy, o.previous_energy; atol = o.atol, rtol = o.rtol) +function ITensorMPS.checkdone!(o::ConvergenceObserver; energy, sweep, psi, outputlevel) + stagnated = isapprox(energy, o.previous_energy; atol = o.atol, rtol = o.rtol) o.previous_energy = energy + var = o.variance_tol > 0 ? variance(o.H, psi) : Inf - return converged || time() - o.init_time > o.time_limit + return stagnated || var < o.variance_tol || time() - o.init_time > o.time_limit end @@ -147,17 +151,19 @@ minimize(Q; device = Metal.mtl) ``` """ function minimize( Q :: AbstractMatrix{T} - , l :: Union{AbstractVector{T}, Nothing} = nothing - , c :: T = zero(T) - ; cutoff :: Float64 = 1e-8 - , atol :: Float64 = cutoff - , rtol :: Float64 = atol - , iterations :: Int = 10 - , time_limit :: Float64 = +Inf - , maxdim = [10, 20, 100, 100, 200] - , device :: Function = identity - , kwargs... - ) where {T} + , l :: Union{AbstractVector{T}, Nothing} = nothing + , c :: T = zero(T) + ; cutoff :: Float64 = 1e-8 # a cutoff of 1E-5 gives sensible accuracy; a cutoff of 1E-8 is high accuracy; and a cutoff of 1E-12 is near exact accuracy. (https://itensor.org/docs.cgi?page=tutorials/dmrg_params) + , atol :: Float64 = cutoff + , rtol :: Float64 = atol + , vtol :: Float64 = 0.0 + , iterations :: Int = 10 + , time_limit :: Float64 = +Inf + , maxdim = [10, 20, 50, 100, 100, 200] + , noise = [1e-5, 1e-6, 1e-7, 1e-8, 1e-10, 1e-12, 0.0] # 1E-5 is a lot of noise and 1E-12 is a minimal amount of noise that can still be considered non-zero. + , device :: Function = identity + , kwargs... + ) where {T} particles = size(Q)[1] # Quantization @@ -167,7 +173,7 @@ function minimize( Q :: AbstractMatrix{T} # Initial product state # Slight entanglement to help DMRG avoid local minima psi0 = random_mps(T, sites; linkdims=2) - observer = ConvergenceObserver(atol, rtol, time_limit) + observer = ConvergenceObserver(atol, rtol, vtol, time_limit) energy, tn = dmrg(device(H), device(psi0) ; nsweeps = iterations @@ -186,8 +192,7 @@ end """ minimize(Q::Matrix, c::Number; kwargs...) -Solve the Quadratic Unconstrained Binary Optimization problem -without a linear term. +Solve the Quadratic Unconstrained Binary Optimization problem with no linear term. min b'Qb + c s.t. b_i in {0, 1}