Nonlocal order parameter of pair superfluids

\documentclass{article}
\usepackage{graphicx}
\usepackage{amsmath}
\usepackage[margin=1in]{geometry}
\usepackage{xcolor}
\usepackage{comment}
\usepackage{hyperref}
\usepackage{cite}

\title{Reply Paper Scipost}
\author{Nitya Cuzzuol, Arianna Montorsi and Luca Barbiero}
\date{July 2024}

\begin{document}

\maketitle

\section{Report1}
\textit{it would be useful to core a review on non-standard Hubbard models where pair supefluidity is discussed thoroughly in various context:}\\
\begin{itemize}
\item 1. Omjyoti Dutta, Mariusz Gajda, Philipp Hauke, Maciej Lewenstein, Dirk-Sören Lühmann, Boris A. Malomed, Tomasz Sowiński, Jakub Zakrzewski, Non-standard Hubbard models in optical lattices, Rep. Prog. Phys. 78, 066001 (2015), arXiv:1406.0181.
\item 2. Titas Chanda, Luca Barbiero, Maciej Lewenstein, Manfred J. Mark, and Jakub Zakrzewski, Recent progress on quantum simulations of non-standard Bose-Hubbard modes, arXiv:2405.07775.
\end{itemize}
We thank the referee for the excellent comments and evaluation of our paper, and for recommending it for publication. Following the suggestion, we introduced the citation to both reviews in the introduction of our paper.

\section{Report2}
\textit{I do not find the paper suitable for publication in its present form, because the text seems to lack accuracy. I suggest a major revision along the lines presented below.}\\
We thank the referee for carefully reading our manuscript. From his/her observations we understand that some of the known results obtained in previous literature could be further detailed in the present manuscript to render its reading more self contained. We appreciate the suggestion.
\begin{itemize}
\item (A) \textit{Could the Authors explain in greater detail how their quantity $O_P^{(o)}(j)$, defined in Eq. 8, provides different information from the closely related parity $O_P^{(o)}$, previously defined in Ref. [10:Berg PRB 2008], and explicitly mentioned on p. 4, 3 lines from the bottom of the page ?}\\
Parity and odd parity are related by
\begin{equation}
O_P(j) = (-1)^{j-1}O_P^{(o)}(j)
\end{equation}
As explained in more details in Ref.[11] for the very same relation holding between charge and spin parity in fermionic case, the finiteness of the expectation value of the uniform part of one or the other captures phases with different orders. In particular, $O_P$ turns out to have a finite uniform part $(O_P(j)+O_P(j+1))/2$ in the MI phase, in which the uniform part of $O_P^{(o)}$ is vanishing. Instead the latter becomes finite in the PSF phase, in which the uniform part of $O_P$ goes to zero. We added a sentence on this point in the revised manuscript.
\item (B) \textit{The Authors write (abstract, l. 7; and again p. 2, 4 lines from end) that the odd parity is the "unique order parameter" of the PSF phase. However, Ref. [46: Bonnes PRL 2011] also introduce the order parameter $\langle b_i^2 \rangle$ in a similar context. What are the differences in the roles of these two quantities ?}\\
First of all, in the present context in which the total number of particle is conserved, $\langle b_i^2 \rangle$ cannot be an order parameter. Indeed, $b_i^2$ correlations in PSF phase decay to zero algebraically in the thermodynamic limit [44], unlike our order parameter which correlations remain finite, as expected for a true order parameter. Secondly, while such algebraic decay for $b_i^2$ correlations is observed in both the SF and the PSF phases, $O_P^{(o)}$ is vanishing in the superfluid phase and becomes finite only when entering the paired superfluid phase.
\end{itemize}
\textit{The following remarks (1-4) are more minor, but may help in making the manuscript more accessible to a non-specialised audience.}
\begin{itemize}
\item 1) \textit{Concerning the second considered model (Eq. 11): Figure 6a seems to suggest that the species A and B are trapped in two spatially-separated optical lattices, but this does not explicitly appear in the Hamiltonian:}
\begin{equation}
H = - t \sum_{i, \sigma = \rm {A,B}} (b_{i,\sigma}^{\dagger} b_{i+1,\sigma} + \text{h.c.}) + U \sum_i n_{i,\rm {A}} n_{i,\rm {B}} %+ V \sum_{i,\sigma=A,B} n_{i,\sigma}n_{i+1,\sigma}
\end{equation}
\textit{does this subtlety in the proposed realisation play a role ?}\\
The model we are analyzing here can be described as a couple of chains of hard-core bosons or equivalently as a single chain populated by two species of hard-core bosons. $t$ is the hopping coefficient of hard core bosons and there is no interchain/interspecies hopping. The interaction between the two species/two chains is regulated by $U$. In order to further clarify the notation we added a comment on the index $i$ in the text.\\

