ITensor / ITensors.jl

A Julia library for efficient tensor computations and tensor network calculations
https://itensor.org
Apache License 2.0
541 stars 124 forks source link

[ITensors] [ENHANCEMENT] Expose `truncerr` in calls to `truncate!` #1516

Open NuclearPowerNerd opened 4 months ago

NuclearPowerNerd commented 4 months ago

Hello,

I'd be interested on thoughts here. I think it would be useful to be able to get the truncation error corresponding to a call to truncate!.

Currently truncate! as defined in abstractmps.jl only returns the MPS, but each truncation that is performed is a call to svd which returns the TruncSVD type which contains the truncation error.

Maybe one solution then would be to allocate a small vector with as many elements as there are bonds of the MPS. Then in the right to left sweep just collect the truncation error into each element and have that be optionally returned (could make it a kwarg).

If this functionality is already available then please just disregard and close this feature request.

Here is some psuedocode for what I had in mind (I've not actually tested this to see if it works, just what I was thinking). If this is something that sounds reasonable, I can prepare the pull request. I can make the corresponding updates to truncate as well.

function truncate!(
  ::Algorithm"frobenius", M::AbstractMPS; site_range=1:length(M), return_truncation_error=false, kwargs...
)
  N = length(M)
  truncerr = zeros(N-1)

  # Left-orthogonalize all tensors to make
  # truncations controlled
  orthogonalize!(M, last(site_range))

  # Perform truncations in a right-to-left sweep
  for (idx, j) in enumerate(reverse((first(site_range) + 1):last(site_range)))
    rinds = uniqueinds(M[j], M[j - 1])
    ltags = tags(commonind(M[j], M[j - 1]))
    U, S, V, spec = svd(M[j], rinds; lefttags=ltags, kwargs...)
    truncerr[idx] = spec.truncerr
    M[j] = U
    M[j - 1] *= (S * V)
    setrightlim!(M, j)
  end
  if return_truncation_error
  return M, truncerr 
  else
   return M
   end
end
emstoudenmire commented 4 months ago

I think this is a good feature for us to pursue, in the sense of making this information available through a return value. Probably, though, we should do it by having truncate return the MPS and a NamedTuple with various things in it, including these truncation errors. The nice thing about a NamedTuple is it would let us add other things in the future and also users could pass keyword arguments to opt in to (possibly) more expensive or technical things that could be returned but which we wouldn't want to return by default.

NuclearPowerNerd commented 4 months ago

Hey Miles - sounds good. When I originally posted this I also was thinking of truncate! in isolation. But I can see how it would be useful if you could return this on any call that eventually ends up calling truncate, such as contract or add (when using densitymatrix).

For the named tuple approach, were you thinking of something like this? If not, do you have any examples within the current ITensor code I could take a look at for reference?

function truncate!(
  ::Algorithm"frobenius", M::AbstractMPS; site_range=1:length(M), return_spec=false, kwargs...
)
  N = length(M)
  spec = (;truncerr=zeros(N-1)) # calling it spec to keep it generic? Since other things might be added in the future

  # Left-orthogonalize all tensors to make
  # truncations controlled
  orthogonalize!(M, last(site_range))

  # Perform truncations in a right-to-left sweep
  for (idx, j) in enumerate(reverse((first(site_range) + 1):last(site_range)))
    rinds = uniqueinds(M[j], M[j - 1])
    ltags = tags(commonind(M[j], M[j - 1]))
    U, S, V, spec = svd(M[j], rinds; lefttags=ltags, kwargs...)
    spec.truncerr[idx] = spec.truncerr
    M[j] = U
    M[j - 1] *= (S * V)
    setrightlim!(M, j)
  end
  if return_spec
  return M, spec
  else
   return M
   end
end