YASTN: Yet another symmetric tensor networks; A Python library for abelian symmetric tensor network calculations

1) Paper: interesting and instructive examples
2) Paper: good balance between conciseness and details
3) Software: one of the few packages with full support for GPUs and automatic differentiation
4) Software: one of the few public packages for simulations with (infinite and possibly fermionic) PEPS
5) Software: extensive set of examples and models (lattices and Hamiltonians)

1) Some statements or code samples in the paper would benefit from a more extensive explanation (see report)
2) Somewhat unclear code structure regarding the interplay with peps-torch (see report)
3) No benchmark comparisons to 'competing' software packages
4) Grammar (use of articles) can be improved
5) Examples focus on PEPS (on CPUs) only

This paper presents the open-source Python package YASTN for tensor network simulations of quantum many-body systems. It has built-in support for matrix product states (MPS) and, notably as compared to some of the other main packages in the field, extensive support for projected entangled-pair states (PEPS). At the core is a tensor type that can represent tensors with abelian symmetries and for which the data can be stored using both Numpy and PyTorch backends, where the latter enables support for GPUs and automatic differentation. These features, together with the extensive PEPS support, set it apart from its main competitor in the Python universe, namely TenPy (which arguably has more support for different MPS algorithms and extensions).

The manuscript itself reads well (aside from an occasional missing or ill-used article "a" / "the"). It concisely introduces the design principles, explaining in a good amount of detail how abelian symmetries constrain the number of nonzero tensor entries and how to specify such symmetries and create and manipulate tensors with them. I have a few small questions or clarification requets on this section which are listed below.

About the second half of the manuscript is spend on a number of impressive examples that are instructive regarding the implications of using abelian symmetries in PEPS simulations. However, given the extensive set of features of this package, the benchmarks do only cover a very small fraction of it, and in particular do not cover aspects such as
* CPU versus GPU with and without symmetries, and the effect of the symmetry overhead in that case
* comparison to other packages such as TenPy and ITensors.jl, for simple tensor network contraction or for ... next bullet point
* anything related to MPS (possibly in relation to the previous bullet point, as the main competitors have a strong MPS emphasis)
* memory and runtime cost of automatic differentiation

I will list a number of questions and suggestions that might result in changes here, but they are optional and could also be addressed in the response:

1) Around line104: it is perhaps worth mentioning that the charges for $Z_n$ are elements of $\mathbb{Z}_n$, because in those cases, the angle theta in the symmetry is restricted to take values that are integer multiples of $2\pi/N$. Also, I don't quite understand why the different symbol $Z_n$ and $\mathbb{Z}_n$ are used, but if the first one is supposed to refer to the multiplicative group, it is probably more correct to denote it as the cyclic group $C_n$.

2) Around line 113: it is stated that boldface is used to underline the vector nature, but it seems that the SciPost style makes all math in boldface, so that also scalar values such as $s$, individual dimensions per sector $D_\rho$, and esssentially all other math is typeset in boldface.

3) The final paragraph of section 2.1: can you elaborate on the the benefits of the lazy initalisation and the allocation of nonzero blocks only. Is it often the case that blocks which are allowed by the symmetry, turn out to be zero anyway? What happens if they receive nonzero data afterwards? Does this require to allocate a new tensor or to copy the data?

4) The code examples in lines 155-156 and 169-170 could perhaps use a bit more explanation, especially for non-Python users (although they might not be the typical audience of this paper). I was in particular wondering about what exactly the 'axes' argument in the 'tensordot' call specifies, and how the indices of the resulting tensor are sorted/defined. Another question is about the recording of the original structure of fusing legs: what happens if different fusion steps are concatenated, e.g. first 1 with 2 and 3 with 3, then the resulting 1 with the resulting 2. Are these intermediate steps recorded. Will 'unfuse_legs' take one step back, or immediately jump back to the original structure?

5) Given that wall times are reported in seconds, would it be possible to provide some details on the hardware on which these simulations were ran?

6) More benchmarks can provide further interesting insight, and as listed in the report section, there are quite a few benchmarks that could be conceivable. I leave it to the authors to decide what they believe could be useful additions here.

7) I believe it would be good if the interplay/interdependency with peps-torch could be clarified. It is sometimes listed as an independent package (introduction, line 67, as well as conclusion, line 348), but at the same time seems closely tied to YASTN and is mentioned as part of it for doing the PEPS optimisations (e.g. conclusion, line 347).

Publish (easily meets expectations and criteria for this Journal; among top 50%)

1) easy to install
2) clear code structure
3) extensive test suite
4) pytorch backend with support of GPU acceleration and auto differentation
5) interface quite similar to well-known numpy

1) The benchmarks do not compare to other implementations, but only check the expected scaling within the package. It might be interesting to compare the overhead implied by charge conservation to other libraries like TeNPy or ITensor.
If the authors are interested, I'd be happy to help setting up a comparison with TeNPy similar to https://github.com/tenpy/tenpy_benchmarks and https://itensor.github.io/ITensorBenchmarks.jl/dev/tenpy_itensor/index.html

2) The documentation includes installation instructions, examples and benchmarks as required by the acceptance criteria. However, it could be a bit more extensive in some places. E.g. https://yastn.github.io/yastn/examples/fpeps/ctmrg.html mentions a comparison to the Onsager solution. While this is indeed done in the corresponding tests/peps/test_ctmrg.py file, this comparison is not included in the docs.

YASTN is a solid python library for tensor network calculations based on the widely adopted numpy and pytorch as backends. It extends those packages by implementing a tensor class exploiting a block-diagonal structure emerging from abelian symmetries, and the benchmarks in the manuscript demonstrate that using the symmetries can lead to significant speedups. The tensor interface is quite similiar to the well-known numpy package; in contrast to e.g. TeNPy it does not introduce labels, or make the legs globally unique like iTensor.
High-level tensor network algorithms are provided for (finite and infinite) PEPS (the focus of the library) and finite MPS.

While YASTN indeed joins a row of other tensor network libraries as implied by the name, it has sufficiently different focus than the other packages and is an independent implementation valuable for cross-checks. Further, the demand has already been demonstrated by several publications employing it for state-of-the-art calculations. I hence think that the acceptance criteria for Scipost Physcis Codebase are met and recommend a publication.

Optional: benchmark comparsion to other tensor network library, see point 1) under weaknesses.

Publish (meets expectations and criteria for this Journal)