\textit{Why do the paired bosons belong to sites of the two lattices which differ by one lattice spacing rather than to corresponding sites ?}\\
A pair of atoms on different chains and same site index corresponds to the double occupied site of the constrained BH model. The two atoms displaced by one lattice site in the two chains are representative of a broken pair maintaining a finite correlation length. In fig.6a the red ellipses highlight such finite correlation length for the two single bosons, which is consistent with finite odd parity.\\
We have added a comment on this point in the caption of fig.6.
\item 2) \textit{The Authors use two different methods: bosonisation (Sec. 2.1) and DMRG (Sec. 3). Matrix Product States are also mentioned once in Sec. 3.2. Would the Authors consider including a section entitled "Methods", in which they discuss the strengths and limitations of these approaches ?}\\
Since they are both standard methods in the context of low dimensional strongly correlated systems, we don't feel that an introduction in this context is appropriate. Notable examples of Bosonization analysis can be found in \cite{Schulz_1986, Berg_2008, DallaTorre_2006, Timonen_1985, Luther_1977, Orignac_1998, Shelton_1996, Nersesyan_2020, Nakamura_2003}. While DMRG results in \cite{Kuhner_1998, Kuhner_2000, DallaTorre_2006, Berg_2008, Ejima_2015, Ejima_2021, Rossini_2012, Dalmonte_2011, Diehl_2010, DegliEspostiBoschi_2003, DegliEspostiBoschi_2016, Fraxanet_2022}.
\item 3) \textit{Two keywords seem to be used without a definition or reference:}
\begin{itemize}
\item (a) \textit{"snake-like path", on p. 7, paragraph 1, line 6;}\\
This keyword means the way one enumerate the sites of a 2D lattice to obtain a single chain. It can be found in the references cited for the TeNPy libraries Ref.[69]. We will add the citation after this keyword.
\item (b) \textit{"$C_4$ symmetric 2D limit", on p. 7, paragraph 3, line 2.}\\
We are referring here to the rotation symmetry of the square lattice by $\pi/2$ rotations. However, we decided to remove this sentence to avoid any confusion.
\end{itemize}
\item 4) \textit{ I am surprised at the Authors' usage of "normal superfluid" as a single keyword (abstract, l. 13; and again p. 2, paragraph 3, l. 9). Is this standard terminology? If so, could the Authors provide a reference using it?}\\
The authors thank for the suggestion and will replace "normal superfluid" with "atomic superfluid" as in Ref.[43, 46]. %\cite{Diehl_2010, Bonnes_2011}. %[41,44].

\end{itemize}

\section{Report3}
\textit{This paper proposes a non-local order parameter "odd parity", which is a modified version of "parity" operator introduced in earlier literature. The authors argue that it characterizes the pair superfluid (PSF) phase, and demonstrate numerically in 1 and 2 dimensions. Indeed the numerical results seem to confirm the picture, and the odd parity operator looks quite interesting and useful.}\\
We thank the referee for the thorough reading of the paper and the valuable observations and suggestions provided. Following the suggestions and addressing the questions raised, we have integrated and updated the paper as discussed below.
\begin{itemize}
\item 1) Parity and odd parity are related as
\begin{equation}
O_P(j) = (-1)^{j-1}O_P^{(o)}(j)
\end{equation}
\textit{Then how can they detect different phases? I suppose that, the hidden assumption is that the nonlocal order parameter has a definite sign, i.e. it does not change the sign for sufficiently large $j$. So if $O_P(j)$ is non-vanishing but is oscillating and changes sign like $(-1)^j$, we do not regard the parity to be long-ranged but the odd parity is long-ranged (and vice versa). Do you agree? In any case, please clarify.}\\
This is the point. One should look at the asymptotic behavior of the expectation value of the uniform part of the two operators. Indeed the uniform part of $O_P(j)$, namely $(O_P(j)+O_P(j+1))/2$, has vanishing expectation value in the asymptotic limit in the pair superfluid phase, unlike that of $O_P^{(o)}(j)$ . Ultimately, this is rooted in the spin-charge separation of the bosonized Hamiltonian, and is perfectly analogous to the connection between spin and charge parity discussed in Ref.[11]. We added a comment on that along the manuscript.\\

