###### Abstract

We show that the spectral set conjecture by Fuglede [6] holds in the setting
of cyclic groups of order $p^{n}q$, where $p$, $q$ are distinct primes and
$n\geq 1$. This means that a subset $E$ of such a group $G$ tiles the group by
translation ($G$ can be partitioned into translates of $E$) if and only if
there exists an orthogonal basis of $L^{2}(E)$ consisting of group characters.
The main ingredient of the present proof is the structure of vanishing sums of
roots of unity of order $N$, where $N$ has at most two prime divisors; the
extension of this proof to the case of cyclic groups of order $p^{n}q^{m}$
seems therefore feasible. The only previously known infinite family of cyclic
groups, for which Fuglede’s conjecture is verified in both directions, is that
of cyclic $p$-groups, i.e. $\mathbb{Z}_{p^{n}}$.

title = Fuglede’s conjecture on cyclic groups of order $p^{n}q$, author =
Romanos Diogenes Malikiosis and Mihail N. Kolountzakis, plaintextauthor =
Romanos Diogenes Malikiosis and Mihail N. Kolountzakis, plaintexttitle =
Fuglede’s conjecture on cyclic groups of order p**n q, keywords = Spectral
set, tiling in finite abelian groups, Fuglede’s conjecture, vanishing sums of
roots of unity, year=2017, number=12, received=7 December 2016, published=5
September 2017, doi=10.19086/da.2071,

[classification=text]

##  1 Introduction

Let $\Omega$ be a measurable subset of $\mathbb{R}^{n}$ of positive Lebesgue
measure. $\Omega$ is called _spectral_ , if it accepts an orthogonal basis of
exponentials, namely $e^{i\lambda\cdot x}$, where $\lambda$ ranges through
$\Lambda\subseteq\mathbb{R}^{n}$; the set $\Lambda$ is called the _spectrum_
of $\Omega$. On the other hand, $\Omega$ is called a _tile_ of
$\mathbb{R}^{n}$, if there is a set $T\subseteq\mathbb{R}^{n}$ such that
almost every point of $\mathbb{R}^{n}$ can be written uniquely as $\omega+t$,
where $\omega\in\Omega$ and $t\in T$; $T$ is hence called the _tiling
complement_ of $\Omega$.

Fuglede’s spectral set conjecture [6] asserts that $\Omega$ is spectral if and
only if $\Omega$ is a tile. For many years, there was a lot of positive
evidence towards the veracity of this conjecture; Fuglede himself proved this
conjecture when the spectrum or the tiling complement of $\Omega$ is a
lattice. Later on, the conjecture was proved for various families of subsets
by numerous authors, for example, for $2$-dimensional convex bodies [7],
unions of two intervals in $\mathbb{R}$ [11], etc.

Against the current of research on this subject, Tao disproved this conjecture
for dimensions $5$ and above [15], constructing a spectral set in $5$
dimensions that does not tile the space111Any counterexample in $n$ dimensions
could be extrapolated to a counterexample in all higher dimensions..
Subsequently, works by Kolountzakis and Matolcsi [9, 10, 14] and Farkas-
Matolcsi-Mora [5] showed that the conjecture is false in dimensions $3$ and
above in both directions, leaving it open for $\mathbb{R}$ and
$\mathbb{R}^{2}$.

This conjecture can be naturally stated for other groups, for example
$\mathbb{Z}$ or any locally compact abelian group. These cases are not only
interesting on their own, but they also have connections to the original case.
In his disproof of the $5$-dimensional case, Tao constructed a spectral set in
$\mathbb{Z}_{3}^{5}$ containing $6$ elements, hence not a tile, as the
cardinality of any tile of a finite abelian group divides the order of the
group; then he lifted this counterexample to $\mathbb{R}^{5}$. So was done in
the disproof of the other direction.

Some examples of note, where Fuglede’s conjecture holds, include finite cyclic
$p$-groups [12, 4], $\mathbb{Z}_{p}\times\mathbb{Z}_{p}$ [8], and
$\mathbb{Q}_{p}$ [3, 4], the field of $p$-adic numbers. Borrowing the notation
from [2], we write S-T$(G)$, respectively T-S$(G)$, if the
Spectral$\Rightarrow$Tile direction, respectively Tile$\Rightarrow$Spectral,
holds in $G$; when we put $G=\mathbb{Z}_{N}$, it is understood that the
statement holds _for all_ $N$. The connection between the conjecture on
$\mathbb{R}$ and on finite cyclic groups or $\mathbb{Z}$ is summarized below
[2]:

$\text{{\bf T-S}}(\mathbb{R})\Longleftrightarrow\text{{\bf
T-S}}(\mathbb{Z})\Longleftrightarrow\text{{\bf T-S}}(\mathbb{Z}_{N}),$

and

$\text{{\bf S-T}}(\mathbb{R})\Longrightarrow\text{{\bf
S-T}}(\mathbb{Z})\Longrightarrow\text{{\bf S-T}}(\mathbb{Z}_{N}).$

According to the above connections, a counterexample in a finite cyclic group
can be lifted to a counterexample in $\mathbb{R}$; on the other hand, if the
conjecture were true for every cyclic group and $\mathbb{Z}$, this would hold
no meaning for the original conjecture in $\mathbb{R}$, unless it were proven
that every spectral set in $\mathbb{R}$ has a rational spectrum. We may ask
nevertheless, to which extent is Fuglede’s conjecture true for finite abelian
groups, or even cyclic ones. Surprisingly, not much is known for cyclic
groups, apart from cyclic $p$-groups, i.e. $\mathbb{Z}_{p^{n}}$ for $p$ prime.
The direction Tile$\Rightarrow$Spectral is known also for cyclic groups of
order $p^{n}q^{m}$, for $p$, $q$ distinct primes [1, 12]; see Section 3 below.

######  Theorem 1.1 ([1, 12]).

Let $A\subseteq\mathbb{Z}_{N}$ be a set that tiles $\mathbb{Z}_{N}$ by
translations, where $N=p^{n}q^{m}$ for distinct primes $p$, $q$. Then $A$ is
spectral.

The novel contribution of our present work is the proof of
Spectral$\Rightarrow$Tile direction for cyclic groups of order $N=p^{n}q$,
thus establishing the veracity of Fuglede’s conjecture in this setting.

######  Theorem 1.2.

Let $A\subseteq\mathbb{Z}_{N}$ be a spectral set, where $N=p^{n}q$ for
distinct primes $p$, $q$. Then $A$ tiles $\mathbb{Z}_{N}$ by translations.

The proof relies heavily on the structure of vanishing sums of roots of unity
of order $N$, where $N$ is divided by at most two primes. The fact that such
sums are _nonnegative_ linear combinations of $p$\- and $q$-cycles [13] gives
efficient bounds on spectral sets $A\subseteq\mathbb{Z}_{N}$; we believe that
these techniques can be extended to every cyclic group of order
$N=p^{n}q^{m}$. Unfortunately, we have not managed to conclude our proof in
this more general setting so far.

