diff --git a/docs/src/man/states.md b/docs/src/man/states.md index 80b018d10..02cf66dbc 100644 --- a/docs/src/man/states.md +++ b/docs/src/man/states.md @@ -8,31 +8,53 @@ using LinearAlgebra: dot ## FiniteMPS -A [`FiniteMPS`](@ref) is - at its core - a chain of mps tensors. +A [`FiniteMPS`](@ref) is - at its core - a chain of MPS tensors. ```@raw html finite MPS ``` -### Usage +### Construction -A `FiniteMPS` can be created by passing in a vector of tensormaps: +If you already have the state of interest, it is straightforward to convert it to an MPS as follows: ```@example states -L = 10 -data = [rand(ComplexF64, ℂ^1 ⊗ ℂ^2 ← ℂ^1) for _ in 1:L]; -state = FiniteMPS(data) +ψ_dense = rand((ℂ^2)^4) +ψ_mps = FiniteMPS(ψ_dense) ``` -Or alternatively by specifying its structure +However, typically this is not the case, as storing the full state becomes expensive rather quickly. +Then, `FiniteMPS` are best constructed by first specifying a [`FiniteMPSManifold`](@ref) that encodes the physical and (maximal) virtual spaces: ```@example states -max_bond_dimension = ℂ^4 -physical_space = ℂ^2 -state = FiniteMPS(rand, ComplexF64, L, physical_space, max_bond_dimension) +pspaces = fill(ℂ^2, 10) # physical spaces (Vector) +max_virtualspace = ℂ^4 # single max virtual space +manifold = FiniteMPSManifold(pspaces, max_virtualspace) +ψ = rand(manifold) # random normalized FiniteMPS ``` -You can take dot products, renormalize!, expectation values,.... +Finally, it is also possible to build them from explicit MPS tensors, by passing them directly. +Here we construct the + +```@example states +As = [rand(ComplexF64, ℂ^1 ⊗ ℂ^2 ← ℂ^1) for _ in 1:L]; +ψ_from_As = FiniteMPS(As) +``` + +!!! warning "Full rank spaces" + + As the `FiniteMPS` object handles tensors in well-chosen gauges, the virtualspaces, as well as the associated tensors might reduce in size. + This can be achieved in a lossless manner whenever the spaces are not full rank, in the following sense: + ```julia + left_virtualspace(A) ⊗ physicalspace(A) ≿ right_virtualspace(A) && + left_virtualspace(A)' ≾ physicalspace(A) ⊗ right_virtualspace(A)' + ``` + +!!! warning "Edge spaces" + + It is possible for a `FiniteMPS` object to have non-trivial left- and/or right edge spaces. + This can be convenient whenever the state is embedded in a larger system (e.g. as part of a [`WindowMPS`](@ref)), or to allow for non-trivially charged symmetric states. + Therefore, be mindful that when constructing a `FiniteMPS` from tensors directly, you need to handle the edges separately. ### Gauging and canonical forms @@ -102,25 +124,26 @@ The idea behind this construction is that one never has to worry about how the s ## InfiniteMPS -An [`InfiniteMPS`](@ref) can be thought of as being very similar to a finite mps, where the set of tensors is repeated periodically. +An [`InfiniteMPS`](@ref) represents a periodically repeating unit cell of MPS tensors. -It can also be created by passing in a vector of `TensorMap`s: +### Construction + +Similar to `FiniteMPS`, the easiest way of constructing an `InfiniteMPS` is by specifying an [`InfiniteMPSManifold`](@ref) describing one unit cell: ```@example states -data = [rand(ComplexF64, ℂ^4 ⊗ ℂ^2 ← ℂ^4) for _ in 1:2] -state = InfiniteMPS(data) +pspaces = [ℂ^2, ℂ^2] # 2-site unit cell +vspaces = [ℂ^4, ℂ^5] # virtual space to the left of each site +manifold = InfiniteMPSManifold(pspaces, vspaces) +ψinf = rand(manifold) ``` -or by initializing it from given spaces +Alternatively, we may also start from explicit site tensors: ```@example states -phys_spaces = fill(ℂ^2, 2) -virt_spaces = [ℂ^4, ℂ^5] # by convention to the right of a site -state = InfiniteMPS(phys_spaces, virt_spaces) +As = [rand(ComplexF64, imanifold[i]) for i in 1:length(imanifold)] +ψinf2 = InfiniteMPS(As) ``` -Note that the code above creates an `InfiniteMPS` with a two-site unit cell, where the given virtual spaces are located to the right of their respective sites. - ### Gauging and canonical forms Much like for `FiniteMPS`, we can again query the gauged tensors `AL`, `AR`, `C` and `AC`. @@ -142,8 +165,9 @@ It represents a window of mutable tensors (a finite MPS), embedded in an infinit It can therefore be created accordingly, ensuring that the edges match: ```@example states -infinite_state = InfiniteMPS(ℂ^2, ℂ^4) -finite_state = FiniteMPS(5, ℂ^2, ℂ^4; left=ℂ^4, right=ℂ^4) +infinite_state = rand(InfiniteMPSManifold(ℂ^2, ℂ^4)) +finite_manifold = FiniteMPSManifold(fill(ℂ^2, 5), ℂ^4; left_virtualspace=ℂ^4, right_virtualspace=ℂ^4) +finite_state = rand(ComplexF64, finite_manifold) window = WindowMPS(infinite_state, finite_state, infinite_state) ``` diff --git a/src/states/abstractmps.jl b/src/states/abstractmps.jl index 6dd8a1341..d2f2511b4 100644 --- a/src/states/abstractmps.jl +++ b/src/states/abstractmps.jl @@ -85,74 +85,6 @@ function MPSTensor(A::AbstractArray{T}) where {T <: Number} t[] .= A return t end - -""" - isfullrank(A::GenericMPSTensor; side=:both) - -Determine whether the given tensor is full rank, i.e. whether both the map from the left -virtual space and the physical space to the right virtual space, and the map from the right -virtual space and the physical space to the left virtual space are injective. -""" -isfullrank(A::GenericMPSTensor; kwargs...) = isfullrank(space(A); kwargs...) -function isfullrank(V::TensorKit.TensorMapSpace; side = :both) - Vₗ = V[1] - Vᵣ = V[numind(V)] - P = ⊗(getindex.(Ref(V), 2:(numind(V) - 1))...) - return if side === :both - Vₗ ⊗ P ≿ Vᵣ' && Vₗ' ≾ P ⊗ Vᵣ - elseif side === :right - Vₗ ⊗ P ≿ Vᵣ' - elseif side === :left - Vₗ' ≾ P ⊗ Vᵣ - else - throw(ArgumentError("Invalid side: $side")) - end -end - -""" - makefullrank!(A::PeriodicVector{<:GenericMPSTensor}; alg=Defalts.alg_qr()) - -Make the set of MPS tensors full rank by performing a series of orthogonalizations. -""" -function makefullrank!(A::PeriodicVector{<:GenericMPSTensor}; alg = Defaults.alg_qr()) - while true - i = findfirst(!isfullrank, A) - isnothing(i) && break - if !isfullrank(A[i]; side = :left) - L, Q = _right_orth!(_transpose_tail(A[i]); alg) - A[i] = _transpose_front(Q) - A[i - 1] = A[i - 1] * L - else - A[i], R = _left_orth!(A[i]; alg) - A[i + 1] = _transpose_front(R * _transpose_tail(A[i + 1])) - end - end - return A -end - -function makefullrank!(virtualspaces::PeriodicVector{S}, physicalspaces::PeriodicVector{S}) where {S <: ElementarySpace} - haschanged = true - while haschanged - haschanged = false - for i in 1:length(virtualspaces) - Vmax = fuse(virtualspaces[i - 1], physicalspaces[i - 1]) - if !(virtualspaces[i] ≾ Vmax) - virtualspaces[i] = infimum(virtualspaces[i], Vmax) - haschanged = true - end - end - for i in reverse(1:length(virtualspaces)) - Vmax = fuse(dual(physicalspaces[i - 1]), virtualspaces[i]) - if !(virtualspaces[i - 1] ≾ Vmax) - virtualspaces[i - 1] = infimum(virtualspaces[i - 1], Vmax) - haschanged = true - end - end - end - - return virtualspaces -end - # Tensor accessors # ---------------- @doc """ @@ -164,9 +96,9 @@ If this hasn't been computed before, this can be computed as: - `kind=:ALAC` : AL[i] * AC[i+1] """ AC2 -#=========================================================================================== -MPS types -===========================================================================================# +# =========================================================================================== +# MPS types +# =========================================================================================== abstract type AbstractMPS end abstract type AbstractFiniteMPS <: AbstractMPS end @@ -201,7 +133,7 @@ TensorKit.sectortype(ψtype::Type{<:AbstractMPS}) = sectortype(site_type(ψtype) """ left_virtualspace(ψ::AbstractMPS, [pos=1:length(ψ)]) - + Return the virtual space of the bond to the left of sites `pos`. !!! warning @@ -245,3 +177,307 @@ physicalspace(ψ::AbstractMPS) = map(Base.Fix1(physicalspace, ψ), eachsite(ψ)) Return an iterator over the sites of the MPS `state`. """ eachsite(ψ::AbstractMPS) = eachindex(ψ) + +# =========================================================================================== +# MPS manifolds +# =========================================================================================== +""" + abstract type AbstractMPSManifold{S <: ElementarySpace} + +Abstract supertype for the characterization of an MPS manifold, i.e. the sizes of all tensors that are involved. +These types are used mostly as convenient dispatch +""" +abstract type AbstractMPSManifold{S <: ElementarySpace} end + +TensorKit.spacetype(::Type{<:AbstractMPSManifold{S}}) where {S} = S + +physicalspace(manifold::AbstractMPSManifold, i::Integer) = physicalspace(manifold)[i] +left_virtualspace(manifold::AbstractMPSManifold, i::Integer) = left_virtualspace(manifold)[i] +right_virtualspace(manifold::AbstractMPSManifold, i::Integer) = right_virtualspace(manifold)[i] + +Base.getindex(manifold::AbstractMPSManifold, site::Integer) = + left_virtualspace(manifold, site) ⊗ physicalspace(manifold, site) ← right_virtualspace(manifold, site) + +""" + FiniteMPSManifold(pspaces, vspaces) <: AbstractMPSManifold{S} + +Full characterization of all [`FiniteMPS`](@ref) spaces. Both `pspaces` and `vspaces` are `Vector`s that hold the +local physical and virtual spaces, such that we have `length(pspaces) + 1 == length(vspaces)`. +These objects can be used to construct MPS for a given space, for example through [`Base.randn`](@ref) or similar functions. + +## Constructors + + FiniteMPSManifold(physicalspaces::AbstractVector{<:TensorSpace}, max_virtualspaces; [left_virtualspace], [right_virtualspace]) + +To construct a `FiniteMPSManifold`, you should provide the `physicalspaces`, along with the maximal desired virtual space. +This latter can be provided as a single `<:ElementarySpace`, or, if a site-dependent maximum is desired through a vector thereof. +In that case, `length(physicalspaces) == length(max_virtualspaces) + 1`. + +It is also possible to add non-trivial virtual spaces on both the left and the right edge of the MPS. +This might be useful to construct "charged" MPS, or to work with [`WindowMPS`](@ref) that may be embedded in a larger system. + +!!! note + As the supplied virtual spaces are interpreted as a maximum, this function will automatically construct "full-rank" spaces. + These are chosen such that every MPS tensor that would be generated could have full rank. See also [`makefullrank!](@ref). +""" +struct FiniteMPSManifold{S <: ElementarySpace, S′ <: TensorSpace{S}} <: AbstractMPSManifold{S} + pspaces::Vector{S′} + vspaces::Vector{S} + + # disable default constructor + function FiniteMPSManifold{S, S′}( + pspaces::Vector{S′}, vspaces::Vector{S} + ) where {S <: ElementarySpace, S′ <: TensorSpace{S}} + return new{S, S′}(pspaces, vspaces) + end +end + +function FiniteMPSManifold( + physicalspaces::AbstractVector{S′}, max_virtualspaces::AbstractVector{S}; + left_virtualspace::S = oneunit(S), right_virtualspace::S = oneunit(S) + ) where {S <: ElementarySpace, S′ <: TensorSpace{S}} + L₁ = length(physicalspaces) + L₂ = length(max_virtualspaces) + L₁ == L₂ + 1 || + throw(DimensionMismatch(lazy"`|physicalspaces|` ($L₁) should be 1 more than `|max_virtualspaces|` ($L₂)")) + + # copy to avoid side-effects and get correct array type + pspaces = collect(physicalspaces) + vspaces = vcat(left_virtualspace, max_virtualspaces, right_virtualspace) + manifold = FiniteMPSManifold{S, S′}(pspaces, vspaces) + + # ensure all spaces are full rank -- use vspaces as maximum + return makefullrank!(manifold) +end +function FiniteMPSManifold( + physicalspaces::AbstractVector{S′}, max_virtualspace::S; kwargs... + ) where {S <: ElementarySpace, S′ <: TensorSpace{S}} + return FiniteMPSManifold(physicalspaces, fill(max_virtualspace, length(physicalspaces) - 1); kwargs...) +end +function FiniteMPSManifold(mps_tensors::AbstractVector{A}) where {A <: GenericMPSTensor} + numin(A) == 1 || throw(ArgumentError("Not a valid MPS tensor space")) + pspaces = map(physicalspace, mps_tensors) + vspaces = Vector{spacetype(A)}(undef, length(mps_tensors) - 1) + for (i, mps_tensor) in enumerate(mps_tensors) + i == length(mps_tensors) && continue + vspaces[i] = right_virtualspace(mps_tensor) + end + return FiniteMPSManifold( + pspaces, vspaces; + left_virtualspace = left_virtualspace(first(mps_tensors)), + right_virtualspace = right_virtualspace(last(mps_tensors)) + ) +end + +""" + InfiniteMPSManifold(pspaces, vspaces) <: AbstractMPSManifold{S} + +Full characterization of all [`InfiniteMPS`](@ref) spaces. Both `pspaces` and `vspaces` are `PeriodicVector`s that hold the +local physical and virtual spaces, such that we must have `length(pspaces) == length(vspaces)`. +These objects can be used to construct MPS for a given space, for example through [`Base.randn`](@ref) or similar functions. + +## Constructors + + FiniteMPSManifold(physicalspaces::AbstractVector{<:TensorSpace}, max_virtualspaces; [left_virtualspace], [right_virtualspace]) + +To construct a `FiniteMPSManifold`, you should provide the `physicalspaces`, along with the maximal desired virtual space. +This latter can be provided as a single `<:ElementarySpace`, or, if a site-dependent maximum is desired through a vector thereof. +In that case, `length(physicalspaces) == length(max_virtualspaces)`. + +!!! note + As the supplied virtual spaces are interpreted as a maximum, this function will automatically construct "full-rank" spaces. + These are chosen such that every MPS tensor that would be generated could have full rank. See also [`makefullrank!](@ref). +""" +struct InfiniteMPSManifold{S <: ElementarySpace, S′ <: Union{S, CompositeSpace{S}}} <: AbstractMPSManifold{S} + pspaces::PeriodicVector{S′} + vspaces::PeriodicVector{S} +end + +function InfiniteMPSManifold( + physicalspaces::AbstractVector{S′}, virtualspaces::AbstractVector{S} + ) where {S <: ElementarySpace, S′ <: Union{S, CompositeSpace{S}}} + L₁ = length(physicalspaces) + L₂ = length(virtualspaces) + L₁ == L₂ || + throw(DimensionMismatch(lazy"`|physicalspaces|` ($L₁) should be equal to `|virtualspaces|` ($L₂)")) + + # copy to avoid side-effects and get correct array type + pspaces = collect(physicalspaces) + vspaces = collect(virtualspaces) + manifold = InfiniteMPSManifold{S, S′}(pspaces, vspaces) + + # ensure all spaces are full rank -- use vspaces as maximum + return makefullrank!(manifold) +end +function InfiniteMPSManifold(mps_tensors::AbstractVector{A}) where {A <: GenericMPSTensor} + pspaces = PeriodicVector(map(physicalspace, mps_tensors)) + vspaces = PeriodicVector(map(left_virtualspace, mps_tensors)) + for i in eachindex(vspaces) + vspaces[i] == right_virtualspace(mps_tensors[i - 1]) || + throw(SpaceMismatch("incompatible spaces between site $(i - 1) and $i")) + end + return InfiniteMPSManifold(pspaces, vspaces) +end + +Base.length(manifold::AbstractMPSManifold) = length(physicalspace(manifold)) + +physicalspace(manifold::Union{FiniteMPSManifold, InfiniteMPSManifold}) = manifold.pspaces +left_virtualspace(manifold::FiniteMPSManifold) = manifold.vspaces[1:(end - 1)] +left_virtualspace(manifold::InfiniteMPSManifold) = manifold.vspaces +right_virtualspace(manifold::FiniteMPSManifold) = manifold.vspaces[2:end] +right_virtualspace(manifold::InfiniteMPSManifold) = PeriodicVector(circshift(manifold.vspaces, 1)) + +# Utility +function Base.repeat(manifold::FiniteMPSManifold, i::Integer) + last(manifold.vspaces) == first(manifold.vspaces) || throw(SpaceMismatch()) + pspaces = repeat(manifold.pspaces, i) + vspaces = push!(repeat(left_virtualspace(manifold), i), last(manifold.vspaces)) + return FiniteMPSManifold(pspaces, vspaces) +end +function Base.repeat(manifold::InfiniteMPSManifold, i::Integer) + pspaces = repeat(manifold.pspaces, i) + vspaces = repeat(manifold.vspaces, i) + return InfiniteMPSManifold(pspaces, vspaces) +end + +function Base.vcat(manifold1::FiniteMPSManifold, manifold2::FiniteMPSManifold) + last(right_virtualspace(manifold1)) == first(left_virtualspace(manifold2)) || throw(SpaceMismatch()) + pspaces = vcat(manifold1.pspaces, manifold2.pspaces) + vspaces = push!(vcat(left_virtualspace(manifold1), left_virtualspace(manifold2)), last(manifold2.vspaces)) + return FiniteMPSManifold(pspaces, vspaces) +end +Base.vcat(manifold::FiniteMPSManifold, manifolds::FiniteMPSManifold...) = foldl(vcat, (manifold, manifolds...)) +function Base.vcat(manifold::InfiniteMPSManifold, manifolds::InfiniteMPSManifold...) + pspaces = vcat(manifold.pspaces, Base.Fix2(getproperty, :pspaces).(manifolds)...) + vspaces = vcat(manifold.vspaces, Base.Fix2(getproperty, :vspaces).(manifolds)...) + return InfiniteMPSManifold(pspaces, vspaces) +end + +# MPS constructors +# ---------------- +for randf in (:rand, :randn, :randexp, :randisometry) + _docstr = """ + $randf([rng=default_rng()], [T=Float64], manifold::AbstractMPSManifold) -> mps + + Generate an `mps` with tensors generated by `$randf`. + + See also [`Random.$(randf)!`](@ref). + """ + _docstr! = """ + $(randf)!([rng=default_rng()], mps::AbstractMPS) -> mps + + Fill the tensors of `mps` with entries generated by `$(randf)!`. + + See also [`Random.$(randf)`](@ref). + """ + + if randf != :randisometry + randfun = GlobalRef(Random, randf) + randfun! = GlobalRef(Random, Symbol(randf, :!)) + else + randfun = randf + randfun! = Symbol(randf, :!) + end + + @eval begin + @doc $_docstr $randfun(::Type, ::AbstractMPSManifold) + @doc $_docstr! $randfun!(::AbstractMPSManifold) + + # filling in default eltype + $randfun(manifold::AbstractMPSManifold) = $randfun(Defaults.eltype, manifold) + function $randfun(rng::Random.AbstractRNG, manifold::AbstractMPSManifold) + return $randfun(rng, Defaults.eltype, manifold) + end + + # filling in default rng + function $randfun(::Type{T}, manifold::AbstractMPSManifold) where {T} + return $randfun(Random.default_rng(), T, manifold) + end + $randfun!(mps::AbstractMPS) = $randfun!(Random.default_rng(), mps) + end +end + +""" + isfullrank(A::GenericMPSTensor; side=:both) + +Determine whether the given tensor is full rank, i.e. whether both the map from the left +virtual space and the physical space to the right virtual space, and the map from the right +virtual space and the physical space to the left virtual space are injective. +""" +isfullrank(A::GenericMPSTensor; kwargs...) = isfullrank(space(A); kwargs...) +function isfullrank(V::TensorKit.TensorMapSpace; side = :both) + Vₗ = V[1] + Vᵣ = V[numind(V)] + P = ⊗(getindex.(Ref(V), 2:(numind(V) - 1))...) + return if side === :both + Vₗ ⊗ P ≿ Vᵣ' && Vₗ' ≾ P ⊗ Vᵣ + elseif side === :right + Vₗ ⊗ P ≿ Vᵣ' + elseif side === :left + Vₗ' ≾ P ⊗ Vᵣ + else + throw(ArgumentError("Invalid side: $side")) + end +end +isfullrank(manifold::AbstractMPSManifold; kwargs...) = all(i -> isfullrank(manifold[i]; kwargs...), 1:length(manifold)) + +""" + makefullrank!(A::PeriodicVector{<:GenericMPSTensor}; alg=Defalts.alg_qr()) + +Make the set of MPS tensors full rank by performing a series of orthogonalizations. +""" +function makefullrank!(manifold::FiniteMPSManifold) + # left-to-right sweep + for site in 1:length(manifold) + if !isfullrank(manifold[site]; side = :right) + maxspace = fuse(left_virtualspace(manifold, site), fuse(physicalspace(manifold, site))) + manifold.vspaces[site + 1] = infimum(right_virtualspace(manifold, site), maxspace) + end + end + # right-to-left sweep + for site in reverse(1:length(manifold)) + if !isfullrank(manifold[site]; side = :left) + maxspace = fuse(right_virtualspace(manifold, site), dual(fuse(physicalspace(manifold, site)))) + manifold.vspaces[site] = infimum(left_virtualspace(manifold, site), maxspace) + end + end + return manifold +end +function makefullrank!(manifold::InfiniteMPSManifold) + haschanged = true + while haschanged + haschanged = false + # left-to-right sweep + for site in 1:length(manifold) + if !isfullrank(manifold[site]; side = :right) + maxspace = fuse(left_virtualspace(manifold, site), fuse(physicalspace(manifold, site))) + manifold.vspaces[site + 1] = infimum(right_virtualspace(manifold, site), maxspace) + haschanged = true + end + end + # right-to-left sweep + for site in reverse(1:length(manifold)) + if !isfullrank(manifold[site]; side = :left) + maxspace = fuse(right_virtualspace(manifold, site), dual(fuse(physicalspace(manifold, site)))) + manifold.vspaces[site] = infimum(left_virtualspace(manifold, site), maxspace) + haschanged = true + end + end + end + return manifold +end +function makefullrank!(A::PeriodicVector{<:GenericMPSTensor}; alg = Defaults.alg_qr()) + while true + i = findfirst(!isfullrank, A) + isnothing(i) && break + if !isfullrank(A[i]; side = :left) + L, Q = _right_orth!(_transpose_tail(A[i]); alg) + A[i] = _transpose_front(Q) + A[i - 1] = A[i - 1] * L + else + A[i], R = _left_orth!(A[i]; alg) + A[i + 1] = _transpose_front(R * _transpose_tail(A[i + 1])) + end + end + return A +end diff --git a/src/states/finitemps.jl b/src/states/finitemps.jl index 44ac3fa64..015d031e5 100644 --- a/src/states/finitemps.jl +++ b/src/states/finitemps.jl @@ -25,124 +25,126 @@ By convention, we have that: --- ## Constructors - FiniteMPS([f, eltype], physicalspaces::Vector{<:Union{S,CompositeSpace{S}}}, - maxvirtualspaces::Union{S,Vector{S}}; - normalize=true, left=oneunit(S), right=oneunit(S)) where {S<:ElementarySpace} - FiniteMPS([f, eltype], N::Int, physicalspace::Union{S,CompositeSpace{S}}, - maxvirtualspaces::Union{S,Vector{S}}; - normalize=true, left=oneunit(S), right=oneunit(S)) where {S<:ElementarySpace} - FiniteMPS(As::Vector{<:GenericMPSTensor}; normalize=false, overwrite=false) - -Construct an MPS via a specification of physical and virtual spaces, or from a list of -tensors `As`. All cases reduce to the latter. In particular, a state with a non-trivial -total charge can be constructed by passing a non-trivially charged vector space as the -`left` or `right` virtual spaces. - -### Arguments -- `As::Vector{<:GenericMPSTensor}`: vector of site tensors - -- `f::Function=rand`: initializer function for tensor data -- `eltype::Type{<:Number}=ComplexF64`: scalar type of tensors - -- `physicalspaces::Vector{<:Union{S, CompositeSpace{S}}`: list of physical spaces -- `N::Int`: number of sites -- `physicalspace::Union{S,CompositeSpace{S}}`: local physical space - -- `virtualspaces::Vector{<:Union{S, CompositeSpace{S}}`: list of virtual spaces -- `maxvirtualspace::S`: maximum virtual space - -### Keywords -- `normalize=true`: normalize the constructed state -- `overwrite=false`: overwrite the given input tensors -- `left=oneunit(S)`: left-most virtual space -- `right=oneunit(S)`: right-most virtual space + +Recommended ways to construct a finite MPS are: + +- Using an MPS manifold of spaces + + ```julia + rand([rng], [T], manifold::FiniteMPSManifold) + randn([rng], [T], manifold::FiniteMPSManifold) + zeros([T], manifold::FiniteMPSManifold) + ones([T], manifold::FiniteMPSManifold) + ``` + + First build a [`FiniteMPSManifold`](@ref) that fixes the physical and (maximal) virtual spaces, then allocate an MPS on that manifold. + See [`FiniteMPSManifold`](@ref) for how to specify site-dependent virtual-space bounds and nontrivial edge charges via `left_virtualspace` and `right_virtualspace`. + +- From site tensors + + ```julia + FiniteMPS(As::Vector{<:GenericMPSTensor}; normalize = false) + ``` + + Construct an MPS from a vector of already-sized MPS site tensors `As`. + When `normalize = true`, the state is normalized during construction. + +- From a full state tensor + + ```julia + FiniteMPS(ψ::AbstractTensor) + ``` + + Factorizes a full many-body state `ψ` into a finite MPS (using a left-canonical sweep). + +In particular, charged MPS can be created by giving nontrivial left and/or right virtual spaces when constructing the manifold: + +```julia +ps = fill(ℂ^2, N) # physical spaces +m = FiniteMPSManifold(ps, ℂ^D; left_virtualspace = Qₗ, right_virtualspace = Qᵣ) +ψ = rand(ComplexF64, m) # normalized random MPS on the manifold +``` + +!!! warning "Deprecated constructors" + Older constructors of the form `FiniteMPS([f, eltype], physicalspaces, max_virtualspaces; ...)` + or `FiniteMPS([f, eltype], N, physicalspace, max_virtualspaces; ...)` are deprecated. Use + `rand`/`randn`/`zeros`/`ones` with a [`FiniteMPSManifold`](@ref) instead. """ struct FiniteMPS{A <: GenericMPSTensor, B <: MPSBondTensor} <: AbstractFiniteMPS ALs::Vector{Union{Missing, A}} ARs::Vector{Union{Missing, A}} ACs::Vector{Union{Missing, A}} Cs::Vector{Union{Missing, B}} + + # constructor from data function FiniteMPS{A, B}( ALs::Vector{Union{Missing, A}}, ARs::Vector{Union{Missing, A}}, ACs::Vector{Union{Missing, A}}, Cs::Vector{Union{Missing, B}} ) where {A <: GenericMPSTensor, B <: MPSBondTensor} return new{A, B}(ALs, ARs, ACs, Cs) end - function FiniteMPS( - ALs::Vector{MA}, ARs::Vector{MA}, - ACs::Vector{MA}, - Cs::Vector{MB} - ) where {MA <: Union{GenericMPSTensor, Missing}, MB <: Union{MPSBondTensor, Missing}} - A = _not_missing_type(MA) - B = _not_missing_type(MB) - length(ACs) == length(Cs) - 1 == length(ALs) == length(ARs) || - throw(DimensionMismatch("length mismatch of tensors")) - sum(ismissing.(ACs)) + sum(ismissing.(Cs)) < length(ACs) + length(Cs) || - throw(ArgumentError("at least one AC/C should not be missing")) - - S = spacetype(A) - left_virt_spaces = Vector{Union{Missing, S}}(missing, length(Cs)) - right_virt_spaces = Vector{Union{Missing, S}}(missing, length(Cs)) - - for (i, tup) in enumerate(zip(ALs, ARs, ACs)) - non_missing = filter(!ismissing, tup) - isempty(non_missing) && throw(ArgumentError("missing site tensor")) - (al, ar, ac) = tup - - if !ismissing(al) - !ismissing(left_virt_spaces[i]) && - ( - left_virt_spaces[i] == _firstspace(al) || - throw(SpaceMismatch("Virtual space of AL on site $(i) doesn't match")) - ) - - left_virt_spaces[i + 1] = _lastspace(al)' - left_virt_spaces[i] = _firstspace(al) - end - - if !ismissing(ar) - !ismissing(right_virt_spaces[i]) && - ( - right_virt_spaces[i] == _firstspace(ar) || - throw(SpaceMismatch("Virtual space of AR on site $(i) doesn't match")) - ) - - right_virt_spaces[i + 1] = _lastspace(ar)' - right_virt_spaces[i] = _firstspace(ar) - end - - if !ismissing(ac) - !ismissing(left_virt_spaces[i]) && - ( - left_virt_spaces[i] == _firstspace(ac) || - throw(SpaceMismatch("Left virtual space of AC on site $(i) doesn't match")) - ) - !ismissing(right_virt_spaces[i + 1]) && - ( - right_virt_spaces[i + 1] == _lastspace(ac)' || - throw(SpaceMismatch("Right virtual space of AC on site $(i) doesn't match")) - ) - - right_virt_spaces[i + 1] = _lastspace(ac)' - left_virt_spaces[i] = _firstspace(ac) - end - end - for (i, c) in enumerate(Cs) - ismissing(c) && continue - !ismissing(left_virt_spaces[i]) && ( - left_virt_spaces[i] == _firstspace(c) || - throw(SpaceMismatch("Left virtual space of C on site $(i - 1) doesn't match")) - ) - !ismissing(right_virt_spaces[i]) && ( - right_virt_spaces[i] == _lastspace(c)' || - throw(SpaceMismatch("Right virtual space of C on site $(i - 1) doesn't match")) - ) - end + # undef constructor + function FiniteMPS{A, B}(::UndefInitializer, L::Integer) where {A, B} + ALs = Vector{Union{Missing, A}}(missing, L) + ARs = Vector{Union{Missing, A}}(missing, L) + ACs = Vector{Union{Missing, A}}(missing, L) + Cs = Vector{Union{Missing, B}}(missing, L + 1) return new{A, B}(ALs, ARs, ACs, Cs) end end +function FiniteMPS( + ALs::Vector{MA}, ARs::Vector{MA}, ACs::Vector{MA}, Cs::Vector{MB} + ) where {MA <: Union{GenericMPSTensor, Missing}, MB <: Union{MPSBondTensor, Missing}} + (L = length(ACs)) == length(Cs) - 1 == length(ALs) == length(ARs) || + throw(DimensionMismatch("length mismatch of tensors")) + sum(ismissing.(ACs)) + sum(ismissing.(Cs)) < length(ACs) + length(Cs) || + throw(ArgumentError("at least one AC/C should not be missing")) + + # determine the MPS manifold from the input and instantiate the MPS + mps_tensors = map(ALs, ARs, ACs) do AL, AR, AC + mps_tensor = coalesce(AL, AR, AC) + ismissing(mps_tensor) && throw(ArgumentError("missing site tensor at site $i")) + return mps_tensor + end + manifold = FiniteMPSManifold(mps_tensors) + A = _not_missing_type(MA) + B = _not_missing_type(MB) + mps = FiniteMPS{A, B}(undef, length(manifold)) + + # populate the non-missing tensors into the MPS + for i in 1:L + V_mps = manifold[i] + if !ismissing(ALs[i]) + space(ALs[i]) == V_mps || throw(SpaceMismatch("incompatible space for AL[$i]")) + getfield(mps, :ALs)[i] = ALs[i] + end + if !ismissing(ARs[i]) + space(ARs[i]) == V_mps || throw(SpaceMismatch("incompatible space for AR[$i]")) + getfield(mps, :ARs)[i] = ARs[i] + end + if !ismissing(ACs[i]) + space(ACs[i]) == V_mps || throw(SpaceMismatch("incompatible space for AC[$i]")) + getfield(mps, :ACs)[i] = ACs[i] + end + if !ismissing(Cs[i]) + V = left_virtualspace(manifold, i) + space(Cs[i]) == (V ← V) || + throw(SpaceMismatch("incompatible space for C[$i]")) + getfield(mps, :Cs)[i] = Cs[i] + end + end + if !ismissing(Cs[end]) + V = right_virtualspace(manifold, L) + space(Cs[end]) == (V ← V) || + throw(SpaceMismatch("incompatible space for C[$(L + 1)]")) + getfield(mps, :Cs)[L + 1] = Cs[L + 1] + end + + return mps +end + _not_missing_type(::Type{Missing}) = throw(ArgumentError("Only missing type present")) function _not_missing_type(::Type{T}) where {T} if T isa Union @@ -202,91 +204,33 @@ function _gaugecenter(ψ::FiniteMPS)::HalfInt isnothing(center) && throw(ArgumentError("No center found, invalid state")) return center end + #=========================================================================================== Constructors ===========================================================================================# -function FiniteMPS(As::Vector{<:GenericMPSTensor}; normalize = false, overwrite = false) - # TODO: copying the input vector is probably not necessary, as we are constructing new - # vectors anyways, maybe deprecate `overwrite`. - As = overwrite ? As : copy(As) - N = length(As) - As[1] = MatrixAlgebraKit.copy_input(qr_compact, As[1]) - local C - for i in eachindex(As) - As[i], C = qr_compact!(As[i]; positive = true) - normalize && normalize!(C) - i == N || (As[i + 1] = _transpose_front(C * _transpose_tail(As[i + 1]))) - end +function FiniteMPS(As::Vector{<:GenericMPSTensor}; normalize::Bool = false) + mps = FiniteMPS - A = eltype(As) - B = typeof(C) - - Cs = Vector{Union{Missing, B}}(missing, N + 1) - ALs = Vector{Union{Missing, A}}(missing, N) - ARs = Vector{Union{Missing, A}}(missing, N) - ACs = Vector{Union{Missing, A}}(missing, N) - - ALs .= As - Cs[end] = C - - return FiniteMPS(ALs, ARs, ACs, Cs) -end + # start with first in case eltype changes + AL1, C = qr_compact(first(As)) + normalize && normalize!(C) -function FiniteMPS( - f, elt, Pspaces::Vector{<:Union{S, CompositeSpace{S}}}, maxVspaces::Vector{S}; - normalize = true, left::S = oneunit(S), right::S = oneunit(S) - ) where {S <: ElementarySpace} - N = length(Pspaces) - length(maxVspaces) == N - 1 || - throw(DimensionMismatch("length of physical spaces ($N) and virtual spaces $(length(maxVspaces)) should differ by 1")) - - # limit the maximum virtual dimension such that result is full rank - fusedPspaces = fuse.(Pspaces) # for working with multiple physical spaces - Vspaces = similar(maxVspaces, N + 1) - - Vspaces[1] = left - for k in 2:N - Vspaces[k] = infimum(fuse(Vspaces[k - 1], fusedPspaces[k - 1]), maxVspaces[k - 1]) - dim(Vspaces[k]) > 0 || @warn "no fusion channels available at site $k" - end + # instantiate the destination + A = typeof(AL1) + B = typeof(C) + mps = FiniteMPS{A, B}(undef, length(As)) - Vspaces[end] = right - for k in reverse(2:N) - Vspaces[k] = infimum(Vspaces[k], fuse(Vspaces[k + 1], dual(fusedPspaces[k]))) - dim(Vspaces[k]) > 0 || @warn "no fusion channels available at site $k" + # left-gauge and fill the ALs and last C + ALs = getfield(mps, :ALs) + for i in eachindex(ALs) + i == 1 && (ALs[i] = AL1; continue) + ALs[i], C = qr_compact!(_mul_front(C, As[i])) + normalize && normalize!(C) end + getfield(mps, :Cs)[end] = C - # construct MPS - tensors = MPSTensor.(f, elt, Pspaces, Vspaces[1:(end - 1)], Vspaces[2:end]) - return FiniteMPS(tensors; normalize, overwrite = true) -end -function FiniteMPS( - f, elt, Pspaces::Vector{<:Union{S, CompositeSpace{S}}}, maxVspace::S; - kwargs... - ) where {S <: ElementarySpace} - maxVspaces = fill(maxVspace, length(Pspaces) - 1) - return FiniteMPS(f, elt, Pspaces, maxVspaces; kwargs...) -end -function FiniteMPS( - Pspaces::Vector{<:Union{S, CompositeSpace{S}}}, maxVspaces::Union{S, Vector{S}}; - kwargs... - ) where {S <: ElementarySpace} - return FiniteMPS(rand, Defaults.eltype, Pspaces, maxVspaces; kwargs...) -end - -# Also accept single physical space and length -function FiniteMPS(N::Int, V::VectorSpace, args...; kwargs...) - return FiniteMPS(fill(V, N), args...; kwargs...) -end -function FiniteMPS(f, elt, N::Int, V::VectorSpace, args...; kwargs...) - return FiniteMPS(f, elt, fill(V, N), args...; kwargs...) -end - -# Also accept ProductSpace of physical spaces -FiniteMPS(P::ProductSpace, args...; kwargs...) = FiniteMPS(collect(P), args...; kwargs...) -function FiniteMPS(f, elt, P::ProductSpace, args...; kwargs...) - return FiniteMPS(f, elt, collect(P), args...; kwargs...) + return mps end # construct from dense state @@ -296,9 +240,74 @@ function FiniteMPS(ψ::AbstractTensor) A = _transpose_front( U * transpose(ψ * U', ((), reverse(ntuple(identity, numind(ψ) + 1)))) ) - return FiniteMPS(decompose_localmps(A); normalize = false, overwrite = true) + return FiniteMPS(decompose_localmps(A); normalize = false) end +for f in (:zeros, :ones) + @eval begin + Base.$f(manifold::FiniteMPSManifold) = $f(Defaults.eltype, manifold) + function Base.$f(::Type{T}, manifold::FiniteMPSManifold) where {T} + As = map(i -> $f(T, manifold[i]), 1:length(manifold)) + return FiniteMPS(As) + end + end +end + +for randfun in (:rand, :randn) + randfun! = Symbol(randfun, :!) + @eval function Random.$randfun(rng::Random.AbstractRNG, ::Type{T}, manifold::FiniteMPSManifold) where {T} + As = map(i -> $randfun(rng, T, manifold[i]), 1:length(manifold)) + return normalize!(FiniteMPS(As)) + end + @eval function Random.$randfun!(rng::Random.AbstractRNG, mps::FiniteMPS) + for i in 1:length(mps) + mps.AC[i] = $randfun!(rng, mps.AC[i]) + end + return mps + end +end + +# Deprecate old constructor syntaxes +# ---------------------------------- +Base.@deprecate( + FiniteMPS( + f, elt, Pspaces::Vector{<:TensorSpace{S}}, maxVspaces::Vector{S}; left::S = oneunit(S), right::S = oneunit(S) + ) where {S <: ElementarySpace}, + f(elt, FiniteMPSManifold(Pspaces, maxVspaces; left_virtualspace = left, right_virtualspace = right)) +) +Base.@deprecate( + FiniteMPS( + f, elt, Pspaces::Vector{<:TensorSpace{S}}, maxVspace::S; left::S = oneunit(S), right::S = oneunit(S) + ) where {S <: ElementarySpace}, + f(elt, FiniteMPSManifold(Pspaces, maxVspace; left_virtualspace = left, right_virtualspace = right)) +) +Base.@deprecate( + FiniteMPS( + Pspaces::Vector{<:TensorSpace{S}}, maxVspaces::Union{S, Vector{S}}; left::S = oneunit(S), right::S = oneunit(S) + ) where {S <: ElementarySpace}, + rand(FiniteMPSManifold(Pspaces, maxVspaces; left_virtualspace = left, right_virtualspace = right)) +) +Base.@deprecate( + FiniteMPS( + f, elt, N::Int, V::TensorSpace{S}, args...; left::S = oneunit(S), right::S = oneunit(S) + ) where {S <: ElementarySpace}, + f(elt, FiniteMPSManifold(fill(V, N), args...; left_virtualspace = left, right_virtualspace = right)) +) +Base.@deprecate( + FiniteMPS( + N::Int, V::TensorSpace{S}, args...; left::S = oneunit(S), right::S = oneunit(S) + ) where {S <: ElementarySpace}, + rand(FiniteMPSManifold(fill(V, N), args...; left_virtualspace = left, right_virtualspace = right)) +) +Base.@deprecate( + FiniteMPS(P::ProductSpace, args...; kwargs...), + rand(FiniteMPSManifold(collect(P), args...; kwargs...)) +) +Base.@deprecate( + FiniteMPS(f, elt, P::ProductSpace, args...; kwargs...), + f(elt, FiniteMPSManifold(collect(P), args...; kwargs...)) +) + #=========================================================================================== Utility ===========================================================================================# diff --git a/src/states/infinitemps.jl b/src/states/infinitemps.jl index 62c278afb..968c7506b 100644 --- a/src/states/infinitemps.jl +++ b/src/states/infinitemps.jl @@ -18,36 +18,66 @@ By convention, we have that: --- ## Constructors - InfiniteMPS([f, eltype], physicalspaces::Vector{<:Union{S, CompositeSpace{S}}, - virtualspaces::Vector{<:Union{S, CompositeSpace{S}}; - kwargs...) where {S<:ElementarySpace} - InfiniteMPS(As::AbstractVector{<:GenericMPSTensor}; kwargs...) - InfiniteMPS(ALs::AbstractVector{<:GenericMPSTensor}, C₀::MPSBondTensor; - kwargs...) - -Construct an MPS via a specification of physical and virtual spaces, or from a list of -tensors `As`, or a list of left-gauged tensors `ALs`. - -### Arguments -- `As::AbstractVector{<:GenericMPSTensor}`: vector of site tensors -- `ALs::AbstractVector{<:GenericMPSTensor}`: vector of left-gauged site tensors -- `C₀::MPSBondTensor`: initial gauge tensor - -- `f::Function=rand`: initializer function for tensor data -- `eltype::Type{<:Number}=ComplexF64`: scalar type of tensors - -- `physicalspaces::AbstractVector{<:Union{S, CompositeSpace{S}}`: list of physical spaces -- `virtualspaces::AbstractVector{<:Union{S, CompositeSpace{S}}`: list of virtual spaces - -### Keywords -- `tol`: gauge fixing tolerance -- `maxiter`: gauge fixing maximum iterations + +Recommended ways to construct an infinite (periodic unit-cell) MPS are: + +- Using an MPS manifold of spaces + + ```julia + rand([rng], [T], manifold::InfiniteMPSManifold; tol, maxiter) + randn([rng], [T], manifold::InfiniteMPSManifold; tol, maxiter) + ``` + + Build an [`InfiniteMPSManifold`](@ref) with physical spaces and (maximal) virtual spaces for the unit cell. + +- From unit-cell site tensors + + ```julia + InfiniteMPS(As::AbstractVector{<:GenericMPSTensor}; tol, maxiter) + ``` + + Takes a vector `As` of (full-rank preferred) site tensors defining one unit cell. + The tensors are gauge-fixed (left/right) and internal bond tensors `C` are produced. + If any tensor isn't full rank, a warning is emitted and `makefullrank!` is applied. + +- From left- or right- gauged tensors and an initial gauge tensors + + ```julia + InfiniteMPS(ALs::AbstractVector{<:GenericMPSTensor}, C₀::MPSBondTensor; tol, maxiter) + InfiniteMPS(C₀::MPSBondTensor, ARs::AbstractVector{<:GenericMPSTensor}; tol, maxiter) + ``` + + Starts from gauged tensors `ALs` or `ARs` and an initial center bond `C₀` and completes the other gauge. + +### Keywords (passed to [`gaugefix!`](@ref)) +- `tol`: convergence tolerance for gauge fixing. +- `maxiter`: maximum iterations in gauge fixing. +- Additional keyword arguments accepted by `gaugefix!` (e.g. `order = :L | :R`). + +### Examples +```julia +using MPSKit, TensorKit +ps = PeriodicVector([ℂ^2, ℂ^2, ℂ^2]) # physical spaces, 3-site unit cell +vs = PeriodicVector([ℂ^4, ℂ^8, ℂ^4]) # maximal virtual spaces +m = InfiniteMPSManifold(ps, vs) +ψ = rand(ComplexF64, m; tol=1e-10) + +# Construct from pre-built site tensors +As = map(i -> rand(ComplexF64, m[i]), 1:length(m)) +ψ2 = InfiniteMPS(As; tol=1e-10) +``` + +!!! warning "Deprecated constructors" + Older constructors like `InfiniteMPS(pspaces, Dspaces)` or with `[f, T]` are deprecated. + Use `rand`/`randn` with an [`InfiniteMPSManifold`](@ref) instead. """ struct InfiniteMPS{A <: GenericMPSTensor, B <: MPSBondTensor} <: AbstractMPS AL::PeriodicVector{A} AR::PeriodicVector{A} C::PeriodicVector{B} AC::PeriodicVector{A} + + # constructor from data function InfiniteMPS{A, B}( AL::PeriodicVector{A}, AR::PeriodicVector{A}, C::PeriodicVector{B}, AC::PeriodicVector{A} = AL .* C @@ -62,50 +92,67 @@ struct InfiniteMPS{A <: GenericMPSTensor, B <: MPSBondTensor} <: AbstractMPS return new{A, B}(AL, AR, C, AC) end - function InfiniteMPS( - AL::PeriodicVector{A}, AR::PeriodicVector{A}, - C::PeriodicVector{B}, AC::PeriodicVector{A} = AL .* C - ) where {A <: GenericMPSTensor, B <: MPSBondTensor} - # verify lengths are compatible - L = length(AL) - L == length(AR) == length(C) == length(AC) || - throw(ArgumentError("incompatible lengths of AL, AR, C, and AC")) - # verify tensors are compatible - spacetype(A) == spacetype(B) || - throw(SpaceMismatch("incompatible space types of AL and C")) - for i in 1:L - N = numind(AL[i]) - N == numind(AR[i]) == numind(AC[i]) || - throw(SpaceMismatch("incompatible spaces at site $i")) - - # verify that the physical spaces are compatible - phys_ind = 2:(N - 1) - all( - space.(Ref(AL[i]), phys_ind) .== space.(Ref(AR[i]), phys_ind) .== - space.(Ref(AC[i]), phys_ind) - ) || - throw(SpaceMismatch("incompatible physical spaces at site $i")) - - # verify that the virtual spaces are compatible - space(AL[i], 1) == dual(space(AL[i - 1], N)) && - space(AR[i], 1) == dual(space(AR[i - 1], N)) && - space(AC[i], 1) == space(AL[i], 1) && - space(AC[i], N) == space(AR[i], N) && - space(C[i], 1) == dual(space(AL[i], N)) && - space(AR[i], 1) == dual(space(C[i - 1], 2)) || - throw(SpaceMismatch("incompatible virtual spaces at site $i")) - # verify that the spaces are non-zero - dim(space(AL[i])) > 0 && dim(space(C[i])) > 0 || - @warn "no fusion channels available at site $i" - end + # undef constructors + function InfiniteMPS{A, B}( + ::UndefInitializer, L::Integer + ) where {A <: GenericMPSTensor, B <: MPSBondTensor} + AL = PeriodicVector(Vector{A}(undef, L)) + AR = PeriodicVector(Vector{A}(undef, L)) + C = PeriodicVector(Vector{B}(undef, L)) + AC = PeriodicVector(Vector{A}(undef, L)) return new{A, B}(AL, AR, C, AC) end + end -#=========================================================================================== -Constructors -===========================================================================================# +function InfiniteMPS{A, B}( + ::UndefInitializer, manifold::InfiniteMPSManifold + ) where {A <: GenericMPSTensor, B <: MPSBondTensor} + psi = InfiniteMPS{A, B}(undef, length(manifold)) + for i in 1:length(manifold) + V = manifold[i] + psi.AL[i] = A(undef, V) + psi.AR[i] = A(undef, V) + psi.AC[i] = A(undef, V) + Vr = right_virtualspace(manifold, i) + psi.C[i] = B(undef, Vr ← Vr) + end + return psi +end +function InfiniteMPS{A}( + ::UndefInitializer, manifold::InfiniteMPSManifold + ) where {A <: GenericMPSTensor} + B = tensormaptype(spacetype(A), 1, 1, storagetype(A)) + return InfiniteMPS{A, B}(undef, manifold) + +end +function InfiniteMPS( + AL::PeriodicVector{A}, AR::PeriodicVector{A}, + C::PeriodicVector{B}, AC::PeriodicVector{A} = AL .* C + ) where {A <: GenericMPSTensor, B <: MPSBondTensor} + (L = length(AL)) == length(AR) == length(C) == length(AC) || + throw(ArgumentError("incompatible lengths of AL, AR, C, and AC")) + spacetype(A) == spacetype(B) || + throw(SpaceMismatch("incompatible space types of AL and C")) + + for i in 1:L + numind(AL[i]) == numind(AR[i]) == numind(AC[i]) || + throw(SpaceMismatch("incompatible spaces at site $i")) + space(AL[i]) == space(AR[i]) == space(AC[i]) || + throw(SpaceMismatch("incompatible spaces at site $i")) + domain(C[i]) == codomain(C[i]) || + throw(SpaceMismatch("Non-square C at site $i")) + right_virtualspace(AL[i]) == left_virtualspace(AR[i + 1]) || + throw(SpaceMismatch("incompatible spaces between site $i and site $(i + 1)")) + + # verify that the spaces are non-zero + (dim(space(AL[i])) > 0 && dim(space(C[i])) > 0) || + @warn "no fusion channels available at site $i" + end + + return InfiniteMPS{A, B}(AL, AR, C, AC) +end function InfiniteMPS( AL::AbstractVector{A}, AR::AbstractVector{A}, C::AbstractVector{B}, @@ -117,40 +164,14 @@ function InfiniteMPS( ) end -function InfiniteMPS( - pspaces::AbstractVector{S}, Dspaces::AbstractVector{S}; - kwargs... - ) where {S <: IndexSpace} - return InfiniteMPS(MPSTensor.(pspaces, circshift(Dspaces, 1), Dspaces); kwargs...) -end -function InfiniteMPS( - f, elt::Type{<:Number}, pspaces::AbstractVector{S}, Dspaces::AbstractVector{S}; - kwargs... - ) where {S <: IndexSpace} - return InfiniteMPS( - MPSTensor.(f, elt, pspaces, circshift(Dspaces, 1), Dspaces); - kwargs... - ) -end -InfiniteMPS(d::S, D::S) where {S <: Union{Int, <:IndexSpace}} = InfiniteMPS([d], [D]) -function InfiniteMPS( - f, elt::Type{<:Number}, d::S, D::S - ) where {S <: Union{Int, <:IndexSpace}} - return InfiniteMPS(f, elt, [d], [D]) -end -function InfiniteMPS(ds::AbstractVector{Int}, Ds::AbstractVector{Int}) - return InfiniteMPS(ComplexSpace.(ds), ComplexSpace.(Ds)) -end -function InfiniteMPS( - f, elt::Type{<:Number}, ds::AbstractVector{Int}, Ds::AbstractVector{Int}, kwargs... - ) - return InfiniteMPS(f, elt, ComplexSpace.(ds), ComplexSpace.(Ds); kwargs...) -end +#=========================================================================================== +Constructors +===========================================================================================# function InfiniteMPS(A::AbstractVector{<:GenericMPSTensor}; kwargs...) # check spaces - leftvspaces = circshift(_firstspace.(A), -1) - rightvspaces = conj.(_lastspace.(A)) + leftvspaces = circshift(left_virtualspace.(A), -1) + rightvspaces = right_virtualspace.(A) isnothing(findfirst(leftvspaces .!= rightvspaces)) || throw(SpaceMismatch("incompatible virtual spaces $leftvspaces and $rightvspaces")) @@ -160,24 +181,12 @@ function InfiniteMPS(A::AbstractVector{<:GenericMPSTensor}; kwargs...) @warn "Constructing an MPS from tensors that are not full rank" makefullrank!(A_copy) - AR = A_copy - - leftvspaces = circshift(_firstspace.(AR), -1) - rightvspaces = conj.(_lastspace.(AR)) - isnothing(findfirst(leftvspaces .!= rightvspaces)) || - throw(SpaceMismatch("incompatible virtual spaces $leftvspaces and $rightvspaces")) - - # initial guess for the gauge tensor - V = _firstspace(A_copy[1]) - C₀ = isomorphism(storagetype(eltype(A_copy)), V, V) - - # initialize tensor storage - AL = similar.(AR) - AC = similar.(AR) - C = similar(AR, typeof(C₀)) - ψ = InfiniteMPS{eltype(AL), eltype(C)}(AL, AR, C, AC) + manifold = InfiniteMPSManifold(A_copy) + ψ = InfiniteMPS{eltype(A)}(undef, manifold) # gaugefix the MPS + V = left_virtualspace(ψ, 1) + C₀ = isomorphism(storagetype(eltype(A_copy)), V, V) gaugefix!(ψ, A_copy, C₀; kwargs...) mul!.(ψ.AC, ψ.AL, ψ.C) @@ -185,18 +194,11 @@ function InfiniteMPS(A::AbstractVector{<:GenericMPSTensor}; kwargs...) end function InfiniteMPS(AL::AbstractVector{<:GenericMPSTensor}, C₀::MPSBondTensor; kwargs...) - AL = PeriodicArray(copy.(AL)) - - all(isfullrank, AL) || - @warn "Constructing an MPS from tensors that are not full rank" + AL = PeriodicArray(AL) - # initialize tensor storage - AC = similar.(AL) - AR = similar.(AL) - T = TensorOperations.promote_contract(scalartype(AL), scalartype(C₀)) - TC = TensorOperations.tensoradd_type(T, C₀, ((1,), (2,)), false) - C = similar(AR, TC) - ψ = InfiniteMPS{eltype(AL), TC}(AL, AR, C, AC) + all(isfullrank, AL) || @error "Constructing an MPS from tensors that are not full rank" + ψ = InfiniteMPS{eltype(AL)}(undef, InfiniteMPSManifold(AL)) + ψ.AL .= AL # gaugefix the MPS gaugefix!(ψ, AL, C₀; order = :R, kwargs...) @@ -206,15 +208,11 @@ function InfiniteMPS(AL::AbstractVector{<:GenericMPSTensor}, C₀::MPSBondTensor end function InfiniteMPS(C₀::MPSBondTensor, AR::AbstractVector{<:GenericMPSTensor}; kwargs...) - AR = PeriodicArray(copy.(AR)) + AR = PeriodicArray(AR) - # initialize tensor storage - AC = similar.(AR) - AL = similar.(AR) - T = TensorOperations.promote_contract(eltype(AR), eltype(C₀)) - TC = TensorOperations.tensoradd_type(T, C₀, ((1,), (2,)), false) - C = similar(AR, TC) - ψ = InfiniteMPS{eltype(AL), TC}(AL, AR, C, AC) + all(isfullrank, AR) || @error "Constructing an MPS from tensors that are not full rank" + ψ = InfiniteMPS{eltype(AR)}(undef, InfiniteMPSManifold(AR)) + ψ.AR .= AR # gaugefix the MPS gaugefix!(ψ, AR, C₀; order = :L, kwargs...) @@ -223,6 +221,49 @@ function InfiniteMPS(C₀::MPSBondTensor, AR::AbstractVector{<:GenericMPSTensor} return ψ end +for randfun in (:rand, :randn) + randfun! = Symbol(randfun, :!) + @eval function Random.$randfun(rng::Random.AbstractRNG, T::Type, manifold::InfiniteMPSManifold; kwargs...) + As = map(i -> $randfun(rng, T, manifold[i]), 1:length(manifold)) + return InfiniteMPS(As; kwargs...) + end + @eval function Random.$randfun!(rng::Random.AbstractRNG, mps::InfiniteMPS) + foreach(Base.Fix1(rng, $randfun!), mps.AC) + C₀ = $randfun!(rng, mps.C[0]) + gaugefix!(mps, mps.AC, C₀) + mul!.(mps.AC, mps.AL, mps.C) + return mps + end +end + +# Deprecate old constructor syntaxes +# ---------------------------------- +Base.@deprecate( + InfiniteMPS(pspace::S, Dspace::S; kwargs...) where {S <: IndexSpace}, + rand(InfiniteMPSManifold(pspace, Dspace); kwargs...) +) +Base.@deprecate( + InfiniteMPS(pspaces::AbstractVector{S}, Dspaces::AbstractVector{S}; kwargs...) where {S <: IndexSpace}, + rand(InfiniteMPSManifold(pspaces, Dspaces); kwargs...) +) +Base.@deprecate( + InfiniteMPS(f, T::Type, pspace::S, Dspace::S; kwargs...) where {S <: IndexSpace}, + f(T, InfiniteMPSManifold(pspace, Dspace); kwargs...) +) +Base.@deprecate( + InfiniteMPS(f, T::Type, pspaces::AbstractVector{S}, Dspaces::AbstractVector{S}; kwargs...) where {S <: IndexSpace}, + f(T, InfiniteMPSManifold(pspaces, Dspaces); kwargs...) +) + +Base.@deprecate( + InfiniteMPS(ds::AbstractVector{Int}, Ds::AbstractVector{Int}), + rand(InfiniteMPSManifold(ComplexSpace.(ds), ComplexSpace.(Ds))) +) +Base.@deprecate( + InfiniteMPS(f, T::Type, ds::AbstractVector{Int}, Ds::AbstractVector{Int}), + f(T, InfiniteMPSManifold(ComplexSpace.(ds), ComplexSpace.(Ds))) +) + #=========================================================================================== Utility ===========================================================================================#