Add CTMRG algorithm for layers of PEPO tensors
#184
Conversation
Codecov ReportAttention: Patch coverage is
... and 2 files with indirect coverage changes 🚀 New features to boost your workflow:
|
this includes a small bug fix
Co-authored-by: Lander Burgelman <[email protected]>
Co-authored-by: Lander Burgelman <[email protected]>
Co-authored-by: Lander Burgelman <[email protected]>
Co-authored-by: Lander Burgelman <[email protected]>
Co-authored-by: Lander Burgelman <[email protected]>
This includes a `_dagger` function on a `PEPOTensor`, which is probably far from ideal
This changes the name from PEPOLayers to PEPOTrace. Other names of variables and files are changed accordingly. This also includes a small bug fix on the scalartype in the testfile
|
I think everything should be fine now. And if @lkdvos is okay with the naming convention, I think we can merge this. |
Co-authored-by: Lander Burgelman <[email protected]>
Co-authored-by: Lander Burgelman <[email protected]>
Co-authored-by: Lander Burgelman <[email protected]>
Co-authored-by: Lander Burgelman <[email protected]>
There was a problem hiding this comment.
This looks good to me, @lkdvos what do you think? The functionality is a bit experimental, but I think it's mostly well defined what one would do with a stack of PEPOs (similar to the density operator support in MPSKit I think). There's no changes to the public interface and the new features are tested, so good enough for me at this point. Once the things you'll use is it for are better established, it would be great to also add an example at a later point @sanderdemeyer.
As a side not, but since it's again relevant here: I think at some point we should add more 'public' support for daggering and transposing transfer operators, both 3D ones like Sander uses here and 2D ones like the MPSKit MPOs. Now every time someone needs a left (i.e. bottom in our view) fixed point people usually do some shady manual intervention. Even though it's inherently a bit ill-defined, I think a properly documented dagger that exactly states what it does (and does not do) would be a valuable addition to the public interface in the future.
There was a problem hiding this comment.
Looks more or less okay to me, left some small comments.
I am also not entirely sold on PEPOTrace, but I also cannot think of that many better options.
If this is only used for finite-temperature things I'd say to call it something based of Projected Entangled Pair Density Operator (PEPDO), in analogy to the Matrix Product Density Operator (MPDO), but I've only ever really seen the MPDO in the literature so I don't think it matters too much.
The finite-temperature stuff is indeed why I am adding this. In principle, this code is not restricted to this case and could be applied to PEPOs, so I prefer |
Change dagger to TensorKit.adjoint Fix docstring of adjoint Remove rrule for InfinitePEPO
include test on adjoint operator change test on mpotensor make the tests fermionic slightly modify the docstring of adjoint
|
Several comments from above are now included. To try to resolve the discussion about the naming convention, i.e. to make sure that what we're doing is truly the (mathematically rigorous) adjoint, I added a test to verify that To make sure that all the twists are correctly implemented, I also changed the test to the fermionic case. |
| adjoint(O::InfinitePEPO) | ||
|
|
||
| Create the adjoint of an InfinitePEPO. | ||
| With this definition, <Oψ₁, ψ₂> = <ψ₁, adjoint(O)ψ₂>, with Oψ the physical action of the PEPO O on the PEPS ψ. |
There was a problem hiding this comment.
I think this is a bit of a strange notation, since this seems to be closer to the mathematicians notation than the physicist notation (which is also what Julia follows).
Note that dot(v1, v2) takes the complex conjugate of the first argument, in accordance to the first argument being the bra, and dot(v1, A, v2) = dot(v1, A * v2), so I think it is more natural to say that adjoint(O) = O' is defined such that dot(psi, O, phi) = dot(psi, O * phi) = dot(O' * psi, phi), or in other words:
Given that the definition in terms of dot is actually really what we are doing, I'd also say that this might be a better notation to add to the docstring, since otherwise it is not even that clear what you mean with the notation. It also avoids questions about unicode and math etc.
| With this definition, <Oψ₁, ψ₂> = <ψ₁, adjoint(O)ψ₂>, with Oψ the physical action of the PEPO O on the PEPS ψ. | |
| This is defined such that `dot(psi, O, phi) == dot(psi, O * phi) == dot(O' * psi, phi)` for any two states `psi` and `phi`. |
There was a problem hiding this comment.
There is quite a lot of code in the global scope here, which effectively means that most of it is type-unstable. This isn't necessarily a big problem, but given that our tests are quite computationally demanding this might not be that great.
It could be improved by turning some of it into functions, but I'm also slightly confused about what the tests are actually doing. How come there are no twists in any of the fusing contractions?
| @testset "mpotensor for PEPOTraceSandwich" begin | ||
| # Test whether the mpotensor of a double layer of PEPOs is the same as in the case where the layers are fused | ||
| mpo = InfiniteSquareNetwork(map(PEPSKit.mpotensor, PEPSKit.unitcell(network))) | ||
| mpo_fused = InfiniteSquareNetwork( | ||
| map(PEPSKit.mpotensor, PEPSKit.unitcell(network_fused)) | ||
| ) | ||
| @test mpo ≈ mpo_fused | ||
| end |
There was a problem hiding this comment.
Would it not make more sense to check that these networks have the same values, rather than checking the tensors?
| nrm = PEPSKit._contract_site((1, 1), network, env) | ||
| nrm_fused = PEPSKit._contract_site((1, 1), network_fused, env_fused) |
There was a problem hiding this comment.
If this is the actual norm, why can't you call norm explicitly? Why do we have to manually divide by the norm, is this not all equivalent to expectation_value?
| projector_alg = projector_algs[1] # only use :halfinfinite for this test due to convergence issues | ||
| @testset "Test adjoint of an InfinitePEPO using $alg with $projector_alg" for alg in | ||
| ctm_styles | ||
| # Test the definition of the adjoint of an operator, i.e. <Oψ₁, ψ₂> = <ψ₁, adjoint(O)ψ₂>, with Oψ the physical action of the PEPO O on the PEPS ψ. |
There was a problem hiding this comment.
Same comment about this being the mathematicians notation for adjoints. Also, could you add the case where this is equal to the peps-pepo-peps contraction?
There was a problem hiding this comment.
only now noticed that this is probably done below, it might make more sense to combine all of these?
This builds upon #134, where the CTMRG algorithm now includes
InfiniteSquareNetworks without a top and bottom layer, where the top physical leg of the topPEPOTensoris contracted with the bottom physical leg of the bottomPEPOTensor. This is based on the typePEPOTraceSandwich(The equivalent ofPEPOSandwichin the previous PR).This PR can be summarised as
PEPOTraceSandwichInfiniteSquareNetworkofPEPOTraceSandwiches based on anInfinitePEPOdaggerfor anInfinitePEPOPEPOTraceSandwichmpotensorfunction forPEPSSandwichto bothPEPOSandwiches andPEPOTraceSandwiches.PEPOis compared to the case where this double layer is now squashed into a single layerPEPOon aPEPSis compared to the case where this application is squashed into a single layerPEPSAny suggestions on the implementation and naming of things in this PR is very much appreciated. What anyway needs to be done, is
daggermpotensorfunctions that now exist