\textit{It is also very important to describe how the authors extracted the non-local order parameter $O_P^{(o)}$ from the numerical data $O_P^{(o)}(j)$ for various different $j$'s. This should be one of the essential information the authors are obliged to present.}\\
The nonlocal order parameter $O_P^{(o)}$ has been extracted from $O_P^{(o)}(j)$ by fixing $j$ to a sufficiently large value to ensure that, within the region of parameter space where the model exhibits the phase defined by this parameter (i.e., not close to the transition), the nonlocal order parameter has already converged to a finite constant. We verified that in the same region the parameter $O_P(j)$ has an oscillating behaviour with vanishing average.
\item 2) \textit{As the authors note, under the three-body constraint, the Bose-Hubbard is mapped to the S=1 system. In the 1D context (S=1 chain), in my understanding, the PSF is nothing but the "XY2" phase identified by Schulz in the seminal paper Ref. [55]. Although the authors do cite Ref. [55], I believe that it would be fair to mention explicitly that the 1D PSF phase was essentially discovered in Ref.[55].}\\
The authors thank for the suggestion. We have added the information.
\item 3) \textit{Related to the previous point, in Ref. [55] the Antiferromagnetic (Neel) phase was also identified. I believe that the odd parity operator is also non-vanishing in the antiferromagnetic phase, which corresponds to a charge-density wave state in the Bose-Hubbard context. So the non-vanishing odd parity alone does not imply that the system is in PSF phase. For 1D Bose-Hubbard model, according to Fig. 2(b), the authors numerically confirmed the power-law decay of $D(r)$. Moreover, Fig. 3 (b) indicates that $\Delta_2$ always vanishes for the range of $U/t$ studied. So I agree that the system is indeed in the PSF phase when the odd parity is non-vanishing within the model studied by the authors. (However, one should also be able to find the antiferromagnetic phase by modifying the Hamiltonian.)}\\
\begin{figure}[ht]
%
\centering
\includegraphics[scale=1]{decayNN.png}
%
\caption{\label{NN}
Analysis of $(-1)^{|i-j|}\times\langle (n_i-\Bar{n})(n_j-\Bar{n}) \rangle$ for different values of the onsite interaction $U$ and nearest neighbour interaction $V$ for 1D and 2D. The correlator has long range order only for finite positive $V$ corresponding to CDW \cite{Cuzzuol_2024}.}
%
\end{figure}
The intuition that the odd parity is different from zero also in the charge density wave phase equivalent to the cited antiferromagnetic phase is correct. Indeed, in ref.\cite{Cuzzuol_2024} this was observed for the extended BH model for positive values of nearest neighbour interaction. However, contrary to PSF, the CDW phase is SSB ordered phase, and, besides being characterized by a finite value of both string and odd parity nonlocal orders, can also be characterized by the local order parameter $(-1)^{|i-j|}\times\langle (n_i-\Bar{n})(n_j-\Bar{n}) \rangle$ \cite{Baldelli_2024}, where $n$ is the number operator and $\Bar{n}$ the filling. This local order parameter has long range order only in CDW as represented in fig.\ref{NN}.\\

\textit{For 2D, Fig. 4(b) might suggest that the system is PSF and not antiferromagnetic/CDW, but I do not feel the evidence is sufficient. How does $D(r)$ look like, for example at $U= - 20$ ?}\\
Also for more negative $U$ values $D(r)$ is observed to have the same power law decay. For clarity we reported in Fig.(\ref{Dr}) the decay of $D(r)$ for $U = -17, -18, -20$. As a counterproof, in 1D the transition to CDW state was observed in \cite{Cuzzuol_2024} only for positive values of nn interactions. Moreover one can verify that $(-1)^{|i-j|}\times\langle (n_i-\Bar{n})(n_j-\Bar{n}) \rangle$ doesn't have long range order for negative $U$ in Fig.\ref{NN}.
\begin{figure}[ht]
%
\centering
\includegraphics[scale=1]{decaybd2b2_N6legs.png}
%
\caption{\label{Dr}
Analysis of $D(r)$ as defined in the paper for $U = -17, -18, -20$. For all three values of interaction $D(r)$ is always decaying power law.}
%
\end{figure}

\item \textit{4) Please give the detailed derivation of Eq. (9) (bosonized expression of the odd parity operator), which is rather important for this paper.}\\
A detailed derivation has been added to the section on bosonization, together with appropriate references to the literature on the related derivation for spin parity.

\item 5) \textit{I have mostly considered the simplest case especially in 1D. Please also carefully re-examine the cases of the arbitrary filling and the two-species system, following the above comments.}\\
Yes, we have also added a comment on the 2D limit.
\end{itemize}

\bibliographystyle{iopart-num.bst}
\bibliography{myRef_reply.bib}
\end{document}