The proof for the direction Tile$\Rightarrow$Spectral for $N=p^{n}q^{m}$ seems
easier for the following reasons: firstly, it is facilitated greatly by the
conditions (T1) and (T2) to be defined in Section 2, that are equivalent to
the tiling condition [1] when $N$ has at most two distinct prime factors and
thus impose a very strong structure on the tiling set $A$, yielding also the
spectrality thereof (not to mention that the tiling condition alone yields
$\\#A\mid N$, restricting the size of $A$ down to the divisors of $N$);
secondly, although spectrality gives some structure on $A$, connecting the set
of roots of the mask polynomial of $A$ with the difference set of a spectrum
$B$, it is difficult to obtain new information about $A$ when $A(X)$ vanishes
on primitive root of unity whose order is not prime power. We overcome this
obstacle by virtue of the structure of the vanishing sums of roots of unity,
obtaining strong lower bounds for the size of the spectrum, leading to a
contradiction whenever one of the conditions (T1) and (T2) is not met. The
latter method is most efficient when

$\widehat{\chi}_{A}(d_{1})=\widehat{\chi}_{A}(d_{2})=\dotsb=\widehat{\chi}_{A}(d_{k})=0,$

for $d_{1}\mid d_{2}\mid\dotsb\mid d_{k}\mid N$ and $k$ relatively large. This
successive divisibility property occurs much more frequently when the power of
$p$ or $q$ in the factorization of $N$ is small, and this limitation did not
allow us to go further than the case $N=p^{n}q$, for the time being.

Lastly, we have to emphasize that these techniques can be no further extended.
Consider the following vanishing sum of roots of unity of order $N=pqr$, where
$p$, $q$, $r$ are distinct primes:

$(\omega_{p}+\omega_{p}^{2}+\dotsb+\omega_{p}^{p-1})(\omega_{q}+\omega_{q}^{2}+\dotsb+\omega_{q}^{q-1})+(\omega_{r}+\omega_{r}^{2}+\dotsb+\omega_{r}^{r-1})=(-1)(-1)+(-1)=0,$

where $\omega_{n}=e^{2\pi i/n}$. It is known that this sum cannot be expressed
as a nonnegative linear combination of $p$-, $q$-, or $r$-cycles [13],
therefore we cannot obtain strong lower bounds on the size of a spectral
subset of $\mathbb{Z}_{N}$. Whether sums such as the above lead to a
counterexample in some finite cyclic group and eventually in $\mathbb{R}$,
remains to be seen.

##  2 Preliminaries

Let $\mathbb{Z}_{N}$ denote the ring of integers modulo $N$. With every
(multi)set $A$ with elements from $\mathbb{Z}_{N}$, we associate a polynomial
in the quotient ring $R=\mathbb{Z}[X]/(X^{N}-1)$, say

$A(X)=\sum_{a\in A}m_{a}X^{a},$

where $m_{a}$ is the multiplicity of $a$ in the multi-set $A$. $A$ is a proper
set, if and only if $A(X)$ has coefficients $0$ and $1$ (it is understood that
we write any element in $R$ as a linear combination of $1,X,\dotsc,X^{N-1}$).
$A(X)$ is called the _mask polynomial_ of $A$; it has the following connection
with the Fourier transform of the characteristic function of $A$:

$\widehat{\chi}_{A}(n)=A(\omega^{-n}),$

where $\omega=e^{2\pi i/N}$ throughout this paper.

A subset $A\subseteq\mathbb{Z}_{N}$ is called _spectral_ if there is a set $B$
with $\\#A=\\#B$, such that the set of exponentials

$x\mapsto e^{\frac{2\pi ibx}{N}},\ b\in B$

is orthogonal on $A$ with the usual inner product.

######  Theorem 2.1.

Let $A\subseteq\mathbb{Z}_{N}$ be spectral. Then

$B-B\subseteq{\left\\{{0}\right\\}}\cup{\left\\{{n:A(\omega^{n})=0}\right\\}}.$

###### Proof.

By definition, we get

$\sum_{a\in A}e^{\frac{2\pi i(b-b^{\prime})a}{N}}=0,$

whenever $b,b^{\prime}$ are distinct elements of $B$. This is equivalent to
the condition

$B-B\subseteq{\left\\{{0}\right\\}}\cup{\left\\{{n:A(\omega^{n})=0}\right\\}}.\qed$

There is a natural action of the Galois group

$G=\textnormal{Gal}(\mathbb{Q}(\omega)/\mathbb{Q})\cong\mathbb{Z}_{N}^{\star}$

on the values of $A(X)$, given by

$\sigma(A(\omega^{a}))=A(\omega^{ag}),$ (2.1)

where $\sigma\in G$ is determined by $\sigma(\omega)=\omega^{g}$, for some
$g\in\mathbb{Z}_{N}^{\star}$. Therefore, in order to determine the support of
$A(\omega^{n})$ we only need to evaluate at the divisors of $N$. For any
integer $N$, the _divisor class_ of $N$ with respect to $n>0$, a divisor of
$N$, is simply $n\mathbb{Z}_{N}^{\star}$, which is the set of residues $\bmod
N$ whose greatest common divisor with $N$ is equal to $n$. The set
${\left\\{{n:A(\omega^{n})=0}\right\\}}$ is just a union of divisor classes.
This is equivalent to the fact that if $A(\omega^{d})=0$ for some $d\mid N$,
then $\Phi_{N/d}(X)\mid A(X)$, where $\Phi_{n}(X)$ denotes the $n$-th
cyclotomic polynomial.

We also denote

$D:={\left\\{{n\in\mathbb{N}:n\mid N,A(\omega^{m})\neq 0,\text{ for all
}m\text{ with }n\mid m\text{ and }m\mid N}\right\\}}.$

######  Proposition 2.2.

If $A\subseteq\mathbb{Z}_{N}$ is spectral and $n\in D$, then $\\#A\leq n$.

###### Proof.

By Theorem 2.1 and the hypothesis we get

$(B-B)\cap n\mathbb{Z}_{N}={\left\\{{0}\right\\}}$

for a spectrum $B$ of $A$. Hence, no two distinct elements of $B$ can have the
same residue $\bmod n$, thus obtaining $\\#A=\\#B\leq n$. ∎

The following properties on polynomials were introduced by Coven and
Meyerowitz [1] in their effort to characterize finite sets that tile the
integers by translations. We will adapt this definition for subsets of
$\mathbb{Z}_{N}$.

######  Definition 2.3.

Let $A\subseteq\mathbb{Z}_{N}$, and let $A(X)$ be its mask polynomial. The set
of prime powers dividing $N$ is denoted by $S$; define

$S_{A}={\left\\{{s\in S:\Phi_{s}(X)\mid A(X)}\right\\}}.$

We say that $A$ satisfies the property (T1), if

$\\#A=A(1)=\prod_{s\in S_{A}}\Phi_{s}(1),$

and that it satisfies (T2), if for every distinct elements
$s_{1},\dotsc,s_{m}$ of $S_{A}$, $\Phi_{s_{1}\dotsm s_{m}}(x)$ divides $A(x)$.

######  Theorem 2.4 (Theorem A [1]).

If $A\subseteq\mathbb{Z}_{N}$ satisfies (T1) and (T2), then $A$ tiles
$\mathbb{Z}_{N}$ by translations.

Let $d\mid N$. A $d$-cycle is a coset by the cyclic subgroup of $N$ with $d$
elements, that is, a set of the form

${\left\\{{j,j+N/d,j+2N/d,\dotsc,j+(d-1)N/d}\right\\}}.$

Especially in the case when $N$ has only two prime divisors, the following
Lemma allows us to discern the structure of $n\cdot A:={\left\\{{na:a\in
A}\right\\}}$ (a multi-set), whenever $A(\omega^{n})=0$. In particular, it
says that $n\cdot A$ must be the union of $p$\- and $q$-cycles.

######  Lemma 2.5.

Let $n\mid N$ be such that $N/n$ has at most two prime divisors, say $p$ and
$q$. If $A(\omega^{n})=0$, then

$A(X^{n})\equiv
P_{n}(X^{n})\Phi_{p}(X^{N/p})+Q_{n}(X^{n})\Phi_{q}(X^{N/q})\bmod X^{N}-1,$
(2.2)

where $P_{n}$ and $Q_{n}$ have nonnegative coefficients.

###### Proof.

By definition, $A(\omega^{n})$ is a vanishing sum of roots of unity of order
$N/n$, in particular

$0=A(\omega^{n})=\sum_{a\in A}\omega^{na}.$

As $N/n$ has at most two prime divisors, $p$ and $q$, the above sum can be
written as linear combination with _nonnegative integer coefficients_ of sums
of the form

$\omega^{k}(1+\omega_{p}+\omega_{p}^{2}+\dotsb+\omega_{p}^{p-1})$

or

$\omega^{k}(1+\omega_{q}+\omega_{q}^{2}+\dotsb+\omega_{q}^{q-1}),$

due to Theorem 3.3 from [13], which shows that $n\cdot A$ is the union of
$p$\- and $q$-cycles, as multi-sets. Every $p$-cycle has a mask polynomial of
the form $X^{k}\Phi_{p}(X^{N/p})$; union of multi-sets corresponds to addition
of the mask polynomials, thus obtaining (2.2). We note that the argument of
$P_{n},Q_{n}$ is $X^{n}$, simply because $A(X^{n})$ can be expressed in terms
of powers of $X^{n}$, as $n\mid N$. ∎

######  Remark.

If $N/n$ has only one prime divisor, say $p$, then it is understood that
$Q_{n}\equiv 0$.

It is also useful to find conditions under which $n\cdot A$ cannot be written
as a union of $p$-cycles or $q$-cycles only, or equivalently,
$P_{n}Q_{n}\not\equiv 0$, for every such possible choice of $P_{n}$ and
$Q_{n}$.

######  Proposition 2.6.

Let $N$ have only two prime divisors, say $p$ and $q$, and $A(\omega^{n})=0$,
for some $n\mid N$, so that

$A(X^{n})\equiv
P_{n}(X^{n})\Phi_{p}(X^{N/p})+Q_{n}(X^{n})\Phi_{q}(X^{N/q})\bmod X^{N}-1.$

If $A(\omega^{np^{a}})\neq 0$ for some $a>0$, then $P_{n}\not\equiv 0$, and if
$A(\omega^{nq^{b}})\neq 0$ for some $b>0$, then $Q_{n}\not\equiv 0$.

###### Proof.

If $Q_{n}\equiv 0$, then

$A(X^{nq^{b}})\equiv P_{n}(X^{nq^{b}})\Phi_{p}(X^{Nq^{b}/p})\equiv
P_{n}(X^{nq^{b}})\Phi_{p}(X^{N/p}),$

thus obtaining $A(\omega^{nq^{b}})=0$. The other case is proven similarly. ∎

(T1) and (T2) are conjectured to hold if and only if $A$ tiles
$\mathbb{Z}_{N}$ [1] (this conjecture was initially formulated in
$\mathbb{Z}$). For every set $A\subseteq\mathbb{Z}_{N}$, however, a weaker
property than (T1) holds, that is very useful for bounding $\\#A$.

######  Proposition 2.7.

Let $A\subseteq\mathbb{Z}_{N}$. Then $\prod_{s\in S_{A}}\Phi_{s}(1)$ divides
$\\#A$. In particular, if $p$ is a prime divisor of $N$, and there are $m$
integers $0<a_{1}<a_{2}<\dotsb<a_{m}$ such that $A(\omega^{N/p^{a_{i}}})=0$
for all $1\leq i\leq m$, then $p^{m}\mid\\#A$, where $p^{a_{m}}\mid N$.

###### Proof.

By definition, $A(X)$ is divided by $\prod_{s\in S_{A}}\Phi_{s}(X)$ in
$\mathbb{Z}[X]$. Putting $X=1$, we get the desired result. ∎

In the rest of this paper we will prove Fuglede’s conjecture on
$\mathbb{Z}_{N}$, where $N=p^{n}q$, for $p\neq q$ primes. The direction
Tile$\Rightarrow$Spectral can be deduced from the work of Coven-Meyerowitz [1]
and Łaba [12] in the more general case when $N$ has at most two prime
divisors, and is shown in Section 3. The direction Spectral$\Rightarrow$Tile
is proven in section 5. For completeness, we will also present the
Spectral$\Rightarrow$Tile proof for $N=p^{n}$ in section 4 due to its elegance
and brevity, although a proof for this case follows from [12]; a different
proof also appeared in [4].

One final tool that will be very useful in this note, is the following:

######  Lemma 2.8.

Let $m,n>0$ be two relatively prime integers, and $0<k<mn$ another integer.
Then, there is at most one pair $(s,t)$ of nonnegative integers, such that
$k=sm+tn$. If $k=mn$, then there are exactly two such pairs, namely $(n,0)$
and $(0,m)$.

###### Proof.

Assume that $0<k<mn$ and there is a pair $(s,t)$ such that $k=sm+tn$, with
$s,t\geq 0$. All other pairs of integer solutions $(s^{\prime},t^{\prime})$ to
the Diophantine equation $k=s^{\prime}m+t^{\prime}n$, satisfy
$s^{\prime}=s-nx$, $t^{\prime}=t+mx$, for $x\in\mathbb{Z}$. If $x>0$, then
$s^{\prime}\leq s-n\leq k/m-n<0$, whereas if $x<0$ we get $t^{\prime}\leq
t-m\leq k/n-m<0$. The case $k=mn$ is proven similarly. ∎

##  3 Tile$\Rightarrow$Spectral, $N=p^{n}q^{m}$

In this section, we will review the proof of Theorem 1.1. The proof is not
new, and is based on combined arguments from [1] and [12], even though the
case for finite cyclic groups is not explicitly mentioned. We will need the
following tools from [1] and [12].

######  Lemma 3.1 (Lemma 1.3, [1]).

Let $N$ be an integer and let $A$ and $B$ be finite multisets of nonnegative
integers with mask polynomials $A(X)$ and $B(X)$. Then the following
statements are equivalent. Each forces $A$ and $B$ to be sets such that
$\\#A\\#B=A(1)B(1)=N$.

  1. 1.

$A\oplus(B\oplus N\mathbb{Z})=\mathbb{Z}$ is a tiling.

  2. 2.

$A\oplus B$ is a complete set of residues modulo $N$.

  3. 3.

$A(X)B(X)\equiv 1+X+\dotsb X^{N-1}\bmod X^{N}-1$.

  4. 4.

$N=A(1)B(1)$ and for every factor $t>1$ of $N$, the cyclotomic polynomial
$\Phi_{t}(X)$ is a divisor of $A(X)$ or $B(X)$.

######  Lemma 3.2 (Lemma 2.1, [1]).

Let $A(X)$ and $B(X)$ be polynomials with coefficients 0 and 1, $N=A(1)B(1)$,
and $S$ the set of prime power factors of $N$. If $\Phi_{t}(X)$ divides $A(X)$
or $B(X)$ for every factor $t>1$ of $N$, then

  1. 1.

$A(1)=\prod_{s\in S_{A}}\Phi_{s}(1)$ and $B(1)=\prod_{s\in S_{B}}\Phi_{s}(1)$.

  2. 2.

$S_{A}$ and $S_{B}$ are disjoint sets whose union is $S$.

Now, assume that $A$ tiles $\mathbb{Z}_{N}$ by translations, and let $B$ be
the complementary tile.We may assume that
$A\subseteq{\left\\{{0,1,2,\dotsc,N-1}\right\\}}$, and also assume that $A$
tiles $\mathbb{Z}$ by translations; furthermore, this tiling has period $N$,
i. e. $A\oplus(B\oplus N\mathbb{Z})=\mathbb{Z}$. We warn the reader, that not
only do we have to prove that $A$, as a subset of $\mathbb{Z}$, is spectral,
but also that the spectrum is a subset of $N^{-1}\mathbb{Z}$, in order to
claim that $A$, _as a subset of $\mathbb{Z}_{N}$_, is spectral.

We see that conditions $(1)$ and $(2)$ of Lemma 3.1 are satisfied, hence $(4)$
is satisfied as well, which is just the hypothesis of Lemma 3.2. By $(2)$ of
Lemma 3.2, we get that $S_{A}\subseteq S$. Next, we will use the following two
theorems:

######  Theorem 3.3 (Theorems B1 and B2, [1]).

Let $A$ be a finite set of nonnegative integers with corresponding polynomial
$A(x)=\sum_{a\in A}x^{a}$. If $A$ tiles the integers, and $\\#{A}$ is divided
by at most two primes, then $A$ satisfies (T1) and (T2).

######  Theorem 3.4 (Theorem 1.5(i), [12]).

If $A\subseteq\mathbb{Z}$ satisfies (T1) and (T2), then $A$ has a spectrum.

######  Remark.

The important part of the proof, is that the spectrum is explicitly
constructed with respect to $S_{A}$, namely, the set of all

$\sum_{s\in S_{A}}\frac{k_{s}}{s},$

where $k_{s}\in{\left\\{{0,1,\dotsc,p-1}\right\\}}$, $s\in S_{A}$ and
$s=p^{\alpha}$, is proven to be a spectrum of $A$, when it satisfies (T1) and
(T2) (see the beginning of Section 2 [12]).

So, if $A\subseteq{\left\\{{0,1,\dotsc,N-1}\right\\}}$ tiles $\mathbb{Z}$ by
translations, the tiling having period $N$, which has at most two prime
divisors, then $A$ satisfies (T1) and (T2) by Theorem 3.3, as $\\#{A}$ divides
$N$. Next, by Theorem 3.4 we get that $A$ is spectral; by Lemma 3.2 and the
Remark above, we get that the spectrum belongs to $N^{-1}\mathbb{Z}$, hence
$A\subseteq\mathbb{Z}_{N}$ is spectral, completing the proof.

##  4 Spectral$\Rightarrow$Tile, $N=p^{n}$

Suppose $\Lambda\subseteq\mathbb{Z}_{N}$ is a spectrum of $A$ and let
$p^{\nu_{1}},\ldots,p^{\nu_{k}}$ be the divisors $d$ of $N$ such that

$\Phi_{d}(X)\mid A(X)=\sum_{a\in A}X^{a}.$

We have

$\Phi_{p^{\nu_{i}}}(X)=1+X^{p^{\nu_{i}-1}}+X^{2p^{\nu_{i}-1}}+\cdots+X^{(p-1)p^{\nu_{i}-1}}.$

Write $E_{\nu}=\\{0,p^{\nu-1},2p^{\nu-1},\ldots,(p-1)p^{\nu-1}\\}$ so that
$E_{\nu_{i}}(X)=\Phi_{p^{\nu_{i}}}(X)$. Define next the set

$E=E_{\nu_{1}}+E_{\nu_{2}}+\cdots E_{\nu_{k}},$ (4.1)

and notice the sum is direct as $e_{1}+\cdots+e_{k}\in E$ is determines
$e_{i}\in E_{\nu_{i}}$ from the ${\nu_{i}}$-th digit in its expansion to base
$p$. This observation implies that $|E|=p^{k}$.

Notice also that $A(X)$ and $E(X)$ have the same zeros at the $N$-th roots of
unity. This implies that $\Lambda$ is also orthogonal on $E$ as this is
determined by the zeros of $E(X)$ at the $N$-th roots of unity. From the
orthogonality we obtain

$|A|=|\Lambda|\leq|E|=p^{k}.$

Let $B\subseteq\mathbb{Z}_{N}$ be the sum of those $E_{\nu}$,
$\nu=1,2,\ldots,n$, not appearing in (4.1). This sum is again direct, as with
the sum (4.1), so we obtain

$|B|=p^{n-k}.$

It follows that $A(X)B(X)$ vanishes on all $N$-th roots of unity except 1,
which implies that $A+B$ is a tiling of $\mathbb{Z}_{N}$ at some level $\ell$.
Then

$\ell p^{n}=|A|\cdot|B|=p^{k}p^{n-k}=p^{n},$

so that $\ell=1$ and $A+B$ is a tiling of $\mathbb{Z}_{N}$ at level 1.

##  5 Spectral$\Rightarrow$Tile, $N=p^{n}q$

In this Section we prove Theorem 1.2. As mentioned in the Introduction, our
method is to establish a contradiction to the assumption of a spectral set not
satisfying (T1) and (T2). When one of these properties fails, the spectrality
of $A$ induces the existence of many roots of the mask polynomial $A(X)$ of
the form

$A(\omega^{d_{1}})=\dotsb=A(\omega^{d_{k}})=0,$ (5.1)

where $d_{1}\mid\dotsb\mid d_{k}\mid N$, which further gives

$\Phi_{N/d_{1}}(1)\dotsm\Phi_{N/d_{k}}(1)\mid\\#A$ (5.2)

and this produces a strong lower bound for $\\#A$ via the vanishing sums of
roots of unity which eventually leads to a contradiction. Obtaining such lower
bounds were a lot easier for $N=p^{n}$, as (5.2) would yield $p^{k}\mid\\#A$.
This follows from the fact that $\Phi_{p^{s}}(1)=p$, for a prime $p$, whereas
$\Phi_{M}(1)=1$, when $M$ is not a prime power, and this explains the added
complexity of this proof, especially when the fractions $N/d_{i}$ are not
prime powers. Lemma 2.5 is the main tool here; although it does not produce
some divisibility condition such as (5.2), it gives some good lower bounds in
terms of $k$ when (5.1) holds. See for example, the proofs of the Claims 1 and
2 below.

We distinguish two cases, depending on whether $q$ belongs to $D$.

$q\in D$ We have $\\#A\leq q$ from Proposition 2.2. Furthermore, the property
(T2) holds vacuously, so we only need to prove (T1) due to Theorem 2.4. If
$\\#A=1$, $A(X)$ is a monomial and has no root of the form $\omega^{d}$, in
particular, $A(\omega^{p^{n}})\neq 0$ and (T1) holds.

If $\\#A>1$, and $A$ is spectral, then $B-B\neq{\left\\{{0}\right\\}}$ for a
spectrum $B$, so $A(\omega^{d})$ must vanish somewhere. Since $q\in D$, there
must be some nonnegative $a\leq n$ such that $A(\omega^{p^{a}})=0$, so that

$A(X^{p^{a}})\equiv
P(X^{p^{a}})\Phi_{p}(X^{N/p})+Q(X^{p^{a}})\Phi_{q}(X^{N/q})\mod X^{N}-1,$

by Lemma 2.5. $q\in D$ yields $A(\omega^{p^{a}q})\neq 0$, hence by Proposition
2.6 we get $Q\not\equiv 0$, so that $\\#A=A(1)=pP(1)+qQ(1)\geq q$, leading to
$P\equiv 0$ and $Q(1)=1$. This certainly implies that $A(\omega^{p^{n}})=0$,
$A(1)=q$, hence (T1) holds in any case where $A$ is spectral and $q\in D$. By
Theorem 2.4, $A$ tiles $\mathbb{Z}_{N}$ by translations.

$q\notin D$ The size of a spectral set $A\subseteq\mathbb{Z}_{N}$ depends on
the number of roots of $A(X)$ of the form $\omega^{p^{a}q}$.

######  Definition 5.1.

The numbers $0\leq a_{1}<\dotsb<a_{m}\leq n-1$ are such that
${\left\\{{a:A(\omega^{p^{a}q})=0}\right\\}}={\left\\{{a_{1},\dotsc,a_{m}}\right\\}}$.
If no such root exists, we simply put $m=0$. Also, if $A$ is spectral and $B$
a spectrum of $A$, we denote $B_{i}={\left\\{{b\in B:b\equiv i\bmod
q}\right\\}}$. The $p$-adic expansion of the least nonnegative residue $\bmod
p^{n}$ of $b\in B$ has the form

$b\equiv b_{0}+b_{1}p+\dotsb+b_{n-1}p^{n-1}\bmod p^{n},$

where $0\leq b_{i}\leq p-1$ for all $i$. We say that the $p$-adic expansions
of $b$ and $b^{\prime}$ coincide at $a_{1},\dotsc,a_{m}$, if
$b_{a_{i}}=b^{\prime}_{a_{i}}$ for $1\leq i\leq m$, that is, they have the
same $p$-adic digits at those places. Finally, $p^{i}\mid\\!\mid b$ exactly
when $b_{i}$ is the smallest nonzero $p$-adic digit of $b$.

We have seen in the previous subsection that Spectral$\Rightarrow$Tile holds
when $m=0$. By induction, we may assume that Spectral$\Rightarrow$Tile holds
for all nonnegative integers up to $m-1$, for some $m>0$, where
${\left\\{{a:A(\omega^{p^{a}q})=0}\right\\}}={\left\\{{a_{1},\dotsc,a_{m}}\right\\}}$,
exactly as in the Definition above.

######  Proposition 5.2.

Suppose that $A\subseteq\mathbb{Z}_{N}$ is spectral and $\\#A>p^{m}$. Then

  1. 1.

$\\#A=p^{m}q$.

  2. 2.

$\\#B_{i}=p^{m}$ for all $i$, and the elements of $B_{i}$ have all the
possible $p$-adic expansions at $a_{1},\dotsc,a_{m}$, each appearing exactly
once.

  3. 3.

For every $i,j$ and every $b\in B$, there is $b^{\prime}$ such that
$p^{a_{j}}\mid\\!\mid b-b^{\prime}$.

  4. 4.

There is some $a\neq a_{i}$ for all $i$, such that $b-b^{\prime}\in
p^{a}\mathbb{Z}_{N}^{\star}$ for some $b,b^{\prime}\in B$, and
$A(\omega^{p^{a}})=0$.

  5. 5.

If $a>a_{m}$, then $A(\omega^{p^{a_{j}}})=0$ for all $j$, as well as
$A(\omega^{p^{a^{\prime}}})=0$ for all $a^{\prime}\geq a$; in particular,
$A(\omega^{p^{n}})=0$, or equivalently, $p^{n}\notin D$.

  6. 6.

If $a<a_{m}$, then $A(X^{p^{a_{m}}q})\equiv
p^{m-1}qX^{k}\Phi_{p}(X^{N/p})\bmod X^{N}-1$, for some $k$.

###### Proof.

As $\\#B>p^{m}$, there are at least two distinct elements of $B$, say
$b,b^{\prime}$, that have the same digits at places
$a_{1},a_{2},\dotsc,a_{m}$. If $q\mid b-b^{\prime}$, then
$b-b^{\prime}\notin\bigcup_{1\leq i\leq m}p^{a_{i}}q\mathbb{Z}_{N}^{\star}$,
contradicting the fact that $A(X)$ has exactly $m$ roots of the form
$\omega^{p^{a}q}$, due to Theorem 2.1. Hence $q\nmid b-b^{\prime}$, and
$b-b^{\prime}\in p^{a}\mathbb{Z}_{N}^{\star}$, so that $A(\omega^{p^{a}})=0$
by Theorem 2.1, where $a$ is the first place where the $p$-adic expansions of
$b$ and $b^{\prime}$ differ, therefore $a\neq a_{i}$ for all $i$, proving (4).

With the same argument we can show that $\\#B_{i}\leq p^{m}$ for all $i$; in
this case any two elements $b,b^{\prime}\in B_{i}$ satisfy $q\mid
b-b^{\prime}$, so by Theorem 2.1 and the hypothesis, they must have at least
one different digit at places $a_{1},\dotsc,a_{m}$, yielding $\\#B_{i}\leq
p^{m}$, and $\\#B=\\#A\leq p^{m}q$. For convenience, put $n_{i}=p^{a_{i}}q$
and $d=p^{a}$. By Lemma 2.5 we get

$A(X^{d})\equiv
P_{d}(X^{d})\Phi_{p}(X^{N/p})+Q_{d}(X^{d})\Phi_{q}(X^{N/q})\bmod X^{N}-1,$
(5.3)

where $Q_{d}\not\equiv 0$ due to Proposition 2.6, as $A(\omega^{p^{a}q})\neq
0$.

Now, consider the largest index $i$, such that $a>a_{i}$, assuming first that
there is such an index, i. e. $a>a_{1}$; otherwise, we put $i=0$. Put
$u=dq=p^{a}q$, and denote by ${\left\|{A(X)}\right\|}_{\infty}$ the largest
coefficient of $A(X)$ in $R=\mathbb{Z}[X]/(X^{N}-1)$ written as a linear
combination of $1,X,\dotsc,X^{N-1}$.

######  Claim 1.

${\left\|{A(X^{u})}\right\|}_{\infty}\geq p^{i}q$.

###### Proof of Claim.

From (5.3) we obtain

$A(X^{u})\equiv P_{d}(X^{u})\Phi_{p}(X^{N/p})+qQ_{d}(X^{u})\bmod X^{N}-1.$
(5.4)

If $i=0$, then from $Q_{d}\not\equiv 0$ we get that some coefficient of
$A(X^{u})$ is at least as large as $q$, as desired. Suppose that $i>0$. By
repeated application of Lemma 2.5, we have

$A(X^{n_{j}})\equiv P_{n_{j}}(X^{n_{j}})\Phi_{p}(X^{N/p}),$

for all $j$. For $j=1$, if we replace $X$ by $X^{n_{2}/n_{1}}$, we also get
$A(X^{n_{2}})\equiv pP_{n_{1}}(X^{n_{2}})$, and comparing this with
$A(X^{n_{2}})\equiv P_{n_{2}}(X^{n_{2}})\Phi_{p}(X^{N/p})$, we deduce that
$\Phi_{p}(X^{N/p})$ divides $P_{n_{1}}(X^{n_{2}})$, thus obtaining

$A(X^{n_{2}})\equiv pP_{n_{1}}^{1}(X^{n_{2}})\Phi_{p}(X^{N/p}),$

for some polynomial $P_{n_{1}}^{1}$ with positive integer coefficients.
Proceeding inductively, we can get

$A(X^{n_{i}})\equiv p^{i-1}P_{n_{1}}^{i-1}(X^{n_{i}})\Phi_{p}(X^{N/p}),$

hence

$A(X^{u})\equiv p^{i}P_{n_{1}}^{i-1}(X^{u}),$ (5.5)

since $a_{i}<a<a_{i+1}$, where $P_{n_{1}}^{i-1}$ has positive integer
coefficients. Comparing (5.4) and (5.5), we get that $p^{i}$ divides all
coefficients of $A(X^{u})$. Furthermore, since $A(\omega^{u})\neq 0$, we
deduce that there is at least a $p$-cycle on which the elements of $n_{m}\cdot
A$ do not have the same multiplicity; using (5.4) again, we deduce that two
such multiplicities must differ by a multiple of $q$. Therefore, by (5.5) we
conclude that their difference is a multiple of $p^{i}q$, and since they are
both nonnegative and distinct, we finally get
${\left\|{A(X^{u})}\right\|}_{\infty}\geq p^{i}q$. ∎

If $a>a_{m}$ the claim yields ${\left\|{A(X^{u})}\right\|}_{\infty}\geq
p^{m}q$, while on the other hand $A(1)\leq p^{m}q$, therefore, $A(X^{u})\equiv
p^{m}qX^{uk}$, for some $k$. Then, (5.4) yields $P_{d}\equiv 0$, so for all
$a^{\prime}\geq a$ (5.4) gives

$A(X^{p^{a^{\prime}}})\equiv Q_{d}(X^{p^{a^{\prime}}})\Phi_{q}(X^{N/q})\bmod
X^{N}-1$

establishing (5) (and (1) when $a>a_{m}$).

For the rest of the proof, we suppose $a<a_{m}$, and we will estimate
${\left\|{A(X^{n_{m}})}\right\|}_{\infty}$.

######  Claim 2.

${\left\|{A(X^{n_{j+1}})}\right\|}_{\infty}\geq
p{\left\|{A(X^{n_{j}})}\right\|}_{\infty}$.

###### Proof of Claim.

Since $A(\omega^{n_{j}})=0$, we obtain from Lemma 2.5 $A(X^{n_{j}})\equiv
P_{n_{j}}(X^{n_{j}})\Phi_{p}(X^{N/p})$. The largest coefficient of
$A(X^{n_{j}})$ would appear in a $p$-cycle, hence
${\left\|{A(X^{n_{j}+1})}\right\|}_{\infty}\geq
p{\left\|{A(X^{n_{j}})}\right\|}_{\infty}$, and finally,
${\left\|{A(X^{n_{j+1}})}\right\|}_{\infty}\geq{\left\|{A(X^{n_{j}+1})}\right\|}_{\infty}\geq
p{\left\|{A(X^{n_{j}})}\right\|}_{\infty}$, as desired. ∎

Applying Claims 1 and 2, we obtain
${\left\|{A(X^{n_{m}})}\right\|}_{\infty}\geq p^{m-1}q$. The largest
coefficient of $A(X^{n_{m}})$ would appear on a $p$-cycle, since
$A(\omega^{n_{m}})=0$, so we would get $\\#A=A(1)\geq p^{m}q$; but we have
already shown that $\\#A\leq p^{m}q$, so $A(1)=p^{m}q$ and $A(X^{n_{m}})\equiv
p^{m-1}qX^{k}\Phi_{p}(X^{N/p})$ thus proving (6) and (1) in all cases. (2)
follows immediately from (1), as $\\#B=p^{m}q$, hence $\\#B_{i}=p^{m}$ for all
$i$. Finally, (3) is a direct consequence of (2); let $i$ and $b\in B_{i}$ be
arbitrary. For every $j$, there is some $b^{\prime}\in B_{i}$ whose $p$-adic
expansion is the same on $a_{1},\dotsc,a_{j-1}$ but differs on $a_{j}$, due to
(2). Then, $b-b^{\prime}\in p^{c}q\mathbb{Z}_{N}^{\star}$, for some $c\leq
a_{j}$. By Theorem 2.1 and the hypothesis, $c=a_{l}$, for some $l$, but
$c=a_{l}$ cannot hold for $l<j$, thus $l=j$, completing the proof. ∎

Next, we will make use of the induction assumption.

######  Proposition 5.3.

Let $A\subseteq\mathbb{Z}_{N}$ spectral and
${\left\\{{a:A(\omega^{p^{a}q})=0}\right\\}}={\left\\{{a_{1},\dotsc,a_{m}}\right\\}}$.
Define a partition of $A$ into sets $A_{0},A_{1},\dotsc,A_{p-1}$ such that

$a\in A_{j}\Longleftrightarrow p^{a_{m}}qa\in[\frac{N}{p}j,\frac{N}{p}(j+1)).$

Then, all $A_{j}$ have the same cardinality, and if $d\mid p^{a_{m}-1}q$
(assuming $a_{m}>1$) with $A(\omega^{d})=0$, then $A_{j}(\omega^{d})=0$ for
all $j$. Furthermore, $A_{j}(\omega^{p^{a_{m}}q})\neq 0$.

###### Proof.

By Proposition 2.6 we have $A(X^{p^{a_{m}}q})\equiv
P(X^{p^{a_{m}}q})\Phi_{p}(X^{N/p})\bmod X^{N}-1$, and without loss of
generality we have $\deg P(X^{p^{a_{m}}q})<N/p$. The hypothesis implies that

$A_{j}(X^{p^{a_{m}}q})\equiv P(X^{p^{a_{m}}q})X^{\frac{N}{p}j}\bmod X^{N}-1,$

so that $A_{j}(1)=P(1)$ for all $j$, proving the first part.

Now, let $d=p^{k}q^{c}$ with $A(\omega^{d})=0$, where $k<a_{m}$ and $c=0$ or
$1$. We will use the language of multi-sets for this part; $d\cdot A$ is a
union of $p$\- and $q$-cycles, as multi-sets. We will show that every such
cycle must belong exclusively to one of the mutli-sets $d\cdot A_{j}$, hence
each $d\cdot A_{j}$ is also a union of $p$\- and $q$-cycles, leading to
$A_{j}(\omega^{d})=0$, for all $j$.

Indeed, suppose that $da$ is a part of such a cycle; the other member of said
cycle are $da+lN/r$, where $r=p$ or $q$, and $0\leq l\leq r-1$. Define
$d^{\prime}$ by $dd^{\prime}=p^{a_{m}}q$. Then, for every $l$,

$d^{\prime}da\equiv d^{\prime}(da+lN/r)\bmod N,$

since $d^{\prime}l\equiv 0\bmod r$. This is true in case where $r=p$, because
$p\mid d^{\prime}$. If $r=q$, we can have $q$-cycles only if $c=0$, or
equivalently $q\mid d^{\prime}$. The above congruence clearly shows that any
such cycle belongs to one of the multi-sets $d\cdot A_{j}$ (we remark that
these multi-sets are mutually exclusive).

For the last part, we just note that $P(\omega^{p^{a_{m}}q})\neq 0$, otherwise
$\Phi_{p}(X^{N/p})$ would be a factor of $P(X^{p^{a_{m}}q})$, an
impossibility, since the degree of the latter does not exceed $N/p$. ∎

######  Proposition 5.4.

Let $A\subseteq\mathbb{Z}_{N}$ be spectral with
${\left\\{{a:A(\omega^{p^{a}m})=0}\right\\}}={\left\\{{a_{1},\dotsc,a_{m}}\right\\}}$
and $B$ a spectrum. Consider the partition of $A$ into
$A_{0},A_{1},\dotsc,A_{p-1}$ as in Proposition 5.3. Suppose that the maximal
$a\notin{\left\\{{a_{1},\dotsc,a_{m}}\right\\}}$ with $0\leq a\leq n$ and
$(B-B)\cap p^{a}\mathbb{Z}_{N}^{\star}\neq\varnothing$ (if it exists!)
satisfies $a<a_{m}$. Then, each $A_{j}$ is spectral (if such $a$ does not
exist, we have the same conclusion).

###### Proof.

Let $B^{i}$ be the subset of elements $B$ whose $p$-adic digit at $a_{m}$ is
equal to $i$. By hypothesis and Proposition 5.3,

$(B^{i}-B^{i})\subseteq(B-B)\setminus(p^{a_{m}}\mathbb{Z}_{N}^{\star}\cup
p^{a_{m}}q\mathbb{Z}_{N}^{\star})\subseteq{\left\\{{d:A(\omega^{d})=0}\right\\}}\setminus
p^{a_{m}}\mathbb{Z}_{N}\subseteq{d:A_{j}(\omega^{d})=0}.$

By pigeonhole principle, we may select one $B^{i}$ such that
$\\#B^{i}\geq\frac{1}{p}\\#B=A_{j}(1)$. This can only be possible if we have
equality, thus showing that each $A_{j}$ is spectral by Theorem 2.1, having
the same spectrum (actually, any $B^{i}$ would serve as such). ∎

If $A(1)=p^{m}q$ and $A(X)$ has at least $m$ roots of the form $\omega^{d}$,
where $d$ is a power of $p$, then $A$ has a special structure.

######  Proposition 5.5.

Let $A\subseteq\mathbb{Z}_{N}$ with $A(1)=p^{m}q$. Suppose that
$A(\omega^{d_{i}})=0$, where $d_{i}$ are increasing powers of $p$, $1\leq
i\leq m$ and $d_{i}\leq p^{n}$. Then either

$A(X^{d_{j}})\equiv P_{d_{j}}(X^{d_{j}})\Phi_{p}(X^{N/p})\bmod X^{N}-1$

for all $j$, where $P_{d_{j}}\not\equiv 0$ have nonnegative coefficients, or

$A(X^{d_{m}})\equiv Q_{d_{m}}(X^{d_{m}})\Phi_{q}(X^{N/q})\bmod X^{N}-1,$

where $Q_{d_{m}}\not\equiv 0$ has nonnegative coefficients.

###### Proof.

By Lemma 2.5 we obtain

$A(X^{d_{j}})\equiv
P_{d_{j}}(X^{d_{j}})\Phi_{p}(X^{N/p})+Q_{d_{j}}(X^{d_{j}})\Phi_{q}(X^{N/q})\bmod
X^{N}-1,$

while on the other hand we can show inductively

$A(X^{d_{j+1}})\equiv
p^{j}P_{d_{1}}^{j}(X^{d_{j+1}})\Phi_{p}(X^{N/p})+\sum_{i=0}^{j}p^{i}Q_{d_{1}}^{i}(X^{d_{j+1}})\Phi_{q}(X^{N/q})\bmod
X^{N}-1,$

for all $j$, where all the polynomials appearing in these two formulae have
nonnegative coefficients, and satisfy the reccurence relations

$P_{d_{1}}^{j}(X^{d_{j+2}})\equiv
P_{d_{1}}^{j+1}(X^{d_{j+2}})\Phi_{p}(X^{N/p})+Q_{d_{1}}^{j+1}(X^{d_{j+2}})\Phi_{q}(X^{N/q})\bmod
X^{N}-1.$

Without loss of generality we can write

$P_{d_{j}}(X^{d_{j}})\equiv
p^{j-1}P_{d_{1}}^{j-1}(X^{d_{j}}),Q_{d_{j}}(X^{d_{j}})\equiv\sum_{i=0}^{j-1}p^{i}Q_{d_{1}}^{i}(X^{d_{j}}),$

for all $j$. Putting $j=m-1$ and $X=1$ we get

$p^{m}q=A(1)=p^{m}P_{d_{1}}^{m-1}(1)+q\sum_{i=0}^{m-1}p^{i}Q_{d_{1}}^{i}(1),$

so by Lemma 2.8 we have either

$P_{d_{1}}^{m-1}(1)=q\ \text{ and }\ \sum_{i=0}^{m-1}p^{i}Q_{d_{1}}^{i}(1)=0,$
(5.6)

or

$P_{d_{1}}^{m-1}(1)=0\ \text{ and }\
\sum_{i=0}^{m-1}p^{i}Q_{d_{1}}^{i}(1)=p^{m}.$ (5.7)

If (5.6) holds, then we have $Q_{d_{1}}^{i}\equiv 0$ for all $i$, hence
$Q_{d_{j}}\equiv 0$ for all $j$, and $A(X^{d_{j}})\equiv
P_{d_{j}}(X^{d_{j}})\Phi_{p}(X^{N/p})$. Otherwise, if (5.7) holds, then
$P_{d_{1}}^{m-1}\equiv 0$, so in this case we get $A(X^{d_{m}})\equiv
Q_{d_{m}}(X^{d_{m}})\Phi_{q}(X^{N/q})$, as desired. ∎

Now, we can proceed with the conclusion of the Spectral$\Rightarrow$Tile
proof. If $p^{n}\in D$, (T2) holds vacuously, so we need only prove (T1),
namely $A(1)=p^{m}$. Suppose on the contrary that $A(1)>p^{m}$, hence by
Proposition 5.2(1) we have $A(1)=p^{m}q$. If a spectrum $B$ satisfies the
hypothesis of Proposition 5.4, then each $A_{j}$ is spectral, so by induction
they satisfy (T1). Since $A_{j}(1)=p^{m-1}q$, we must have
$A_{j}(\omega^{p^{n}})=0$ for all $j$, yielding $A(\omega^{p^{n}})=0$,
contradicting our assumption that $p^{n}\in D$. If the maximal $a$ such that
$(B-B)\cap p^{a}\mathbb{Z}_{N}^{\star}\neq\varnothing$ satisfies $a>a_{m}$,
then by Proposition 5.2(5) we get that $A(\omega^{p^{a_{i}}})=0$ for all
$1\leq i\leq m$. Taking
${\left\\{{d_{1},\dotsc,d_{m}}\right\\}}={\left\\{{p^{a_{1}},\dotsc,p^{a_{m-1}},p^{a}}\right\\}}$
and applying Proposition 5.5 we get either

$A(X^{d_{m}})\equiv P_{d_{m}}(X^{d_{m}})\Phi_{p}(X^{N/p})\bmod X^{N}-1$

or

$A(X^{d_{m}})\equiv Q_{d_{m}}(X^{d_{m}})\Phi_{q}(X^{N/q})\bmod X^{N}-1,$

where in each case $P_{d_{m}},Q_{d_{m}}\not\equiv 0$ have nonnegative
coefficients. However, this contradicts Proposition 2.6, since
$A(\omega^{p^{n}})A(\omega^{p^{a}q})\neq 0$. We conclude that $A$ must satisfy
(T1) as well, thus tiling $\mathbb{Z}_{N}$ by translations due to Theorem 2.4.

Lastly, we suppose that $p^{n}\notin D$, so that $A(\omega^{p^{n}})=0$. In
this case, $q\mid A(1)$ by Proposition 2.7, hence $A(1)>p^{m}$ and
$A(1)=p^{m}q$ by Proposition 5.2(1). Therefore, (T1) holds and it remains to
prove (T2), namely $A(\omega^{p^{a_{i}}})=0$ for all $i$. If the maximal $a$
for which $(B-B)\cap p^{a}\mathbb{Z}_{N}^{\star}\neq\varnothing$ holds,
satisfies $a>a_{m}$, for a spectrum $B$, then by applying Proposition 5.2(5)
we deduce that (T2) holds. Otherwise, $B$ satisfies the conditions of
Proposition 5.4. Therefore, each $A_{j}$ is spectral and satisfies (T2),
yielding $A_{j}(\omega^{p^{a_{i}}})=0$ for all $j$ and $1\leq i\leq m-1$,
hence $A(\omega^{p^{a_{i}}})=0$ for $1\leq i\leq m-1$. It remains to show that
$A(\omega^{p^{a_{m}}})=0$. We remark that an $a$ as described in Proposition
5.4 _actually exists_ in this case due to Proposition 5.2(4), and $a<a_{m}$.
We apply Proposition 5.5 for
${\left\\{{d_{1},\dotsc,d_{m}}\right\\}}={\left\\{{p^{a_{1}},\dotsc,p^{a_{m-1}},p^{a}}\right\\}}$;
if

$A(X^{d_{j}})\equiv P_{d_{j}}(X^{d_{j}})\Phi_{p}(X^{N/p})\bmod X^{N}-1$

for all $j$, where $P_{d_{j}}\not\equiv 0$ have nonnegative coefficients, then
$A(\omega^{p^{a}q})=0$, a contradiction. Therefore,

$A(X^{d_{m}})\equiv Q_{d_{m}}(X^{d_{m}})\Phi_{q}(X^{N/q})\bmod X^{N}-1,$

where $Q_{d_{m}}\not\equiv 0$ has nonnegative coefficients. Substituting $X$
by $X^{p^{a_{m}}/d_{m}}$, we get

$A(X^{p^{a_{m}}})\equiv Q_{d_{m}}(X^{p^{a_{m}}})\Phi_{q}(X^{N/q}),$

yielding $A(\omega^{p^{a_{m}}})=0$, completing the proof.

## Acknowledgments

The authors are grateful to the anonymous reviewer for their valuable remarks.

## References

  * [1] Ethan M. Coven and Aaron Meyerowitz, _Tiling the integers with translates of one finite set_ , J. Algebra 212 (1999), no. 1, 161–174.
  * [2] Dorin Ervin Dutkay and Chun-Kit Lai, _Some reductions of the spectral set conjecture to integers_ , Math. Proc. Cambridge Philos. Soc. 156 (2014), no. 1, 123–135.
  * [3] Aihua Fan, Shilei Fan, Lingmin Liao, and Ruxi Shi, _Fuglede’s conjecture holds in $\mathbb{Q}_{p}$_, 24pp., https://arxiv.org/abs/1512.08904.
  * [4] Aihua Fan, Shilei Fan, and Ruxi Shi, _Compact open spectral sets in $\mathbb{Q}_{p}$_, J. Funct. Anal. 271 (2016), no. 12, 3628–3661.
  * [5] Bálint Farkas, Máté Matolcsi, and Péter Móra, _On Fuglede’s conjecture and the existence of universal spectra_ , J. Fourier Anal. Appl. 12 (2006), no. 5, 483–494.
  * [6] Bent Fuglede, _Commuting self-adjoint partial differential operators and a group theoretic problem_ , J. Functional Analysis 16 (1974), 101–121.
  * [7] Alex Iosevich, Nets Katz, and Terence Tao, _The Fuglede spectral conjecture holds for convex planar domains_ , Math. Res. Lett. 10 (2003), no. 5-6, 559–569.
  * [8] Alex Iosevich, Azita Mayeli, and Jonathan Pakianathan, _The Fuglede conjecture holds in $\mathbb{Z}_{p}\times\mathbb{Z}_{p}$_, Anal. PDE 10 (2017), no. 4, 757–764.
  * [9] Mihail N. Kolountzakis and Máté Matolcsi, _Complex Hadamard matrices and the spectral set conjecture_ , Collect. Math. (2006), Vol. Extra, 281–291.
  * [10] Mihail N. Kolountzakis and Máté Matolcsi, _Tiles with no spectra_ , Forum Math. 18 (2006), no. 3, 519–528.
  * [11] I. Łaba, _Fuglede’s conjecture for a union of two intervals_ , Proc. Amer. Math. Soc. 129 (2001), no. 10, 2965–2972 (electronic).
  * [12] I. Łaba, _The spectral set conjecture and multiplicative properties of roots of polynomials_ , J. London Math. Soc. (2) 65 (2002), no. 3, 661–671.
  * [13] T. Y. Lam and K. H. Leung, _On vanishing sums of roots of unity_ , J. Algebra 224 (2000), no. 1, 91–109.
  * [14] Máté Matolcsi, _Fuglede’s conjecture fails in dimension 4_ , Proc. Amer. Math. Soc. 133 (2005), no. 10, 3021–3026 (electronic).
  * [15] Terence Tao, _Fuglede’s conjecture is false in 5 and higher dimensions_ , Math. Res. Lett. 11 (2004), no. 2-3, 251–258.

[RDM] Romanos Diogenes Malikiosis
Technische Universität Berlin, Institut für Mathematik
Sekretariat MA 4-1
Straße des 17. Juni 136
D-10623 Berlin, Germany
malikiosmathtu-berlinde
https://sites.google.com/site/romanosdiogenesmalikiosis/ [MNK] Mihail N.
Kolountzakis
Department of Mathematics and Applied Mathematics, University of Crete
Voutes Campus
700 13 Heraklion, Greece
kolountgmailcom
http://mk.eigen-space.org/


# Symmetric projected entangled-pair states analysis of a phase transition in
coupled spin-1/2 ladders

Juraj Hasik  Institute for Theoretical Physics, University of Amsterdam,
Science Park 904, 1098 XH Amsterdam, The Netherlands  Laboratoire de Physique
Théorique, Université de Toulouse, CNRS, UPS, France  Glen Bigan Mbeng
Institut für Theoretische Physik, Universität Innsbruck, A-6020 Innsbruck,
Austria  Sylvain Capponi  Laboratoire de Physique Théorique, Université de
Toulouse, CNRS, UPS, France  Federico Becca  Dipartimento di Fisica,
Università di Trieste, Strada Costiera 11, I-34151 Trieste, Italy  Andreas M.
Läuchli  Laboratory for Theoretical and Computational Physics, Paul Scherrer
Institute, 5232 Villigen, Switzerland  Institute of Physics, École
Polytechnique Fédérale de Lausanne (EPFL), 1015 Lausanne, Switzerland
Institut für Theoretische Physik, Universität Innsbruck, A-6020 Innsbruck,
Austria

###### Abstract

Infinite projected entangled-pair states (iPEPS) have been introduced to
accurately describe many-body wave functions on two-dimensional lattices. In
this context, two aspects are crucial: the systematic improvement of the
Ansatz by the optimization of its building blocks, i.e., tensors characterized
by bond dimension $D$, and the extrapolation scheme to reach the
“thermodynamic” limit $D\to\infty$. Recent advances in variational
optimization and scaling based on correlation lengths demonstrated the ability
of iPEPS to capture the spontaneous breaking of a continuous symmetry in
phases such as the antiferromagnetic (Néel) phase with high fidelity, in
addition to valence-bond solids which are already well described by finite-$D$
iPEPS. In contrast, systems in the vicinity of continuous quantum phase
transitions still present a challenge for iPEPS, especially when non-abelian
symmetries are involved. Here, we consider the iPEPS Ansatz to describe the
continuous transition between the (gapless) antiferromagnet and the (gapped)
paramagnet that exists in the $S=1/2$ Heisenberg model on coupled two-leg
ladders. In particular, we show how accurate iPEPS results can be obtained
down to a narrow interval around criticality and analyze the scaling of the
order parameter in the Néel phase in a spatially anisotropic situation.

##  I Introduction

One of the most challenging problem in condensed matter theory is to obtain
sufficiently accurate approximations of the ground state and low-energy
excitations of generic many-body Hamiltonians. For this reason, several
numerical approaches have been devised in the last 30 years, including
density-matrix renormalization group (DMRG) [1], dynamical mean-field theory
[2], and quantum Monte Carlo (QMC) techniques [3]. Among them, the natural
extensions of DMRG, such as matrix-product states [4, 5] in one dimension and
projected entangled-pair states (PEPS) [6, 7] in two dimensions, represent a
promising computational framework to get excellent approximations of strongly-
correlated lattice models. Within this approach, the ground-state wave
function is represented by means of local tensors, typically associated to the
sites of the underlying lattice. These tensors have two kinds of indices: a
single physical index specifying the local physical configuration (e.g.,
$S^{z}=\pm 1/2$ in a spin-$1/2$ model) and a collection of bond indices (whose
number, usually, equals the the coordination number of the lattice), each one
having dimension $D$. This latter quantity controls the accuracy of the
Ansatz, the exact ground-state wave function being achieved for $D\to\infty$.
The bond indices on nearest-neighbor sites are contracted together, thus
defining a tensor network. Furthermore, infinite systems may be considered by
embedding the tensor network within a suitable environment. Such tensor
networks, dubbed iPEPS, are the main focus of this work.

In the last decade, several algorithms to optimize iPEPS Ansätze have been
developed and used to address challenging problems in different correlated
systems, ranging from the characterization of unconventional states of matter
(e.g., spin liquids in highly-frustrated spin model) to the competition
between stripes and superconductivity in Hubbard or $t{-}J$ models. For
example, the ground state of the spin-$1/2$ Heisenberg model has been analysed
on the kagome [8] and Shastry-Sutherland [9] lattices. Hamiltonians with an
enlarged $SU(N)$ “spin” symmetry have been also considered, to assess both
magnetically ordered and disordered phases [10, 11]. In this regard, a
particular attention has been devoted to the possibility to stabilize chiral
spin liquids [12, 13]. As far as electronic models are concerned, the evidence
in favor of stripes has been pushed for $t{-}J$ [14] and Hubbard [15] models,
in the vicinity of the hole doping $\delta=1/8$.

By construction, the iPEPS Ansatz has been devised to describe both gapped and
gapless phases, satisfying the area law [16]; however, within the optimization
procedure in generic systems, gapless states are not easily obtained.
Similarly, it is also not easy to correctly describe states in the vicinity of
continuous quantum phase transitions. This is in a stark contrast with what
happens in one spatial dimension, where tensor-network states, equipped with
an appropriate scaling theory, may capture very well gapless phases and
critical phenomena [17, 18, 19, 20]. Recently, a progress in the iPEPS
analysis of critical and/or gapless systems has been made thanks to two key
developments. The first one is the introduction of optimization techniques,
based either on diagram summations [21, 22] or the so-called automatic
differentiation [23], which substantially improve the accuracy with respect to
the commonly used imaginary time evolution methods. The second one is the
development of finite-correlation-length scaling (FCLS) [24, 25] that allowed
to leverage the well-established finite-size scaling approach to iPEPS states.
In fact, the accuracy of thermodynamic estimates based on finite-$D$ iPEPS
calculations can be considerably improved. The efficiency of these advances
was recently demonstrated on a paradigmatic problem of the spin-$1/2$
Heisenberg $J_{1}-J_{2}$ model on a square lattice by estimating the
magnetization curve within the Néel phase, even in the vicinity of the
transition to the quantum spin-liquid state [26].

In this work, we pursue the idea of describing a continuous phase transition
within the iPEPS formalism. In particular, as for the $J_{1}-J_{2}$ model, we
focus on a case where the Néel order is melted by quantum fluctuations, namely
a two-dimensional system of coupled spin-$1/2$ Heisenberg ladders defined by

${\cal H}=J\sum_{R}{\bf S}_{R}\cdot{\bf S}_{R+\hat{x}}+\sum_{R}J_{R}{\bf
S}_{R}\cdot{\bf S}_{R+\hat{y}},$ (1)

where ${\bf S}_{R}=(S^{x}_{R},S^{y}_{R},S^{z}_{R})$ is the $S=1/2$ operator on
the site $R=(x,y)$ of a square lattice, $\hat{x}$ and ${\hat{y}}$ are unit
vectors in $x$ and $y$ direction, and $J_{R}=J$ or $J_{R}=\alpha J$, depending
on the parity of $y$. By varying $\alpha$, this model interpolates between the
Heisenberg model on the square lattice at $\alpha=1$ and a system of decoupled
two-leg ladders at $\alpha=0$. In the following, we will take $J$ as the
energy scale. Owing to the absence of the sign problem, this system can be
studied by unbiased QMC techniques and, therefore, offers an excellent
benchmark for the accuracy of iPEPS to describe a non-trivial quantum phase
transition, e.g., beyond the simplest case of the Ising model in transverse
field. Indeed, the Hamiltonian (1) displays a quantum critical point at
$\alpha_{c}=0.31407(5)$, separating a gapless antiferromagnet and a gapped
paramagnet [27]. The critical exponents are compatible with the ones of the
classical three-dimensional $O(3)$ Heisenberg model, as expected.

A recent iPEPS investigation of this Heisenberg model [28] highlighted the
difficulties faced by the unrestricted optimizations of iPEPS tensors across
the phase transition. In fact, due to the unbroken $SU(2)$ symmetry in the
paramagnetic side, strong finite-$D$ effects are present, thus impeding a
systematic analysis. Moreover, the optimization procedure is problematic also
within the antiferromagnetic phase, where the expected $U(1)$ symmetry around
the direction of the staggered magnetization is usually broken. In this work,
we want to constrain the iPEPS Ansatz by imposing this $U(1)$ symmetry, thus
limiting finite-$D$ effect, and optimize such symmetric tensors with the
automatic differentiation [23]. In addition, we extend the FCLS analysis [24,
25] in situations with a spatial anisotropy, due to the presence of two length
scales, since the correlations along the spatial $x$ and $y$ directions will
be generically different in the ground state of Eq. (1). These two
improvements allow us to get accurate results for the antiferromagnetic order
parameter up to the critical point, as confirmed by a direct comparison with
QMC calculations. Finally, we included an external staggered magnetic field in
the Hamiltonian (1), which allows us to inspect the response in the order
parameter, even very close to the critical point.

##  II Methods

Here, we briefly describe the structure of the iPEPS Ansatz that is used to
investigate the Hamiltonian (1). In order to account for both the
antiferromagnetic order in the gapless regime and the short-range valence-bond
correlations in the gapped one, we consider four rank-5 tensors
$t=\\{a,b,c,d\\}$, each one with a physical index $s$ (with dimension $2$,
suitable for $S=1/2$ degrees of freedom) and four auxilliary indices $u,l,d,r$
(with dimension $D$). These tensors are arranged in a $2\times 2$ unit cell to
tile the square lattice:

$|\textrm{iPEPS}(a,b,c,d)\rangle=\vbox{\hbox{\includegraphics[scale={0.25}]{fig1.pdf}}}.$
(2)

Observables are obtained by computing effective environments within the
corner-transfer matrix method, generalized for extended unit cells [14].

The gapped phase retains all the symmetries of the Hamiltonian and, in
particular it has the full $SU(2)$ spin symmetry; by contrast, the gapless
phase breaks the spin-rotational symmetry, leading to a finite collinear
magnetization and a residual $U(1)$ symmetry about it. Therefore, we impose an
explicit $U(1)$ spin symmetry on the tensors. This can be achieved by
associating integer “charges” to the two physical spin-$1/2$ components and
the $D$ virtual degrees of freedom (on each of the bond indices) [28]. To
preserve point-group symmetries of the underlying lattice, a straightforward
choice is to take the same charges on all the bond indices that are equivalent
under the symmetries (e.g., the $C_{4v}$ symmetry used in Ref. [28] led to
such choice). Here, we do not want to impose any symmetry within the $2\times
2$ unit cell, which allows the possibility that different charges are defined
on different indices. Thus, each of the four tensors will possess a set of
charges ${\vec{v}}^{\gamma}_{j}=(v_{0},\dots,v_{D-1})$ with $\gamma=a,b,c,d$
and $j=u,l,d,r$. Instead, the same charges
${\vec{u}}=(u^{\uparrow},u^{\downarrow})$ for the physical indices are taken
for all tensors. Therefore, the $U(1)$ symmetry is realized by enforcing a
selection rule for the non-zero elements of the tensors $t^{s}_{uldr}$:

$u^{s}+v^{\gamma}_{u}+v^{\gamma}_{l}+v^{\gamma}_{d}+v^{\gamma}_{r}=N,$ (3)

where, without loss of generality, we chose the case with $N=0$ to work with
invariant tensors.

The definition of these charges allows us to improve the numerical efficiency,
since tensors may arranged in blocks with definite values of them, allowing a
block-sparse representation. For a detailed treatment of tensor networks with
$U(1)$ symmetry and linear algebra with block-sparse tensors, see Ref. [29].
Importantly, once the blocks are defined, the gradient optimization changes
the tensors elements without mixing different blocks. The implementation of
linear tensor algebra for abelian-symmetric tensors is provided by open-source
library YAST [30] and the iPEPS algorithms built on top of it are available in
the peps-torch library [31].

From a practical point of view, one possibility to identifying $U(1)$ charges
(for the best variational states) is to perform imaginary-time evolution at
different values of $\alpha$ by using a two-site simple-update scheme [32]
starting from product states with $U(1)$ symmetry. Remarkably, this approach
does not lead to substantial improvements (in terms of physical quantities, as
energy or magnetization) with respect to the charges obtained from the
unrestricted optimizations of single-site iPEPS Ansatz with the $C_{4v}$
lattice symmetry in the Heisenberg model with $\alpha=1$ [26], see Appendix A.

##  III Results

Let us now show the results obtained within the optimization of the iPEPS
Ansätze. In order to perform a careful comparison with QMC calculations (which
are numerically exact, since the spin model does not suffer from the sign
problem), we apply the stochastic series expansion method [33] and perform
extrapolations by decreasing the temperature (to reach ground-state
properties) and increasing the size of the cluster (to reach the thermodynamic
limit). Most of the QMC calculations were done on $L\times L$ lattices with
periodic boundary conditions, at temperature $T=1/(2L)$, by using the ALPS
libraries [34, 35].

###  III.1 Phase diagram

In Fig. 1, we report the energy per site $e$ and the staggered magnetization
$m^{2}$ (averaged over the $2\times 2$ unit cell), for different values of
$\alpha$ across the quantum critical point. The outcome shows that, for small
bond dimension, i.e., $D\leq 3$, the accuracy is quite poor in both the
antiferromagnetic and paramagnetic phases. However, once the bond dimension
becomes large enough to capture the entanglement structure of the ground
state, the tensor network provides an excellent variational description, not
only in the paramagnetic phase, where the iPEPS parametrization is
particularly suitable, but also within the magnetically ordered phase. The
accuracy of the ground-state energy is exemplified by showing the relative
error with respect to the value $e_{QMC}$, obtained within the QMC approach
(extrapolated in the thermodynamic limit). Indeed, the quantity $\Delta
e/|e_{QMC}|$, with $\Delta e=|e-e_{QMC}|$ is vanishing within the statistical
errors, see Fig. 1. Most importantly, the magnetization curves for finite
values of the bond dimension $D$ follow the QMC data faithfully in the bulk of
each phase, still overestimating the order parameter close to the critical
point. We would like to mention that the present QMC results allow us to
locate the critical point at $\alpha_{c}=0.31467(1)$, which is slightly
different from the previous estimate $\alpha_{c}=0.31407(5)$ [27] (see
Appendix B).

Figure 1:  Energy per site (upper panel) and staggered magnetization (lower
panel) for optimized iPEPS with $U(1)$ symmetry and $2\times 2$ unit cell. The
results for $D=3,\ldots,6$ are reported, as well as the ones obained with the
FCLS extrapolation, see text. In addition, QMC results are also shown for the
thermodynamic extrapolation. The inset in the upper panel shows the relative
error of the iPEPS energies (with FCLS extrapolations) with respect to the QMC
ones (with thermodynamic extrapolations). The inset in the bottom panel gives
a sketch of the super-exchanges in Hamiltonian (1), where intra-ladder $J$
(inter-ladder $\alpha J$) couplings are denoted by black solid (gray dashed)
lines.

###  III.2 Finite-correlation-length scaling

Obtaining accurate and reliable estimates of the observables (e.g., the energy
per site or the staggered magnetization) within iPEPS Ansätze represents an
important issue that has been addressed since their definition. Recently, an
interesting step forward has been obtained for magnetic states that break a
continuous symmetry (e.g., the antiferromagnetic one in the Heisenberg model).
In this case, the presence of gapless modes induces a diverging correlation
length. When working on finite clusters, as in standard QMC approaches, this
fact translates into well-defined size-scaling laws of the physical quantities
[36, 37]. Instead, within iPEPS, we directly work in the thermodynamic limit
(that is mimicked by the embedding procedure) and, therefore, it is not
possible to straightforwardly apply the same scaling laws; still, the
correlation length $\xi$ is finite for any finite value of the bond dimension
$D$ and diverges with $D\to\infty$. As a result, it has been suggested [24,
25] that the FCLS analysis can be defined in terms of $\xi$ instead of the
cluster size as

$\displaystyle e(\xi)$ $\displaystyle=$ $\displaystyle
e_{0}+\frac{a}{\xi^{3}}+O\left(\frac{1}{\xi^{4}}\right),$ (4) $\displaystyle
m^{2}(\xi)$ $\displaystyle=$ $\displaystyle
m^{2}_{0}+\frac{b}{\xi}+O\left(\frac{1}{\xi^{2}}\right),$ (5)

where $a$ and $b$ are suitable constants. The system we are interested in,
i.e., the Heisenberg model on coupled two-leg ladders, presents a
complication, as it breaks the $\pi/2$-rotational symmetry for $\alpha\neq 1$,
leading to two distinct length scales. Within finite-$D$ iPEPS calculations,
we have access to both of them by looking at the spectrum of the transfer
matrices along $x$ (i.e., within ladders) and $y$ (i.e., across ladders).
Thus, we obtain two correlation lengths ($\xi_{x}$ and $\xi_{y}$) and perform
FCLS independently for both of them. The extrapolated values $e_{0}$ and
$m^{2}$ are independent of choice of the correlation length direction, while
the parameters $a$ and $b$ depend on the choice. It must be emphasized that
these FCLS relations are expected to be valid exclusively for a Goldstone
phase with a spontaneously broken continuous symmetry, where the correlation
length diverges with increasing $D$. Instead, within the paramagnetic gapped
phase, $\xi$ remains finite in the thermodynamic limit, thus implying that
Eqs. (4) and (5) cannot be applied. The impossibility to fit the numerical
data for a given $\alpha$ within this scheme is then taken as the evidence
that the state is not gapless (antiferromagnetic). At the quantum critical
point, we still expect Eq. (4) to hold, while the order parameter scaling is
replaced with a quantum critical form $m(\xi)\approx\xi^{-\Delta_{m}}$, with
$\Delta_{m}$ the scaling dimension of the order parameter at the critical
point [24].

The results for the energy per site $e$ are reported for selected values of
$\alpha$ in Fig. 2. Within the magnetically ordered phase (but also close to
the quantum critical point), the extrapolations of the FCLS analysis using
either $\xi_{x}$ or $\xi_{y}$ give consistent values, providing an internal
verification of the approach. Most importantly, the staggered magnetization
$m^{2}$ can be also extracted for infinite correlation length, see Fig. 3. The
extrapolated results are in very good agreement with QMC estimates, up to
values of the inter-ladder couplings that are very close to the critical
point. Beyond the critical point, the system develops a gap and thus the FCLS
fails to provide the correct scaling, where correlation lengths remain finite.
Nevertheless, within the paramagnetic regime, finite-$D$ iPEPS calculations
are able to recover accurate results. Indeed, for $\alpha\lesssim 0.28$, a
vanishing magnetization is already obtained with $D=5$. However, upon
approaching criticality, i.e., at $\alpha=0.3$, finite-$D$ iPEPS retain a
small residual magnetization (at least up to $D=6$). In fact, close to the
critical point, the correlation length far surpasses the values of $\xi_{x}$
and $\xi_{y}$ that our finite-$D$ iPEPS can support. Hence, the FCLS analysis
can be used to improve the estimates of the energy. Still, some important
deviations in the scaling laws are clearly visible, see Fig. 3. Finally, let
us remark that the iPEPS optimization in the narrow region near criticality is
technically difficult due to presence of instabilities (see Appendix C).

Figure 2:  FCLS extrapolations for the energy $e$, by using both horizontal
($\xi_{x}$, full symbols and lines) and vertical ($\xi_{y}$, empty symbols and
dashed lines) correlation lengths, for a few values of the inter-ladder
couplings $\alpha$. The symbols correspond to iPEPS with $D=3$ (diamonds),
$D=4$ (squares), $D=5$ (triangles), and $D=6$ (pentagons). The colors follow
the $\alpha$ color scheme of Fig. 3. Figure 3:  FCLS extrapolations for
staggered magnetization $m^{2}$, by using both vertical ($\xi_{y}$, left part)
and horizontal ($\xi_{x}$, right part) correlation lengths for
$0.29\leq\alpha\leq 1$. The symbols correspond to iPEPS with $D=3$ (diamonds),
$D=4$ (squares), $D=5$ (triangles), and $D=6$ (pentagons); thermodynamic
extrapolations from QMC data are also reported (stars). Figure 4:
$U(1)$-symmetric iPEPS with $D=3,\ldots,6$ in the presence of external
staggered field $h$. Upper panels: Running exponent as a function of $h$ for
the gapped phase with $\alpha=0.2$ (left), at criticality $\alpha=0.31407$
(middle), and in the Néel phase $\alpha=0.4$ (right). The value of the
critical exponent $1/\delta$ is reported at the critical point. QMC results
are also reported for comparison, with faded colors marking the data affected
by finite lattice size. Bottom panels: horizontal correlation length $\xi_{x}$
as a function of $h$.

###  III.3 The effect of a staggered magnetic field

In the previous section, we focused on estimating observables of the ground
state of the $SU(2)$-symmetric Hamiltonian (1). Here, we explicitly break this
symmetry to investigate the phase transition with a complementary non-
symmetric approach [38, 39, 24]. In particular, we follow the approach of
Refs. [24, 39], which recently used tensor networks to compute running
exponents in critical quantum systems.

To explicitly break the $SU(2)$, we supplement the Hamiltonian (1) with an
external staggered field $h$, which directly couples to the order parameter:

${\cal H}_{h}={\cal H}-h\sum_{R}(-1)^{x+y}S^{z}_{R}\,,$ (6)

where $R=(x,y)$. Then, for any $h\neq 0$, the ground state has a finite
correlation length and develops a finite staggered magnetization
$m=m(h,\alpha)$ in response to the staggered field $h$. The response for $h\to
0^{+}$ is different within magnetically ordered or disordered phases. In both
cases, the magnetization is an analytic function of $h$, with two distinct
regimes:

$\displaystyle m(\alpha,h)$ $\displaystyle=$ $\displaystyle
a(\alpha)h+O(h^{2})\qquad\mathrm{for}\;\alpha<\alpha_{c},\,h\to 0^{+},$ (7)
$\displaystyle m(\alpha,h)$ $\displaystyle=$ $\displaystyle
m(\alpha,0)+O(h)\qquad\mathrm{for}\;\alpha>\alpha_{c},\,h\to 0^{+},$ (8)

where $a(\alpha)$ is a suitable constant. By contrast, at criticality (i.e.,
for $\alpha=\alpha_{c}$), the magnetization is not analytic and the response
follows a power-law behavior:

$m(\alpha_{c},h)\propto h^{1/\delta}\;,\;\;\mathrm{for}\;h\to 0^{+}\,,$ (9)

where $\delta$ is a critical exponent, which only depends on the universality
class of the phase transition. The best estimate of $\delta$ (within the
expected universality class of the classical three-dimensional $O(3)$
Heisenberg model) $1/\delta=0.20916$ [40].

To probe the system’s response at fixed $\alpha$, we define the logarithmic
derivative

$[1/\delta](\alpha,h)=\frac{\partial\log m(\alpha,h)}{\partial\log h}\,,$ (10)

which is usually referred to as running exponent. For each value of $h$, we
optimize the $U(1)$-symmetric iPEPS Ansatz and compute the average staggered
magnetization $m(\alpha,h)$. Then, we estimate the logarithmic derivative by
using finite differences. The QMC results are based on a direct improved
estimator for the running exponent, see Ref. [39]. The outcomes are shown in
Fig. 4, where the running exponent and the correlation length $\xi_{x}$ of
iPEPS are reported for the three different regimes (with $D=2,\ldots,6$).

In the gapped phase (for $\alpha=0.2$), the magnetization is linear in $h$
(linear response) and, therefore, the running exponent saturates at $1$ for
$h\to 0$, as also obtained numerically from both QMC and iPEPS data (with
$D\geq 4$). The iPEPS correlation length saturates for $h\lesssim 0.005$, with
small finite-$D$ corrections. Indeed, in the gapped phase, even finite-$D$
iPEPS provide precise description of the system as the area law is upheld,
which is supported by a direct comparison between iPEPS (with $D\geq 4$) and
QMC (for $L\geq 32$). The $D=3$ iPEPS cannot faithfully capture gapped phase
at $\alpha=0.2$, since it retains finite $m$ and thus responds more akin to
the Néel phase.

Within the magnetically ordered phase (for $\alpha=0.4$), the magnetization
saturates to a finite value and the running exponent goes to zero for $h\to
0$. Indeed, the logarithmic derivative $[1/\delta](\alpha,h)$ has a non-
monotonic behavior, with a broad peak at intermediate values of the staggered
field $h$. This feature is correctly captured by iPEPS. QMC calculations also
reproduce the broad peak, even though at small values of $h$ a huge upturn is
present due to size effects; therefore, for large enough system size $L$, the
two methods converge to the same curve for sufficiently large values of the
staggered field $h$. Here, it seems as though QMC simulations provide an upper
bound of the thermodynamic limit, while iPEPS provide a lower bound of it.
Still the presence of the broad peak is robust and appears in the region of
$h$ where both $D\geq 5$ iPEPS and $L\geq 32$ QMC data are in excellent
agreement.

At the critical point $\alpha=\alpha_{c}$, the running exponent is expected to
converge towards $1/\delta=0.20916$, as mentioned before. The available iPEPS
and QMC calculations (for $D\leq 6$ and $L\leq 128$, respectively) are still
somewhat away from the saturation regime, however the numerical data
monotonously increase, without the broad peak in the magnetically orderd
phase. While the validity of the QMC results is limited by the finite size
$L$, for iPEPS it is limited by the finite correlation length induced by the
bond dimension $D$. Within iPEPS, as the field $h$ is decreased, the
correlation lengths initially grow, but then flatten out (with large
finite-$D$ corrections), in contrast to the expected diverging behavior of the
gapless regime. From the direct comparison between QMC and iPEPS data, we can
estimate that the iPEPS response for $D=6$ gives a faithful estimation of the
exact result down to $h\approx 10^{-2}$. Then, for smaller fields, the running
exponent data for iPEPS and large-$L$ QMC simulations start to deviate as
iPEPS becomes increasingly biased by the induced finite correlation length.

##  IV Conclusions and outlook

The analysis of quantum phase transitions by variational methods relies on a
few aspects: first, a flexible Ansatz for the ground-state wave fuction,
allowing a description of different kinds of phases by tuning its parameters;
second, a practical optimization scheme of such parameters to get the best
approximate ground state; finally, a way to analyze the results, possibly
extrapolating to the thermodynamic limit. Here, we have shown that symmetric
tensor networks, when implemented with a variational optimization, allow for
accurate description of the transition between a gapped paramagnet and gapless
Néel antiferromagnet, thus overcoming the issues arising within the imaginary-
time evolution approach that has been emphasized in a recent work [28]. The
selection of the correct symmetry structure for tensors by the choice of their
charges is crucial, as they determine the physical properties of the wave
function. Most importantly, once the tensors have been variationally
optimized, our simulations show that the physical observables, such as the
energy or the order parameter, are robust to small variations of these
charges. Although the region close to the critical point remains a challenge,
it is in principle tractable by increasing the bond dimension $D$. Already
with data up to $D=6$ and FCSL scaling, we could locate the critical point
with an accuracy of about $5\%$, demonstrating the applicability of this
analysis even in the case of two different length scales $\xi_{x}$ and
$\xi_{y}$. Finally, inspired by the recent works of Refs. [24, 39], we have
included an external staggered field $h$ in the Hamiltonian and compared the
iPEPS results with the ones obtained by QMC. Away from the critical point
$\alpha>\alpha_{c}$ (i.e., within the Néel phase) the running exponent
$[1/\delta](\alpha,h)$ shows a broad peak at finite field $h$, indicative of
magnetic order [39]. The robustness of the peak can be established by locating
the inflection point in the growing iPEPS correlation length as the external
field $h$ is decreased. Peaks at fields larger than inflection point represent
genuine features of the system, while finite-$D$ effects may generate spurious
peaks at fields smaller than the inflection point. Our findings corroborate
the fact that the analysis of the running exponent is a useful diagnosis
within the antiferromagnetic regime, away from the critical point. This proof
of concept analysis shows the potential of the method put forward in Ref. [39]
for future application of iPEPS, e.g., on questions of the stability of
quantum spin liquids.

###### Acknowledgements.

J.H. thanks D. Poilblanc and P. Corboz for valuable discussions. This project
has received funding from the European Research Council (ERC) under the
European Union’s Horizon 2020 research and innovation programme (grant
agreement No 101001604), TNSTRONG ANR-16-CE30-0025, TNTOP ANR-18-CE30-0026-01
grants awarded from the French Research Council, the Austrian Science Fund FWF
project I-4548, and the SFB BeyondC Project No. F7108-N38. This work was also
granted access to the HPC resources of CALMIP supercomputing center under the
allocations 2017-P1231 and 2021-P0677.

## References

  * White [1992] S. White, Density matrix formulation for quantum renormalization groups, Phys. Rev. Lett. 69, 2863 (1992).
  * Georges and Kotliar [1992] A. Georges and G. Kotliar, Hubbard model in infinite dimensions, Phys. Rev. B 45, 6479 (1992).
  * White  _et al._ [1989] S. White, D. Scalapino, R. Sugar, E. Loh, J. Gubernatis, and R. Scalettar, Numerical study of the two-dimensional Hubbard model, Phys. Rev. B 40, 506 (1989).
  * Östlund and Rommer [1995] S. Östlund and S. Rommer, Thermodynamic limit of density matrix renormalization, Phys. Rev. Lett. 75, 3537 (1995).
  * Perez-Garcia  _et al._ [2007] D. Perez-Garcia, F. Verstraete, M. Wolf, and J. Cirac, Matrix product state representations (2007), arXiv:quant-ph/0608197 [quant-ph] .
  * Verstraete and Cirac [2004a] F. Verstraete and J. Cirac, Renormalization algorithms for quantum-many body systems in two and higher dimensions (2004a), arXiv:cond-mat/0407066 [cond-mat.str-el] .
  * Verstraete and Cirac [2004b] F. Verstraete and J. Cirac, Valence-bond states for quantum computation, Phys. Rev. A 70, 060302 (2004b).
  * Liao  _et al._ [2017] H. Liao, Z. Xie, J. Chen, Z. Liu, H. Xie, R. Huang, B. Normand, and T. Xiang, Gapless spin-liquid ground state in the $s=1/2$ kagome antiferromagnet, Phys. Rev. Lett. 118, 137202 (2017).
  * Corboz and Mila [2013] P. Corboz and F. Mila, Tensor network study of the Shastry-Sutherland model in zero magnetic field, Phys. Rev. B 87, 115144 (2013).
  * Corboz  _et al._ [2011] P. Corboz, A. Läuchli, K. Penc, M. Troyer, and F. Mila, Simultaneous dimerization and su(4) symmetry breaking of 4-color fermions on the square lattice, Phys. Rev. Lett. 107, 215301 (2011).
  * Corboz  _et al._ [2012] P. Corboz, M. Lajkó, A. Läuchli, K. Penc, and F. Mila, Spin-orbital quantum liquid on the honeycomb lattice, Phys. Rev. X 2, 041013 (2012).
  * Chen  _et al._ [2020] J.-Y. Chen, S. Capponi, A. Wietek, M. Mambrini, N. Schuch, and D. Poilblanc, SU(3)1 chiral spin liquid on the square lattice: A view from symmetric projected entangled pair states, Phys. Rev. Lett. 125, 017201 (2020).
  * Chen  _et al._ [2021] J.-Y. Chen, J.-W. Li, P. Nataf, S. Capponi, M. Mambrini, K. Totsuka, H.-H. Tu, A. Weichselbaum, J. von Delft, and D. Poilblanc, Abelian su(n)1 chiral spin liquids on the square lattice, Phys. Rev. B 104, 235104 (2021).
  * Corboz  _et al._ [2014] P. Corboz, T. Rice, and M. Troyer, Competing states in the $t$-$j$ model: Uniform $d$-wave state versus stripe state, Phys. Rev. Lett. 113, 046402 (2014).
  * Zheng  _et al._ [2017] B.-X. Zheng, C.-M. Chung, P. Corboz, G. Ehlers, M.-P. Qin, R. Noack, H. Shi, S. White, S. Zhang, and G.-L. Chan, Stripe order in the underdoped region of the two-dimensional Hubbard model, Science 358, 1155 (2017), https://www.science.org/doi/pdf/10.1126/science.aam7127 .
  * Eisert  _et al._ [2010] J. Eisert, M. Cramer, and M. Plenio, Colloquium: Area laws for the entanglement entropy, Rev. Mod. Phys. 82, 277 (2010).
  * Pollmann  _et al._ [2009] F. Pollmann, S. Mukerjee, A. Turner, and J. Moore, Theory of finite-entanglement scaling at one-dimensional quantum critical points, Phys. Rev. Lett. 102, 255701 (2009).
  * Pirvu  _et al._ [2012] B. Pirvu, G. Vidal, F. Verstraete, and L. Tagliacozzo, Matrix product states for critical spin chains: Finite-size versus finite-entanglement scaling, Phys. Rev. B 86, 075117 (2012).
  * Rams  _et al._ [2018] M. Rams, P. Czarnik, and L. Cincio, Precise extrapolation of the correlation function asymptotics in uniform tensor network states with application to the Bose-Hubbard and XXZ models, Phys. Rev. X 8, 041033 (2018).
  * Vanhecke  _et al._ [2019] B. Vanhecke, J. Haegeman, K. Van Acoleyen, L. Vanderstraeten, and F. Verstraete, Scaling hypothesis for matrix product states, Phys. Rev. Lett. 123, 250604 (2019).
  * Corboz [2016] P. Corboz, Variational optimization with infinite projected entangled-pair states, Phys. Rev. B 94, 035133 (2016).
  * Vanderstraeten  _et al._ [2016] L. Vanderstraeten, J. Haegeman, P. Corboz, and F. Verstraete, Gradient methods for variational optimization of projected entangled-pair states, Phys. Rev. B 94, 155123 (2016).
  * Liao  _et al._ [2019] H.-J. Liao, J.-G. Liu, L. Wang, and T. Xiang, Differentiable programming tensor networks, Phys. Rev. X 9, 031041 (2019).
  * Rader and Läuchli [2018] M. Rader and A. Läuchli, Finite Correlation Length Scaling in Lorentz-Invariant Gapless iPEPS Wave Functions, Phys. Rev. X 8, 031030 (2018).
  * Corboz  _et al._ [2018] P. Corboz, P. Czarnik, G. Kapteijns, and L. Tagliacozzo, Finite Correlation Length Scaling with Infinite Projected Entangled-Pair States, Phys. Rev. X 8, 031031 (2018).
  * Hasik  _et al._ [2021] J. Hasik, D. Poilblanc, and F. Becca, Investigation of the Néel phase of the frustrated Heisenberg antiferromagnet by differentiable symmetric tensor networks, SciPost Phys. 10, 12 (2021).
  * Matsumoto  _et al._ [2001] M. Matsumoto, C. Yasuda, S. Todo, and H. Takayama, Ground-state phase diagram of quantum Heisenberg antiferromagnets on the anisotropic dimerized square lattice, Phys. Rev. B 65, 014407 (2001).
  * Hasik and Becca [2019] J. Hasik and F. Becca, Optimization of infinite projected entangled pair states: The role of multiplets and their breaking, Phys. Rev. B 100, 054429 (2019).
  * Singh  _et al._ [2011] S. Singh, R. Pfeifer, and G. Vidal, Tensor network states and algorithms in the presence of a global u(1) symmetry, Phys. Rev. B 83, 115125 (2011).
  * [30] G. Rams, M.M. Wójtowicz and J. Hasik, YAST - Yet Another Symmetric Tensor: Differentiable linear algebra with block-sparse tensors, supporting abelian symmetries, https://gitlab.com/marekrams/yast.
  * [31] J. Hasik and G. Mbeng, peps-torch: A differentiable tensor network library for two-dimensional lattice models, https://github.com/jurajHasik/peps-torch.
  * Jiang  _et al._ [2008] H. Jiang, Z. Weng, and T. Xiang, Accurate determination of tensor network state of quantum lattice models in two dimensions, Phys. Rev. Lett. 101, 090603 (2008).
  * Sandvik [1999] A. Sandvik, Stochastic series expansion method with operator-loop update, Phys. Rev. B 59, R14157 (1999).
  * Albuquerque  _et al._ [2007] A. Albuquerque, F. Alet, P. Corboz, P. Dayal, A. Feiguin, S. Fuchs, L. Gamper, E. Gull, S. Gürtler, A. Honecker, R. Igarashi, M. Körner, A. Kozhevnikov, A. Läuchli, S. Manmana, M. Matsumoto, I. McCulloch, F. Michel, R. Noack, G. Pawłowski, L. Pollet, T. Pruschke, U. Schollwöck, S. Todo, S. Trebst, M. Troyer, P. Werner, and S. Wessel, The ALPS project release 1.3: Open-source software for strongly correlated systems, Journal of Magnetism and Magnetic Materials 310, 1187 (2007).
  * Bauer  _et al._ [2011] B. Bauer, L. Carr, H. Evertz, A. Feiguin, J. Freire, S. Fuchs, L. Gamper, J. Gukelberger, E. Gull, S. Guertler, A. Hehn, R. Igarashi, S. Isakov, D. Koop, P. Ma, P. Mates, H. Matsuo, O. Parcollet, G. Pawłowski, J. Picon, L. Pollet, E. Santos, V. Scarola, U. Schollwöck, C. Silva, B. Surer, S. Todo, S. Trebst, M. Troyer, M. Wall, P. Werner, and S. Wessel, The ALPS project release 2.0: open source software for strongly correlated systems, Journal of Statistical Mechanics: Theory and Experiment 2011, P05001 (2011).
  * Neuberger and Ziman [1989] H. Neuberger and T. Ziman, Finite-size effects in Heisenberg antiferromagnets, Phys. Rev. B 39, 2608 (1989).
  * Fisher [1989] D. Fisher, Universality, low-temperature properties, and finite-size scaling in quantum antiferromagnets, Phys. Rev. B 39, 11783 (1989).
  * Binder and Landau [1984] K. Binder and D. Landau, Finite-size scaling at first-order phase transitions, Phys. Rev. B 30, 1477 (1984).
  * D’Emidio  _et al._ [2021] J. D’Emidio, A. A. Eberharter, and A. M. Läuchli, Diagnosing weakly first-order phase transitions by coupling to order parameters (2021).
  * Chester  _et al._ [2021] S. Chester, W. Landry, J. Liu, D. Poland, D. Simmons-Duffin, N. Su, and A. Vichi, Bootstrapping Heisenberg magnets and their cubic instability, Phys. Rev. D 104, 105013 (2021).
  * Press  _et al._ [2007] W. Press, S. Teukolsky, W. Vetterling, and B. Flannery,  _Numerical Recipes, The Art of Scientific Computing, 3rd Edition_ (Cambridge University Press, Cambridge, 2007).

##  Appendix A Comparison among $U(1)$ classes

First of all, we discuss the preparation of the initial $U(1)$-symmetric
states for the optimization by imaginary-time evolution. We consider two
scenarios of evolution by starting from simple $U(1)$-symmetric product
states: a classical Néel antiferromagnet (NEEL) and a product state of
singlets placed on the rungs of ladders (VBS). The imaginary-time evolutions
were performed by using two-site simple update (SU) with the second-order
Trotter approximation. To track the convergence of SU, we evaluate energy
using corner-transfer matrix environments, with modest environment dimension
$\chi\approx D^{2}$. Starting with the time step $\delta\tau=0.02$, the state
is evolved until an energy increase is observed, in which case the SU step is
not performed and the time step is halved $\delta\tau\rightarrow\delta\tau/2$.
The evolution ends once the time step becomes smaller than
$\delta\tau=10^{-8}$.

The $U(1)$-structure of the evolved NEEL and VBS states show slight variations
in charges $\vec{v}^{\gamma}_{j}$. Moreover, depending on bond dimension $D$
and $\alpha$, the charges do not necessarily respect spatial symmetries of the
model (1). For example, at $D=5$ and $6$ the charges on left ($j=l$) and right
($j=r$) bonds of tensors come out different. Then, we perform a further
optimization of the NEEL and VBS states by using gradient descent and compare
the resulting physical observables with the ones obtained from optimal states
with identical charges $\vec{v}^{\gamma}_{j}$ for all bonds $j$), here denoted
as U1B class. The results for the energy per site $e$ and the local
magnetization $m$ are shown in Fig. 5. The data for $D=3$ show that optimized
VBS states have substantially worse energies than NEEL and U1B states.
Instead, for $D>3$, all these states have very similar energies, differences
being at most of order $2\times 10^{-4}$. The order parameter displays a
similar behavior. A minor exception is represented by the $D=4$ VBS case,
which displays slight breaking of translation symmetry along $y$-axis even in
the $\alpha=1$ limit. In particular, the nearest-neighbor spin-spin
correlations show staggered pattern with
$\langle\vec{S}_{x,y}\cdot\vec{S}_{x,y+1}\rangle-\langle\vec{S}_{x,y+1}\cdot\vec{S}_{x,y+2}\rangle\approx
3\times 10^{-4}$. Overall, the SU-evolved states provide reasonable accurate
initial states for the variational optimization of $U(1)$-symmetric iPEPS.
Irrespective of the differences in the $U(1)$ charges, the final variational
minima for different classes give quantitatively the similar physical picture.

Figure 5:  Left panels: energy per site obtained optimizing NEEL and VBS
states with respect to the one obtained within the U1B class, i.e., $\delta
e=|e({\it NEEL})-e({\it U1B})|$ (circles) and $\delta e=|e({\it VBS})-e({\it
U1B})|$ (squares) for different values of the bond dimension $D$. Right
panels: magnetization $m$ obtained from NEEL, VBS, and U1B states. Observables
are evaluated at finite $\chi$, i.e., $\chi=71$ for $D=3$, $\chi=64$ for
$D=4$, $\chi=50$ for $D=5$ and $\chi=72$. Figure 6:  QMC results for the
scaling analysis of the spin stiffness $\rho_{s}$ as a function of $\alpha$
for system with $L_{x}=L_{y}$ (upper panel) and $L_{x}=2L_{y}$ (lower panel).
The temperature has been chosen as $T=1/(2L_{x})$ to ensure convergence to the
ground-state properties. The vertical black line and shaded area indicate
$\alpha_{c}=0.31475(10)$ estimate and its uncertainty.

Figure 7:  Upper panel: example of optimizations without line search at
$\alpha=0.28$ (with $D=5$ and $\chi=75$) both with and without explicit $U(1)$
symmetry. Optimizations are initialized from the $\alpha=0.3$. The red line
shows the comparison with QMC extrapolation. Lower panel: difference between
thermodynamic estimate of the energy ($\chi\to\infty$ limit) and its minimal
value for fixed environment dimension $\chi$ for optimized $U(1)$ symmetric
states. The optimizations were done for $D=3$, $4$, $5$, and $6$ at
$\chi_{opt}=54$, $64$, $100$, and $108$, respectively. In all cases, the line
search was employed.

##  Appendix B Size scaling analysis of the transition by the QMC technique

Here, we provide some standard details on locating the transition point
$\alpha_{c}$ from QMC simulations. To estimate $\alpha_{c}$, we have computed
the spin stiffness $\rho_{s}$, defined using the winding number:

$\rho_{s}=\frac{T}{2}\langle\sum_{i}(W_{i}^{x})^{2}+(W_{i}^{y})^{2}\rangle$
(11)

where $i$ runs over loops in the space-time SSE configurations and the
brackets denote the Monte Carlo average.

The results for different system sizes $L_{x}\times L_{y}$ are shown in Fig 6.
At criticality, the scaled spin stiffness $L_{x}\rho_{s}$ is constant (taking
the dynamical exponent $z=1$ at the transition), thus allowing to estimate the
location of the transition at $\alpha_{c}=0.31467(1)$. Note that the Binder
cumulant crossing also leads to a similar critical value (data not shown). In
addition, we have also checked that data collapse of $\rho_{s}$ and
magnetization are compatible with known exponents $\nu=0.707$ and
$\beta=0.3689$.

##  Appendix C Instability in optimization

The optimization of the iPEPS wave function is performed by using the
steepest-descent algorithm, implemented with the line search method (i.e., by
using a step length that is adapted along the optimization procedure) [41]. In
the proximity of the critical point, performing the line search is crucial to
obtain stable variational optimizations, especially for the cases with $D=5$
and $6$. In Fig. 7 (upper panel), we show typical cases done without line
search, such that unstable regimes, with an erratic behavior of the energy,
are obtained. The same kind of problem appears in simulations with or without
the $U(1)$ symmetry.

