Refactor keyword argument processing, remove applyexp

Project.toml

@@ -1,7 +1,7 @@
 name = "ITensorTDVP"
 uuid = "25707e16-a4db-4a07-99d9-4d67b7af0342"
 authors = ["Matthew Fishman <mfishman@flatironinstitute.org> and contributors"]
-version = "0.1.4"
+version = "0.2"
 
 [deps]
 ITensors = "9136182c-28ba-11e9-034c-db9fb085ebd5"

README.md

@@ -6,10 +6,7 @@
 [![Coverage](https://codecov.io/gh/mtfishman/ITensorTDVP.jl/branch/main/graph/badge.svg)](https://codecov.io/gh/mtfishman/ITensorTDVP.jl)
 [![Code Style: Blue](https://img.shields.io/badge/code%20style-blue-4495d1.svg)](https://github.com/invenia/BlueStyle)
 
-| :information_source: NOTE          |
-|:---------------------------|
-| This package is effectively a generalization of the DMRG code in [ITensors.jl](https://github.com/ITensor/ITensors.jl), using the MPS/MPO types from that package. It provides a general MPS "solver" interface which allows us to implement other MPS/MPO optimization/solver functionality like DMRG (`ITensorTDVP.dmrg` or `eigsolve`), TDVP (`tdvp` or `exponentiate`), linear solving (`linsolve`), DMRG-X (`dmrg_x`), etc. while sharing most of the code across those different functions. Therefore, it effectively supercedes the DMRG functionality in ITensors.jl (`dmrg`), and provides its own `ITensorTDVP.dmrg`/`eigsolve` function that is essentially the same as the `dmrg` function from ITensors.jl. This package is fairly stable and appropriate for general use. The primary missing feature is a lack of modern subspace expansion tools for methods like TDVP. However, 2-site TDVP or TEBD is often sufficient for performing subspace expansion (except when [it's not](https://arxiv.org/abs/2005.06104)). |
-| However, note that future developments, including modern subspace expansion tools, are being developed in our next-generation tensor network library [ITensorNetworks.jl](https://github.com/mtfishman/ITensorNetworks.jl). That package has contraction, optimization, and evolution tools for general tensor networks, as well as methods like DMRG, TDVP, and linear solving for tree tensor networks, which effectively supercedes the functionality of this package which is limited to solvers for just MPS/MPO (linear/path graph) tensor networks. ITensorNetworks.jl is under heavy development and is _not_ meant for general usage at the moment, except for those who are brave enough to handle missing features and breaking interfaces. Basically, for the average user who wants stable and reliable code, if you need to use TDVP or linear solving on MPS, you should use this package, but keep in mind that it will be superceded by ITensorNetworks.jl eventually, once that code catches up on functionality for TTN that we currently have for MPS/MPO in ITensors.jl (like Hamiltonian construction from `OpSum`, general TTN algebra like addition and contraction, etc.). If you need to use functionality like DMRG, TDVP, or linear solving for TTN, you can try out ITensorNetworks.jl, but keep in mind that it is still a work in progress. |
+## Installation
 
 To install this package, you can use the following steps:
 ```
@@ -19,3 +16,22 @@ julia> ]
 
 pkg> add ITensorTDVP
 ```
+
+## News
+
+### ITensorTDVP.jl v0.2 Release Notes
+
+#### Breaking changes
+
+- The `applyexp` Krylov exponentiation solver backend has been removed, and `solver_backend="applyexp"` option for `tdvp` now just calls `exponentiate` from KrylovKit.jl. They were in many ways the same but `exponentiate` has more advanced features like restarts. Using `solver_backend="applyexp"` now prints a warning to that effect. To remove the warning, don't specify the `solver_backend` keyword argument. Alternatively, you can set it to `"exponentiate"`, though specifying it is now superfluous since there is now only one `solver_backend` available and the `solver_backend` keyword argument may be removed in future minor releases.
+
+#### Bug fixes
+
+- `svd_alg` now doesn't specify a default value, so the default value is set by the `svd` function in ITensors/NDTensors. This fixes an issue using ITensorTDVP.jl and GPU backends, where the default value being set in ITensorTDVP.jl wasn't compatible with the options available in some GPU backends like CUDA.
+- More generally, keyword arguments are handled better throughout the package, so default values are handled more systematically and keyword arguments are listed or forwarded more explicitly, so it should catch more mistakes like passing an incorrect keyword argument name.
+
+## About
+
+This package is effectively a generalization of the DMRG code in [ITensors.jl](https://github.com/ITensor/ITensors.jl), using the MPS/MPO types from that package. It provides a general MPS "solver" interface which allows us to implement other MPS/MPO optimization/solver functionality like DMRG (`ITensorTDVP.dmrg` or `eigsolve`), TDVP (`tdvp` or `exponentiate`), linear solving (`linsolve`), DMRG-X (`dmrg_x`), etc. while sharing most of the code across those different functions. Therefore, it effectively supercedes the DMRG functionality in ITensors.jl (`dmrg`), and provides its own `ITensorTDVP.dmrg`/`eigsolve` function that is essentially the same as the `dmrg` function from ITensors.jl. This package is fairly stable and appropriate for general use. The primary missing feature is a lack of modern subspace expansion tools for methods like TDVP. However, 2-site TDVP or TEBD is often sufficient for performing subspace expansion (except when [it's not](https://arxiv.org/abs/2005.06104)).
+
+However, note that future developments, including modern subspace expansion tools, are being developed in our next-generation tensor network library [ITensorNetworks.jl](https://github.com/mtfishman/ITensorNetworks.jl). The goal of that package is to provide contraction, optimization, and evolution tools for general tensor networks, as well as methods like DMRG, TDVP, and linear solving for tree tensor networks, and the eventual goal is to replace this package which is limited to solvers for just MPS/MPO (linear/path graph) tensor networks. However, ITensorNetworks.jl is under heavy development and is _not_ meant for general usage at the moment, except for those who are brave enough to handle missing features and breaking interfaces. Basically, for the average user who wants stable and reliable code, if you need to use MPS-based TDVP or linear solving, you should use this package for the time being.

src/ITensorTDVP.jl

@@ -17,15 +17,9 @@ using ITensors:
   ProjMPS,
   set_nsite!
 
-# Compatibility of ITensor observer and Observers
+include("defaults.jl")
 include("update_observer.jl")
-
-# Utilities for making it easier
-# to define solvers (like ODE solvers)
-# for TDVP
 include("solver_utils.jl")
-
-include("applyexp.jl")
 include("tdvporder.jl")
 include("tdvpinfo.jl")
 include("tdvp_step.jl")

src/defaults.jl

@@ -0,0 +1,47 @@
+default_nsweeps() = nothing
+default_checkdone() = nothing
+default_write_when_maxdim_exceeds() = nothing
+default_nsite() = 2
+default_reverse_step() = true
+default_time_start() = 0
+default_time_step(t) = t
+default_order() = 2
+default_observer!() = NoObserver()
+default_step_observer!() = NoObserver()
+default_outputlevel() = 0
+default_normalize() = false
+default_sweep() = 1
+default_current_time() = 0
+
+# Truncation
+default_maxdim() = typemax(Int)
+default_mindim() = 1
+default_cutoff(type::Type{<:Number}) = eps(real(type))
+default_noise() = 0
+
+# Solvers
+default_tdvp_solver_backend() = "exponentiate"
+default_ishermitian() = true
+default_issymmetric() = true
+default_solver_verbosity() = 0
+
+# Customizable based on solver function
+default_solver_outputlevel(::Function) = 0
+
+default_solver_tol(::Function) = error("Not implemented")
+default_solver_which_eigenvalue(::Function) = error("Not implemented")
+default_solver_krylovdim(::Function) = error("Not implemented")
+default_solver_verbosity(::Function) = error("Not implemented")
+
+## Solver-specific keyword argument defaults
+
+# dmrg/eigsolve
+default_solver_tol(::typeof(eigsolve)) = 1e-14
+default_solver_krylovdim(::typeof(eigsolve)) = 3
+default_solver_maxiter(::typeof(eigsolve)) = 1
+default_solver_which_eigenvalue(::typeof(eigsolve)) = :SR
+
+# tdvp/exponentiate
+default_solver_tol(::typeof(exponentiate)) = 1e-12
+default_solver_krylovdim(::typeof(exponentiate)) = 30
+default_solver_maxiter(::typeof(exponentiate)) = 100

src/dmrg.jl

@@ -1,24 +1,60 @@
-function eigsolve_solver(; kwargs...)
-  function solver(H, t, psi0; kws...)
+function dmrg_solver(
+  f::typeof(eigsolve);
+  solver_which_eigenvalue,
+  ishermitian,
+  solver_tol,
+  solver_krylovdim,
+  solver_maxiter,
+  solver_verbosity,
+)
+  function solver(H, t, psi0; current_time, outputlevel)
     howmany = 1
-    which = get(kwargs, :solver_which_eigenvalue, :SR)
-    solver_kwargs = (;
-      ishermitian=get(kwargs, :ishermitian, true),
-      tol=get(kwargs, :solver_tol, 1E-14),
-      krylovdim=get(kwargs, :solver_krylovdim, 3),
-      maxiter=get(kwargs, :solver_maxiter, 1),
-      verbosity=get(kwargs, :solver_verbosity, 0),
+    which = solver_which_eigenvalue
+    vals, vecs, info = f(
+      H,
+      psi0,
+      howmany,
+      which;
+      ishermitian=default_ishermitian(),
+      tol=solver_tol,
+      krylovdim=solver_krylovdim,
+      maxiter=solver_maxiter,
+      verbosity=solver_verbosity,
     )
-    vals, vecs, info = eigsolve(H, psi0, howmany, which; solver_kwargs...)
     psi = vecs[1]
     return psi, info
   end
   return solver
 end
 
-function dmrg(H, psi0::MPS; kwargs...)
+function dmrg(
+  H,
+  psi0::MPS;
+  ishermitian=default_ishermitian(),
+  solver_which_eigenvalue=default_solver_which_eigenvalue(eigsolve),
+  solver_tol=default_solver_tol(eigsolve),
+  solver_krylovdim=default_solver_krylovdim(eigsolve),
+  solver_maxiter=default_solver_maxiter(eigsolve),
+  solver_verbosity=default_solver_verbosity(),
+  kwargs...,
+)
   t = Inf # DMRG is TDVP with an infinite timestep and no reverse step
   reverse_step = false
-  psi = tdvp(eigsolve_solver(; kwargs...), H, t, psi0; reverse_step, kwargs...)
+  psi = tdvp(
+    dmrg_solver(
+      eigsolve;
+      solver_which_eigenvalue,
+      ishermitian,
+      solver_tol,
+      solver_krylovdim,
+      solver_maxiter,
+      solver_verbosity,
+    ),
+    H,
+    t,
+    psi0;
+    reverse_step,
+    kwargs...,
+  )
   return psi
 end

src/linsolve.jl

@@ -1,4 +1,3 @@
-
 """
 Compute a solution x to the linear system:
 

src/projmps2.jl

@@ -118,7 +118,6 @@ function proj_mps(P::ProjMPS2)
   # Apply the map
   m = ITensor(true)
   for it in itensor_map
-    #@show inds(it)
     m *= it
   end
   return m