Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

[ITensors] Port VUMPS functions #1544

Open
wants to merge 4 commits into
base: main
Choose a base branch
from
Open

Conversation

ryanlevy
Copy link
Contributor

This PR adds a few of the functions living in https://github.com/ITensor/ITensorInfiniteMPS.jl/blob/main/src/ITensors.jl

The one thing missing at the moment is the more efficient sqrt(T::DiagTensor) which I couldn't get working, but the eigen implementation seems sane at the moment so there may be no need

@ryanlevy
Copy link
Contributor Author

Looks like there's some assumption being made that for a sliced mps e.g. psi[2:4] returns a vector of itensors instead of a MPS, not sure if we want to adjust that here or in ITensorInfiniteMPS

@ryanlevy
Copy link
Contributor Author

Re the sliced MPS, seems the assumption is setting/getting a unit range

function deepcopy_ortho_center(M::AbstractMPS)
  M = copy(M)
  c = ortho_lims(M)
  # TODO: define `getindex(::AbstractMPS, I)` to return `AbstractMPS`
  M[c] = deepcopy(typeof(M)(M[c]))
  return M
end

What are your thoughts given the TODO Matt?

@mtfishman
Copy link
Member

mtfishman commented Oct 15, 2024

I think slicing an MPS/MPO with a unit range should return an MPS/MPO, if we make that change we can just update other code like the one you quote accordingly.

However, the definition in ITensorInfiniteMPS.jl: https://github.com/ITensor/ITensorInfiniteMPS.jl/blob/fc1078a60442081010447caf65c9a22c7278e4ce/src/ITensors.jl#L186 is too naive, as the comment states, since it doesn't preserve information about the orthogonality center.

@mtfishman mtfishman changed the title Port VUMPS functions [ITensors] Port VUMPS functions Oct 20, 2024
@ryanlevy
Copy link
Contributor Author

ryanlevy commented Oct 21, 2024

Just pushed a hopefully better performance version of sqrt. The original logic is a bit strange to me, it will fail for any matrix that isn't PSD without a large atol, how do you feel about allowing for conversion to complex like Base.sqrt does?

julia> eigvals(array(A,i,i'))
2-element Vector{Float64}:
 0.031100580805835387
 1.5475689615743486

julia> sqrt(Hermitian(array(A,i,i')))  array(sqrt(A),i,i')
true

julia> eigvals(array(B,i,i'))
2-element Vector{Float64}:
 -0.9688994191941646
  0.5475689615743488

jjulia> sqrt(array(B,i,i'))
2×2 Matrix{ComplexF64}:
 0.638411+0.135107im  0.254641-0.338726im
 0.254641-0.338726im  0.101568+0.84922im

julia> sqrt(B)
ERROR: DomainError with -0.9688994191941646:
sqrt was called with a negative real argument but will only return a complex result if called with a complex argument. Try sqrt(Complex(x)).
Stacktrace:
  [1] throw_complex_domainerror(f::Symbol, x::Float64)
    @ Base.Math ./math.jl:33
  [2] sqrt
    @ ./math.jl:608 [inlined]
...

@ryanlevy
Copy link
Contributor Author

ryanlevy commented Nov 1, 2024

@mtfishman now that ITensorMPS is out of ITensors I think this PR is complete pending your review and I'll move over there to fix the slice issue

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants