Fix Multiline (MultilineMPS/MultilineMPO) correctness and API inconsistencies - #508
Fix Multiline (MultilineMPS/MultilineMPO) correctness and API inconsistencies#508borisdevos wants to merge 32 commits into
Multiline (MultilineMPS/MultilineMPO) correctness and API inconsistencies#508Conversation
…referring to the 2d shape
…o bd/multiline
|
After having discussed with @leburgel, I realised I didn't fully appreciate/understand what multiline actually represents, so recent changes rectify some mistakes I made, and tried to explain these niche things more cleanly in the docs/docstrings.
|
| - `Multiline` (and therefore `MultilineMPS`/`MultilineMPO`) now consistently treats | ||
| `length`/`eltype`/`iterate`/`m[i]` as referring to the individual lines it stores | ||
| (`length(m) == nrows`), while `size`/`axes`/`eachindex` refer to the `(nrows, ncols)` lattice | ||
| shape. |
There was a problem hiding this comment.
I think this is the change I have the most questions about. I like that you state that it is consistent, even though it is clearly inconsistent between most of these methods 😆 (length(x) == prod(size(x)) is now explicitly broken).
Am I understanding this correctly that you are effectively claiming this is interpreted as a vector of lines, with a little bit of syntactic sugar to make m[i, j] work and then have size follow the syntactic sugar?
With that in mind, I think most of the code is actually not calling m[i,j] at all, since for MPS that does not really mean anything (you'd need the .AC or something anyways), and for MPO I think we are mostly doing things by first unpacking the rows and only then indexing the sites.
Would it also work to actually pull this through all the way, just remove the lattice shape from the equation entirely?
There was a problem hiding this comment.
It's consistent for me in the sense that when working with a multiline MPS what really makes it different is the layers of MPSs, so it made sense to me that length and eltype should refer to the layers. Also, m[i, j] is actually not for iMPS because we don't even allow slicing it (as you say, it doesn't make sense), though it is supported for iMPOs.
So yeah in a nutshell, I am indeed trying to make multiline behave like an AbstractVector at the level of iterating and index slicing correctly over its length, but it's not fully a subtype of that because you still want 2d cartesian indexing over something like m.AL. The initial motivation for this was that some map (among others) on a multiline MPS would previously make an nrows*ncols vector with only nrows defined entries.
Do you mean by removing the lattice shape to just treat multiline as a vector of lines? This is certainly doable, though every size(mpsi, 2) would need to be replaced by length(mpsi[1]), among others which rely on the 2d structure. If you think this is better, then might as well make Multiline{T} <: AbstractVector{T}, right?
There was a problem hiding this comment.
I think if you want to change it and make it consistent, might as well make it behave like a vector of lines, indeed, and possibly just add a new name for size(mps, 2). It just really looks iffy to me to have both cartesian 2d indexing as well as linear indexing but have both of them not compatible, and I'm definitely fine with keeping m.AL alive, as that doesn't interfere with this, since it outputs a separate ALView, for which we can then separately decide to make it behave like a Matrix instead of a Vector.
TLDR, I think technically the current and previous approach are piracy, in the sense that we are stealing the AbstractArray interface and then incorrectly implementing it :)
| # Dominant eigenvalues | ||
| # -------------------- | ||
| """ | ||
| dominant_eigenvalue(ψ, O, [environments]) -> λ |
There was a problem hiding this comment.
Can you maybe explain a bit better what is going on here? I am not really following most of the discussion here to be honest. Why is this not ⟨ϕ|O|ϕ⟩?
From what I can tell it really is that, in the sense that because you take the product over the transfer of all of the rows, you get exactly the overlap after taking it through an entire unit cell, so I don't really understand what is going on here.
I think I might get that it is not an expectation value in the sense that it is not the expectation value of an operator in the statmech sense, for which indeed you'd use dot(psi, O, phi) instead, but if expectation_value(psi, O) = dot(psi, O, psi) by definition, I don't immediately see the benefit of making this a separate function.
There was a problem hiding this comment.
The motivation for this was indeed semantics, trying to stress the fact that it's a stat-mech method, and there the way to reason this eigenvalue is really from applying transfer operators onto a state which lives on one side of the partition function, say the top. To really consider an overlap, which is the standard way of interpreting an expectation value, you'd need to sandwich your transfer operators from both sides, and the fixed point in the bottom is generally not that from the top. That's kind of what I tried explaining in the note in this docstring.
Really, from an identity point of view, you're right, but it's simply not what's being calculated, and I wanted to avoid a collision of meanings behind the same "expectation value" by just giving meaning to the more correct "dominant eigenvalue"
There was a problem hiding this comment.
I'm still confused by this, it really is what is being calculated? I think my main point here is that this seems like you are introducing a lot of extra code mostly for pedantics, which I am sympathetic to but I fail to see a lot of benefits here. In general, most of our code uses the operator type to distinguish the meanings, and since expectation_value means dot / norm, which is precisely what we want in both cases, I am not really sure why this is so confusing?
| # ------------ | ||
| const MultilineMPS = Multiline{<:InfiniteMPS} | ||
| #TODO: add support for finite MPS | ||
| const _MPSs = Union{InfiniteMPS, FiniteMPS} |
There was a problem hiding this comment.
No, because I want to leave out window MPS
| return checkbounds(Bool, _lineview(psi), b) | ||
| end | ||
| _lineview(psi::ACView) = ACView(first(parent(psi.parent))) | ||
| _lineview(psi::ALView) = ALView(first(parent(psi.parent))) | ||
| _lineview(psi::ARView) = ARView(first(parent(psi.parent))) | ||
| _lineview(psi::CView) = CView(first(parent(psi.parent))) |
There was a problem hiding this comment.
| return checkbounds(Bool, _lineview(psi), b) | |
| end | |
| _lineview(psi::ACView) = ACView(first(parent(psi.parent))) | |
| _lineview(psi::ALView) = ALView(first(parent(psi.parent))) | |
| _lineview(psi::ARView) = ARView(first(parent(psi.parent))) | |
| _lineview(psi::CView) = CView(first(parent(psi.parent))) | |
| return checkbounds(Bool, parent(psi), a) && | |
| checkbounds(Bool, @inbounds psi[a], b) | |
| end |
Does this also work?
|
|
||
| # the row index is always periodic | ||
| # for infinite lines the column index is periodic as well | ||
| Base.checkbounds(::Type{Bool}, ::MultilineOrthoView{<:Multiline{<:InfiniteMPS}}, a, b) = true |
There was a problem hiding this comment.
See comment below, which then might make this obsolete?
|
|
||
| VectorInterface.add!!(x::Multiline, y::Multiline, α::Number, β::Number) = add!(x, y, α, β) | ||
|
|
||
| # FIXME? is it intentional that a nontrivial multilinemps of normalised rows never has norm 1? |
There was a problem hiding this comment.
I'm actually wondering if this should be prod instead of sum now...
| TensorKit.storagetype(::Type{Multiline{T}}) where {T} = storagetype(T) | ||
| for ftype in (:spacetype, :sectortype, :storagetype) | ||
| @eval TensorKit.$ftype(::Type{Multiline{T}}) where {T} = $ftype(T) | ||
| @eval TensorKit.$ftype(m::Multiline) = $ftype(typeof(m)) |
There was a problem hiding this comment.
I somehow was under the impression that the latest versions of TensorKit already do this forwarding, is that wrong?
[Relevant edit deeper in conversation]
I found a cluster of bugs that all traced back to the same root cause, namely
Multiline{T}was inconsistent about what it is. Sometimes it was treated as a 2D array indexed by(row, col), and sometimes as a 1D sequence ofT-typed lines, this while different parts of the code picked different conventions without checking consistency. The real bugs were:MultilineMPO * MultilineMPShad never worked (straight up method error), and fixing that trivially still left both a shape bug (lengthcountingrows*colswhilemap/collectonly ever producedrowselements) and a silently wrong row-shift bug (fidelity 0where it should be1).expectation_value(::MultilineMPS, ::MultilineMPO, ...)had a fallback method that silently returned something meaningless (prodinstead ofsum, no row shift, discardedenvswhich was noted through a todo) for anyMultilineMPOline type it wasn't explicitly guarded against. Related to this, there was half-support for Hamiltonian or finite MPO multilines even though they're meaningless. Also the environments weren't being correctly converted, so they were being recalculated every row.isfinite,changebonds(::MultilineMPO, ::SvdCut), andaxes(m, 2)all threw outright.A summary of the changes I made:
"Semantics":
Multiline{T}now consistently treatslength/eltype/iterate/m[i]as referring to the lines (matching the periodic rows), whilesize/axes/eachindexrefer to the(nrows, ncols)lattice shape. These two views are intentionally different (length(m) != prod(size(m))). I also documented this in docstrings.expectation_valuefallback. Changed the infinite MPS-MPO expectation value to correctly account for environments.isfinite,changebonds(::MultilineMPO/::MultilineMPS, ::SvdCut), andaxes(m, i).spacetype/sectortype/storagetypefor multiline. I needed this for the custom show (see below).Construction restrictions:
MultilineMPOonly acceptsInfiniteMPOlines now, enforced at constructor level. This is deliberate, sinceMultilineis also used via converts in environments, quasiparticles etc.AbstractMatrixconstructor that silently built finite-lineMultilineMPOs, and the dead code inorthoview.jlthat existed only to support finite MPS multilines.Display:
MultilineMPS/MultilineMPOcustom show: added asummarywhich looks at the 2d shape, and ashowthat renders each row using that row's own existingshowmethod. Here below an example of how it looks now:Here the multiline MPS:
Details
And here a multiline MPO:
Details
Tests:
test/operators/multilinempo.jl: exact-fidelity tests for the row-shift convention in bothMultilineMPO * MultilineMPSandMultilineMPO * MultilineMPO, using a permutation MPO whose fixed point is known exactly.test/algorithms/statmech.jl: multi-rowexpectation_valueregression test. Also a test that the Hamiltonian loophole above still correctly throws.test/states/multilinemps.jl: added instance-levelisfinite(ψ)(only the type-level form was covered before).Documentation:
docs/src/man/states.md: added the row-shift convention, what's currently supported (infinite lines only), and subtleties.docs/src/man/operators.md: newMultilineMPOsection covering the same convention, and emphasising the MPO*MPO non-composability caveat.Multiline,MultilineMPO,MultilineMPS, andexpectation_value.Some things I didn't address because I'm not sure what to do with them:
MultilineMPS, simply becauseInfiniteMPSitself doesn't implementzerovector/scale.norm(ψ) == sqrt(nrows)for aMultilineMPSof individually-normalized rows is unchanged. This goes throughdot/innersum over rows. I don't know if this is what you want.Checklist
julia --project=test test/runtests.jl, or the relevant subset)docs/src/)[Unreleased]indocs/src/changelog.md, if this PR is user-facing (new feature, behavior change, bug fix, deprecation, or removal)