When performing the optimization, we always work at constant environment bond
dimension $\chi_{opt}$. The states obtained in the erratic regime can show
energies $e(\chi_{opt})$, which are lower than the reference QMC results. This
is possible because the iPEPS thermodynamic estimate of energy, for which the
variational principle holds, is obtained only in the limit of $\chi\to\infty$.
When the line search is employed, it suppresses the erratic regime as the
energy is not allowed to increase in course of optimization. In this case, we
asses the severity of this instability by comparing the lowest energies that
are realized in the vicinity of $\chi_{opt}$ with the energy obtained from the
$\chi\rightarrow\infty$ limit. The results for $D=3,\ldots,6$ and
$0.2\leq\alpha\leq 0.4$ are shown in Fig. 7 (lower panel). This behavior
resembles the overtraining, often encountered in the optimization of
artificial neural networks, since the variational optimization of iPEPS does
not directly optimize the thermodynamic estimate of energy
$e(\chi\rightarrow\infty)$, but only its finite-$\chi$ approximation. In the
previous study of the $J_{1}-J_{2}$ model (using highly-constrained single-
site iPEPS Anstaz [26]), $e(\chi)$ was a monotonically decreasing function
already from modest values of $\chi$. Here, we instead observe that $e(\chi)$
attains shallow spurious minima, with depth at most $O(10^{-5})$, compared to
corresponding $e(\chi\rightarrow\infty)$ thermodynamic estimates.