From 3de5685e39d4b8b115b841d425bf4e26905fb337 Mon Sep 17 00:00:00 2001 From: Paul Brehmer Date: Thu, 11 Jul 2024 11:21:51 +0200 Subject: [PATCH] Add :FixedIter and :LeftMoves safeguard, pdate tests --- src/PEPSKit.jl | 2 +- src/algorithms/peps_opt.jl | 66 ++++++++++++++++++++++++++++++-------- test/ctmrg/fixedsvd.jl | 63 +++++++++++++++++++++++++----------- test/heisenberg.jl | 5 +-- 4 files changed, 101 insertions(+), 35 deletions(-) diff --git a/src/PEPSKit.jl b/src/PEPSKit.jl index 1e1ec392..07ee4f40 100644 --- a/src/PEPSKit.jl +++ b/src/PEPSKit.jl @@ -61,7 +61,7 @@ module Defaults const fpgrad_tol = 1e-6 end -export SVDAdjoint, IterSVD, FixedSVD, NonTruncSVDAdjoint +export SVDAdjoint, IterSVD, NonTruncSVDAdjoint export FixedSpaceTruncation, ProjectorAlg, CTMRG, CTMRGEnv, correlation_length export LocalOperator export expectation_value, costfun diff --git a/src/algorithms/peps_opt.jl b/src/algorithms/peps_opt.jl index c6d724ae..507e0924 100644 --- a/src/algorithms/peps_opt.jl +++ b/src/algorithms/peps_opt.jl @@ -2,9 +2,15 @@ abstract type GradMode{F} end """ struct GeomSum(; maxiter=Defaults.fpgrad_maxiter, tol=Defaults.fpgrad_tol, - verbosity=0, iterscheme=:FixedIter) <: GradMode + verbosity=0, iterscheme=:FixedIter) <: GradMode{iterscheme} Gradient mode for CTMRG using explicit evaluation of the geometric sum. + +With `iterscheme` the style of CTMRG iteration which is being differentiated can be chosen. +If set to `:FixedIter`, the differentiated CTMRG iteration is assumed to have a pre-computed +SVD of the environments with a fixed set of gauges. Alternatively, if set to `:DiffGauge`, +the differentiated iteration consists of a CTMRG iteration and a subsequent gauge fixing step, +such that `gauge_fix` will also be differentiated everytime a CTMRG derivative is computed. """ struct GeomSum{F} <: GradMode{F} maxiter::Int @@ -22,9 +28,15 @@ end """ struct ManualIter(; maxiter=Defaults.fpgrad_maxiter, tol=Defaults.fpgrad_tol, - verbosity=0, iterscheme=:FixedIter) <: GradMode + verbosity=0, iterscheme=:FixedIter) <: GradMode{iterscheme} Gradient mode for CTMRG using manual iteration to solve the linear problem. + +With `iterscheme` the style of CTMRG iteration which is being differentiated can be chosen. +If set to `:FixedIter`, the differentiated CTMRG iteration is assumed to have a pre-computed +SVD of the environments with a fixed set of gauges. Alternatively, if set to `:DiffGauge`, +the differentiated iteration consists of a CTMRG iteration and a subsequent gauge fixing step, +such that `gauge_fix` will also be differentiated everytime a CTMRG derivative is computed. """ struct ManualIter{F} <: GradMode{F} maxiter::Int @@ -41,10 +53,16 @@ function ManualIter(; end """ - struct LinSolver(; solver=KrylovKit.GMRES(), iterscheme=:FixedIter) <: GradMode + struct LinSolver(; solver=KrylovKit.GMRES(), iterscheme=:FixedIter) <: GradMode{iterscheme} Gradient mode wrapper around `KrylovKit.LinearSolver` for solving the gradient linear problem using iterative solvers. + +With `iterscheme` the style of CTMRG iteration which is being differentiated can be chosen. +If set to `:FixedIter`, the differentiated CTMRG iteration is assumed to have a pre-computed +SVD of the environments with a fixed set of gauges. Alternatively, if set to `:DiffGauge`, +the differentiated iteration consists of a CTMRG iteration and a subsequent gauge fixing step, +such that `gauge_fix` will also be differentiated everytime a CTMRG derivative is computed. """ struct LinSolver{F} <: GradMode{F} solver::KrylovKit.LinearSolver @@ -57,8 +75,8 @@ function LinSolver(; end """ - PEPSOptimize{G}(; boundary_alg = CTMRG(), optimizer::OptimKit.OptimizationAlgorithm = LBFGS() - reuse_env::Bool = true, gradient_alg::G, verbosity::Int = 0) + PEPSOptimize{G}(; boundary_alg=CTMRG(), optimizer::OptimKit.OptimizationAlgorithm=LBFGS() + reuse_env::Bool=true, gradient_alg::G=LinSolver(), verbosity::Int=0) Algorithm struct that represent PEPS ground-state optimization using AD. Set the algorithm to contract the infinite PEPS in `boundary_alg`; @@ -69,14 +87,36 @@ step by setting `reuse_env` to true. Otherwise a random environment is used at e step. The CTMRG gradient itself is computed using the `gradient_alg` algorithm. Different levels of output verbosity can be activated using `verbosity` (0, 1 or 2). """ -@kwdef struct PEPSOptimize{G} - boundary_alg::CTMRG = CTMRG() # Algorithm to find boundary environment - optimizer::OptimKit.OptimizationAlgorithm = LBFGS( - 4; maxiter=100, gradtol=1e-4, verbosity=2 - ) - reuse_env::Bool = true # Reuse environment of previous optimization as initial guess for next - gradient_alg::G = GeomSum() # Algorithm to solve gradient linear problem - verbosity::Int = 0 +struct PEPSOptimize{G} + boundary_alg::CTMRG + optimizer::OptimKit.OptimizationAlgorithm + reuse_env::Bool + gradient_alg::G + verbosity::Int + + function PEPSOptimize( # Inner constructor to prohibit illegal setting combinations + boundary_alg::CTMRG{S}, + optimizer, + reuse_env, + gradient_alg::G, + verbosity, + ) where {S,G} + if gradient_alg isa GradMode + if S == :LeftMoves && G.parameters[1] == :FixedIter + throw(ArgumentError(":LeftMoves and :FixedIter are not compatible")) + end + end + return new{G}(boundary_alg, optimizer, reuse_env, gradient_alg, verbosity) + end +end +function PEPSOptimize(; + boundary_alg=CTMRG(), + optimizer=LBFGS(4; maxiter=100, gradtol=1e-4, verbosity=2), + reuse_env=true, + gradient_alg=LinSolver(), + verbosity=0, +) + return PEPSOptimize(boundary_alg, optimizer, reuse_env, gradient_alg, verbosity) end """ diff --git a/test/ctmrg/fixedsvd.jl b/test/ctmrg/fixedsvd.jl index f038b8e2..9b0ec4d1 100644 --- a/test/ctmrg/fixedsvd.jl +++ b/test/ctmrg/fixedsvd.jl @@ -3,6 +3,7 @@ using Random using TensorKit using PEPSKit using PEPSKit: + FixedSVD, ctmrg_iter, gauge_fix, fix_relative_phases, @@ -14,22 +15,46 @@ using PEPSKit: χenv = 16 ctm_alg = CTMRG(; tol=1e-12, miniter=4, maxiter=100, verbosity=1, ctmrgscheme=:AllSides) -# initialize states -Random.seed!(91283219347) -psi = InfinitePEPS(2, χbond) -env_conv1 = leading_boundary(CTMRGEnv(psi; Venv=ComplexSpace(χenv)), psi, ctm_alg) - -# do extra iteration to get SVD -env_conv2, info = ctmrg_iter(psi, env_conv1, ctm_alg) -env_fix, signs = gauge_fix(env_conv1, env_conv2) -@test check_elementwise_convergence(env_conv1, env_fix) - -# fix gauge of SVD -U_fix, V_fix = fix_relative_phases(info.U, info.V, signs) -svd_alg_fix = SVDAdjoint(; fwd_alg=FixedSVD(U_fix, info.S, V_fix)) -ctm_alg_fix = CTMRG(; svd_alg=svd_alg_fix, trscheme=notrunc(), ctmrgscheme=:AllSides) - -# do iteration with FixedSVD -env_fixedsvd, = ctmrg_iter(psi, env_conv1, ctm_alg_fix) -env_fixedsvd = fix_global_phases(env_conv1, env_fixedsvd) -@test check_elementwise_convergence(env_conv1, env_fixedsvd) +@testset "(1, 1) unit cell" begin + # initialize states + Random.seed!(91283219347) + psi = InfinitePEPS(2, χbond) + env_conv1 = leading_boundary(CTMRGEnv(psi; Venv=ComplexSpace(χenv)), psi, ctm_alg) + + # do extra iteration to get SVD + env_conv2, info = ctmrg_iter(psi, env_conv1, ctm_alg) + env_fix, signs = gauge_fix(env_conv1, env_conv2) + @test check_elementwise_convergence(env_conv1, env_fix) + + # fix gauge of SVD + U_fix, V_fix = fix_relative_phases(info.U, info.V, signs) + svd_alg_fix = SVDAdjoint(; fwd_alg=FixedSVD(U_fix, info.S, V_fix)) + ctm_alg_fix = CTMRG(; svd_alg=svd_alg_fix, trscheme=notrunc(), ctmrgscheme=:AllSides) + + # do iteration with FixedSVD + env_fixedsvd, = ctmrg_iter(psi, env_conv1, ctm_alg_fix) + env_fixedsvd = fix_global_phases(env_conv1, env_fixedsvd) + @test check_elementwise_convergence(env_conv1, env_fixedsvd) +end + +@testset "(3, 4) unit cell" begin + # initialize states + Random.seed!(91283219347) + psi = InfinitePEPS(2, χbond; unitcell=(3, 4)) + env_conv1 = leading_boundary(CTMRGEnv(psi; Venv=ComplexSpace(χenv)), psi, ctm_alg) + + # do extra iteration to get SVD + env_conv2, info = ctmrg_iter(psi, env_conv1, ctm_alg) + env_fix, signs = gauge_fix(env_conv1, env_conv2) + @test check_elementwise_convergence(env_conv1, env_fix) + + # fix gauge of SVD + U_fix, V_fix = fix_relative_phases(info.U, info.V, signs) + svd_alg_fix = SVDAdjoint(; fwd_alg=FixedSVD(U_fix, info.S, V_fix)) + ctm_alg_fix = CTMRG(; svd_alg=svd_alg_fix, trscheme=notrunc(), ctmrgscheme=:AllSides) + + # do iteration with FixedSVD + env_fixedsvd, = ctmrg_iter(psi, env_conv1, ctm_alg_fix) + env_fixedsvd = fix_global_phases(env_conv1, env_fixedsvd) + @test check_elementwise_convergence(env_conv1, env_fixedsvd) +end diff --git a/test/heisenberg.jl b/test/heisenberg.jl index 54c40972..352caaa9 100644 --- a/test/heisenberg.jl +++ b/test/heisenberg.jl @@ -14,12 +14,13 @@ ctm_alg = CTMRG(; maxiter=100, verbosity=1, trscheme=truncdim(χenv), - svd_alg=SVDAdjoint(; fwd_alg=TensorKit.SVD(), rrule_alg=GMRES(; tol=1e-12)), + svd_alg=SVDAdjoint(; fwd_alg=TensorKit.SVD(), rrule_alg=GMRES(; tol=1e-10)), + ctmrgscheme=:AllSides ) opt_alg = PEPSOptimize(; boundary_alg=ctm_alg, optimizer=LBFGS(4; maxiter=100, gradtol=1e-3, verbosity=2), - gradient_alg=LinSolver(; solver=GMRES(; tol=1e-6, maxiter=100)), + gradient_alg=LinSolver(; solver=GMRES(; tol=1e-6, maxiter=100), iterscheme=:FixedIter), reuse_env=true, verbosity=2, )