Some ground-state expectation values for the free parafermion Z(N) spin chain

TL;DR

This paper derives exact ground-state expectation values for the non-Hermitian Z(N) spin chain using the Hellmann-Feynman theorem, extending known results for the quantum Ising model to general N and analyzing boundary effects.

Contribution

It provides the first exact analytical expressions for certain ground-state expectation values in free parafermion Z(N) chains for arbitrary N.

Findings

Exact expectation values for sites deep inside the chain.

Numerical validation for N=3.

Analysis of boundary effects on expectation values.

Abstract

We consider the calculation of ground-state expectation values for the non-Hermitian Z(N) spin chain described by free parafermions. For N=2 the model reduces to the quantum Ising chain in a transverse field with open boundary conditions. Use is made of the Hellmann-Feynman theorem to obtain exact results for particular single site and nearest-neighbour ground-state expectation values for general N which are valid for sites deep inside the chain. These results are tested numerically for N=3, along with how they change as a function of distance from the boundary.

Equations

$$H_{\mathrm{f}}=-\sum _{j=1}^{L-1}{\sigma }_j^z{\sigma }_{j+1}^z-\lambda \sum _{j=1}^L{\sigma }_j^x.$$

$$H_{\mathrm{pf}}=-\sum _{j=1}^{L-1}{\sigma }_j{\sigma }_{j+1}^{†}-\lambda \sum _{j=1}^L{\tau }_j,$$

$${\sigma }_j$$

$${\tau }_j$$

$${\sigma }_{m,n}={\omega }^{m-1}{\delta }_{m,n},{\tau }_{m,n}={\delta }_{m,n+1},$$

$$\sigma \tau =\omega \tau \sigma ,{\tau }^{†}={\tau }^{N-1},{\sigma }^{†}={\sigma }^{N-1},$$

$$M_k=(1-{\lambda }^2{)}^{k(N-k)/(2N^2)},$$

$$E=\pm {ϵ}_1\pm {ϵ}_2\pm \cdots \pm {ϵ}_L.$$

$$C=\left[\begin{array}{ccccc}\lambda & & & & \\ 1& \lambda & & & \\ & 1& \lambda & & \\ & & \ddots & \ddots & \\ & & & 1& \lambda \end{array}\right],C^{†}=\left[\begin{array}{ccccc}\lambda & 1& & & \\ & \lambda & 1& & \\ & & \ddots & \ddots & \\ & & & \lambda & 1\\ & & & & \lambda \end{array}\right].$$

$${ϵ}_j={\left(1+{\lambda }^2+2\lambda \cos k_j\right)}^{1/2}.$$

$$\sin (L+1)k=-{\lambda }^{-1}\sin Lk.$$

$$C{C}^{†}{\varphi }_k={ϵ}_k^2{\varphi }_k,$$

$$C^{†}C{\psi }_k={ϵ}_k^2{\psi }_k.$$

$${\varphi }_{k_j}(n)=A_{k_j}\sin (L+1-n)k_j,$$

$${\psi }_{k_j}(n)=-(-1{)}^j{A}_{k_j}\sin n{k}_j,$$

$$A_{k_j}^2=\frac{4}{2L+1-\csc k_j\sin (2L+1)k_j}.$$

$$G_{mn}=-\sum _k{\varphi }_k(m){\psi }_k(n),$$

$$⟨{\Psi }_0|{\sigma }_1^x{\sigma }_L^x|{\Psi }_0⟩$$

$$⟨{\Psi }_0|{\sigma }_1^y{\sigma }_L^y|{\Psi }_0⟩$$

$$⟨{\Psi }_0|{\sigma }_1^z{\sigma }_L^z|{\Psi }_0⟩$$

$${\rho }_a=\underset{L\to \infty }{\lim }⟨{\Psi }_0|{\sigma }_1^a{\sigma }_L^a|{\Psi }_0⟩,$$

$${\rho }_z=\{\begin{array}{l}1-{\lambda }^2,\lambda <1\\ 0,\lambda \ge 1.\end{array}$$

$${\rho }_z=\{\begin{array}{l}(1-{\lambda }^2{)}^{1/4},\lambda <1\\ 0,\lambda \ge 1,\end{array}$$

$$k_j=\frac{2j\pi }{2L+1},j=1,\dots ,L.$$

$$k_L=\pi +\mathrm{i}v,$$

$$\sinh (L+1)v={\lambda }^{-1}\sinh Lv.$$

$$k_j\approx \frac{\pi j}{L}-\frac{\lambda }{1+\lambda }\frac{\pi j}{L^2},j=1,\dots ,L.$$

$$L{k}_j=\pi j-\pi {\kappa }_j+O\left(\frac{1}{L}\right),j=1,\dots ,L$$

$$\cot (\pi {\kappa }_j)=\frac{{\lambda }^{-1}+\cos (\pi j/L)}{\sin (\pi j/L)}$$

$$\pi {\kappa }_j=\frac{\pi j}{2L}+{\tan }^{-1}\left[\frac{\lambda -1}{\lambda +1}\tan \left(\frac{\pi j}{2L}\right)\right].$$

$$e_0(\lambda )=-F\left(-\frac{1}{2},-\frac{1}{2};1;{\lambda }^2\right).$$

$$G_{1L}$$

$$=$$

Full text

00footnotetext: Dedicated to the memory of Vladimir Rittenberg

1 Introduction

The Ising model is one of the most widely studied models in physics. The understanding obtained from the exact solution of the two-dimensional Ising model and its one-dimensional quantum counterpart played a pivotal role in the formulation of the theory of phase transitions and critical phenomena [1, 2, 3, 4]. Although describing interacting classical and quantum spin systems, a feature of the Ising model is that the exact solution is described in terms of free fermions. It is this property of free fermions which allows the relatively straightforward calculation of physical properties, given the powerful framework of statistical mechanics.

In terms of the Pauli matrices $\sigma^{x}_{j}$ and $\sigma^{z}_{j}$ acting at site $j$ on a chain of length $L$ the hamiltonian of the one-dimensional quantum Ising model in a transverse field can be defined by

[TABLE]

The field variable $\lambda$ plays a temperature-like role, with a quantum phase transition at $\lambda=\lambda_{c}=1$ between ferromagnetic ($\lambda<1$) and paramagnetic ($\lambda>1$) phases [4]. Here we have imposed open boundary conditions, for reasons which will become apparent.

The quantum Ising chain belongs to several families of exactly solved $N$-state quantum models. Such a more general model of interest here is the $Z(N)$ model defined by the hamiltonian

[TABLE]

where $\sigma_{j}$ and $\tau_{j}$ are $Z(N)$ operators defined in matrix form by

[TABLE]

Here the identity $I$, $\sigma$ and $\tau$ are $N\times N$ matrices, with $\sigma$ and $\tau$ in position $j$. The clock and shift matrices $\sigma$ and $\tau$ have components

[TABLE]

with $\omega={\mathrm{e}}^{2\pi{\mathrm{i}}/N}$ and $\tau_{m,N}=\delta_{m,1}$. For $N=2$ these are the Pauli matrices $\sigma^{z}$ and $\sigma^{x}$, with hamiltonian (2) reducing to hamiltonian (1). For general $N$ the clock and shift matrices satisfy the relations

[TABLE]

with $\sigma^{N}=\tau^{N}=I$.

The $Z(N)$ model defined by hamiltonian (2) was introduced by Baxter [5, 6, 7]. It follows from a version of the $\tau_{2}$ model, also known as the Bazhanov-Stroganov model [8], a two-dimensional classical chiral spin model connecting the six-vertex model and the chiral Potts model in such a way that the chiral Potts model can be viewed as a descendent of the six-vertex model [8, 9].

It is important to note that the $Z(N)$ hamiltonian (2) is non-Hermitian for $N\geq 3$, with a complex eigenspectrum. With the addition of the Hermitian conjugate term this model has been studied recently in the context of parafermionic edge modes [10, 11], mostly for $N=3$. The non-Hermitian $Z(N)$ hamiltonian (2) has been solved exactly only for open boundary conditions. The energy eigenspectrum displays the remarkable property of free parafermions [12]. This free parafermionic structure has also been confirmed in the related $\tau_{2}$ model with open boundaries [13, 14, 15]. The $Z(N)$ hamiltonian (2) and the related $\tau_{2}$ model are the only known models with an entire spectrum described precisely by free parafermions.

Relatively little is known so far about the physical properties of free parafermions, compared to free fermions, for which a wealth of information is known. The free parafermion $Z(N)$ model (2) has been shown [16, 17] to share some critical exponents with the other parafermionic $Z(N)$ models. For example, the specific heat and anisotropic correlation length exponents, $\alpha=1-2/N$, $\nu_{\parallel}=1$ and $\nu_{\perp}=2/N$, are the same as for the superintegrable $N$-state chiral Potts model [6, 18, 9]. This specific heat exponent is also shared by the $Z(N)$ Fateev-Zamolodchikov model [22]. We note also that the free parafermion $Z(N)$ model provides what appears to be the first example of a one-dimensional hamiltonian exhibiting boundary dependent bulk behavior. For example, the bulk ground-state energy per site and various critical exponents differ between open and periodic boundary conditions [17].

Our specific interest here is in the nature of spin correlations. Spin correlations and related order parameters for the Ising model, in both the classical and quantum formulations, are well understood (see, e.g., [2]). Although much work remains to be done for the generalization of the Ising results to the $N$-state chiral Potts model, there is a remarkable result for the order parameters $M_{k}$, $k=1,2,\ldots,N-1$. The general result,

[TABLE]

reduces to the spontaneous magnetization of the Ising model for $N=2$ and $k=1$. The formula (7) was first conjectured by Albertini, McCoy, Perk and Tang in 1989 [18], then finally proven by Baxter in 2005 [19]. This result is the culmination of a very fruitful conjectural path beginning with Howes, Kadanoff and den Nijs [20] for $N=3$ and the discovery of the $N$-state superintegrable chiral Potts chain by von Gehlen and Rittenberg [21].

Here we investigate some ground-state expectation values for the free parafermion $Z(N)$ model. We focus on open boundary conditions for which the $Z(N)$ free parafermion model is exactly solved and for which exact information is readily obtained for the Ising case. As a special case, the $N=2$ results are thus of relevance to the study of the more general free parafermion model. In Section 2 we summarise the known solution for the eigenvalues and eigenvectors of the Ising hamiltonian (1) with open boundary conditions. These results are used in Section 3 to consider various exact spin correlations for the open quantum Ising chain. Exact results are derived for some ground-state expectation values for the open $Z(N)$ chain in Section 4. These results are tested numerically in Section 5 for $N=3$. A concluding discussion of the results is given in Section 6.

2 Eigenvectors and spin correlations for the open Ising chain

The hamiltonian (1) was solved long ago by Pfeuty [23], closely following the similar solution of the XY model with open boundary conditions [24]. The open quantum Ising chain has been discussed by a number of authors in different contexts, including the dynamics (see, e.g., Refs [25, 26, 27] and references therein). Early work was also done on this model in relation to the connection between finite-size corrections and conformal field theory [28]. A pedagogic description of the solution of both models in terms of free fermions has been given by Karevski [25]. The key ingredient is the Jordan-Wigner transformation, as applied in the seminal paper by Lieb, Schultz and Mattis [24] outlining the solution of the XY model in terms of free fermions. Indeed the XY and quantum Ising models are closely related and can be treated under the same banner in the form of the XYh model (see, e.g., Refs [25, 27] for such a treatment for open boundary conditions). In this Section we summarise some results of relevance to the calculation of the spin correlations for the open quantum Ising chain defined by hamiltonian (1).

The energy eigenspectrum of hamiltonian (1) has the simple free fermion form

[TABLE]

This covers all $2^{L}$ eigenvalues in the energy spectrum for a chain of length $L$. The quasi energy levels $\epsilon_{j}$ appearing in (8) are functions of $\lambda$. They are determined by the eigenvalues of the $L\times L$ matrices $C^{\dagger}C$ or $CC^{\dagger}$, where

[TABLE]

These eigenvalues are given by

[TABLE]

The roots $k_{j}$, $j=1,\ldots,L$, satisfy the equation

[TABLE]

The relevant eigenvectors are given in terms of the matrices $C$ and $C^{\dagger}$ by $C\psi_{k}=\epsilon_{k}\phi_{k}$ and $C^{\dagger}\phi_{k}=\epsilon_{k}\psi_{k}$. Thus

[TABLE]

The eigenvector components are given explicitly by

[TABLE]

for $n=1,\ldots,L$. The normalization factor $A_{k_{j}}$ is

[TABLE]

Given the eigenvector components, a key ingredient in the calculation of spin correlations is the quantity

[TABLE]

in terms of which, for example, the end-to-end spin correlations are determined by [24]

[TABLE]

Here $\Psi_{0}$ is the ground-state eigenvector of hamiltonian (1). It is convenient to define the limits

[TABLE]

with $a=x,y,z$.

In this way, Pfeuty [23] obtained the result

[TABLE]

As noted at the time, this differs from the result

[TABLE]

obtained for periodic boundary conditions [29, 30, 23]. To obtain this result one has to consider the correlation $\langle\sigma_{i}^{z}\sigma_{i+R}^{z}\rangle$ in the limit $R\to\infty$. The difference between the results (21) and (22) is attributed to an end-effect [23]. The result (22) leads directly to the bulk magnetization $M_{z}$, with $\rho_{z}=M_{z}^{2}$ and thus $M_{z}=(1-\lambda^{2})^{1/8}$ for $\lambda<1$, with $M_{z}=0$ for $\lambda\geq 1$.

3 Exact results for $N=2$

Before calculating various spin correlations, we need to consider the solutions of equation (11), which differ according to the physical regime described by $\lambda$. For $\lambda=\lambda_{c}=1$ the $L$ roots are given simply by

[TABLE]

In contrast to the periodic case, more work is required for other values of $\lambda$. There are $L$ real roots in the interval $(0,\pi)$ for $\lambda\geq 1+1/L$, while for $\lambda<1+1/L$ there are $L-1$ roots in $(0,\pi)$, the remaining root being complex, of the form

[TABLE]

where $v$ satisfies the equation

[TABLE]

In particular, it follows that $\mathrm{e}^{v}\approx 1/\lambda$ for large $L$, which features in the calculations below.

On the other hand, for $\lambda>1$ the roots can be approximated for large $L$ by

[TABLE]

This result follows by writing

[TABLE]

which gives

[TABLE]

with solution

[TABLE]

This expansion mirrors that for the XY spin chain with open boundary conditions [24]. Here we have taken the $N=2$ value of the corresponding results for the free parafermion chain [16].

The ground-state energy per site in the limit $L\to\infty$ is given simply in terms of the hypergeometric function $F(a,b;c;z)=\,_{2}\!\,F_{1}(a,b;c;z)$ by

[TABLE]

This result is valid for $\lambda\leq 1$, with $e_{0}(\lambda)=\lambda\,e_{0}(1/\lambda)$ for $\lambda>1$.

3.1 End-to-end correlations

First consider the calculation of the correlations (17), (18) and (19) at $\lambda=1$, for which the relevant sums involving the eigenvector components with the exact roots (23) can be evaluated in closed form. We have

[TABLE]

[TABLE]

A similar result for the sum

[TABLE]

is too cumbersome to reproduce here. Nevertheless, it can be readily established from the result that

[TABLE]

As can be seen from the closed form expressions (32) and (34) both $G_{1L}$ and $G_{L1}$ vanish in this limit. Putting these results together gives

[TABLE]

with $\rho_{z}=0$ and $\rho_{y}=0$.

Turning to $\lambda<1$, the sum $G_{1L}$ is dominated with increasing $L$ by the complex root $k_{L}$, with

[TABLE]

from which Pfeuty’s result (21) for $\rho_{z}$ follows. On the other hand, $G_{L1}$ vanishes for large $L$, and so $\rho_{y}=0$ in this regime. It follows that to calculate $\rho_{x}$ from (17), one needs to determine the value of $G_{11}=G_{LL}$. In the limit $L\to\infty$

[TABLE]

giving the result

[TABLE]

This result reduces to (37) at $\lambda=1$.

Unlike the end-to-end correlation $\rho_{z}$ given in (21), the end-to-end correlation $\rho_{x}$ does not vanish for $\lambda>1$. Rather the result for $\lambda>1$ satisfies the duality relation

[TABLE]

with $0\leq\rho_{x}(\lambda)\leq 1$.

3.2 Nearest-neighbor correlations

When periodic boundary conditions are applied, the end-to-end correlations are equivalent to nearest-neighbor correlations. For periodic boundary conditions, the results for various Ising correlations can be obtained from the key results given in McCoy’s book [2]. Using different notation to distinguish the periodic case, the two nearest-neighbour correlations in the limit $L\to\infty$ are

[TABLE]

[TABLE]

At $\lambda=1$ these results reduce to the values

[TABLE]

A related result is the single site magnetization $M_{x}=\langle\sigma_{1}^{x}\rangle$, given by [2]

[TABLE]

with $M_{x}=\frac{2}{\pi}$ at $\lambda=1$. For systems with periodic boundary conditions the correlations $\langle\sigma_{n}^{z}\sigma_{n+1}^{z}\rangle$, $\langle\sigma_{n}^{x}\sigma_{n+1}^{x}\rangle$ and $M_{x}=\langle\sigma_{n}^{x}\rangle$ are independent of site $n$. The Hellmann-Feynman theorem [2] can then be used to relate $\langle\sigma_{n}^{z}\sigma_{n+1}^{z}\rangle$ and $M_{x}$ to the ground-state energy $e_{0}(\lambda)$, with result

[TABLE]

The various relations among these results can be checked by applying the relation

[TABLE]

among other identities [2, 32].

On the other hand, no such translational invariance holds for systems with open boundary conditions. For the Ising model, the nearest-neighbor correlations are site-dependent. At the boundary,

[TABLE]

The limits of interest are

[TABLE]

with $a=x,y,z$. These correlations do not appear to have been calculated explicitly for open boundary conditions. They are straightforward to evaluate at $\lambda=1$. For this value the general nearest-neighbor correlations $\rho^{z}_{n,n+1}$ and $\rho^{x}_{n,n+1}$ are given by

[TABLE]

For $n=1$ these results simplify to $\rho^{z}_{12}=\frac{8}{5\pi}$ and $\rho^{x}_{12}=\frac{4096}{525\pi^{2}}$. The periodic values (44) are recovered in the limit $n\to\infty$. Physically, the periodic results are thus recovered for sites deep inside the open chain.

4 Exact results for general $N$

In this Section, in the absence of results for the eigenvectors of the $Z(N)$ free parafermion chain, we derive the ground-state expectation values which can be obtained by applying the Hellmann-Feynman theorem adapted for non-Hermitian systems (see, e.g., [31]). The quantities considered are defined for finite-size chains by

[TABLE]

For non-Hermitian systems, both the left $\langle\Phi|$ and right $|\Psi\rangle$ eigenvectors need to be considered. For example, the ground-state energy $E_{0}$ is given by

[TABLE]

The eigenvectors are not orthogonal, but rather form a biothonormal basis

[TABLE]

Moreover, for $N>2$ the eigenvectors are no longer real. Nevertheless, the expectation values (55) and (56) are observed to be real. As is the ground-state energy $E_{0}$.

More specifically, the quantities of interest here are defined by

[TABLE]

In analogy with applying this approach for the Ising case with open boundary conditions, the values of $n$ are such that $1\ll n\ll L$, i.e., for sites $n$ far from the ends of the chain. We apply the Hellmann-Feynman theorem under this condition, and assuming that the contributions from sites near the boundaries become negligible compared to the contributions from the bulk. The expected results are thus

[TABLE]

We can now make use of the exact result for the ground-state energy per site, given in the $L\to\infty$ limit by [16]

[TABLE]

Like the $N=2$ result (30), this result is valid for $\lambda\leq 1$, with $e_{0}(\lambda)=\lambda\,e_{0}(1/\lambda)$ for $\lambda>1$.

Subject to the above proviso for the values of $n$, and making use of the relation (48), the results for general $N$ follow in terms of hypergeometric functions as

[TABLE]

and

[TABLE]

Here we have made use of the identity [32]

[TABLE]

As to be expected, when $N=2$ the above results reduce to the Ising results (45) and (42) given in Section 3.2. At $\lambda=1$ they simplify to the values

[TABLE]

In general $\rho^{X}(\lambda)=\rho^{Z}(1/\lambda)$.

In the next Section we test the results (64) and (65) numerically.

5 Numerical tests

In testing the above results numerically the left eigenvectors appearing in the expectation values (55) and (56) can be calculated from $H_{\mathrm{pf}}^{\dagger}|\Phi_{0}\rangle=E_{0}|\Phi_{0}\rangle$. For this particular non-Hermitian model, this result is equivalent to that of a spatial symmetry observed by Baxter [5], see also [12]. The numerical results are obtained using DMRG implemented by the ITensor library [33]. Some modifications were required to support non-Hermitian systems; primarily, the Davidson eigensolver was replaced with Arnoldi. This results in a considerable slowdown, so the chain lengths typical for Hermitian systems could not be achieved. For this model we can take advantage of the ground-state energy being known exactly for finite size $L$, from which the difference with the numerical energy determined using DMRG was used as the convergence parameter. The numerical results presented here were converged to a tolerance of $10^{-6}$ relative to the exact result for $E_{0}$.

Specifically, the ground-state energy of the $N^{L}\times N^{L}$ matrix $H_{\mathrm{pf}}$ is given by

[TABLE]

In terms of the coupling parameter $g=1/\lambda^{N/2}$, the free parafermion quasi-energies $\epsilon_{j}$ appearing in (68) are the eigenvalues of the $L\times L$ matrix [17]

[TABLE]

As already mentioned, we directly diagonalise this smaller matrix numerically to test the convergence. Alternatively, the value $E_{0}$ used in the convergence test may be obtained by numerically solving an equation for the quasi-momenta values determining $\epsilon_{j}$ [16]. In the numerical computations it is convenient to work in terms of the coupling parameter

[TABLE]

The critical point $\lambda_{c}=1$ thus corresponds to $\gamma_{c}=\frac{1}{2}$, with the range of coupling parameters now restricted to $0\leq\gamma\leq 1$.

Using DMRG as described above, the eigenvectors $\langle\Phi_{0}|$ and $|\Psi_{0}\rangle$ of $H_{\mathrm{pf}}$ are constructed for finite $L$ to evaluate the ground-state expectation values (55) and (56). These quantities are shown in figures 1 and 2 for $N=3$ and $L=60$ with site $n$ varying over the length of the chain. Very similar behaviour is observed for the $N=2$ Ising case. In particular, away from the ends of the chain, the curves flatten, indicative of a ‘bulk’ limiting result. For the Ising case, the flattened curves agree with the values obtained for periodic boundary conditions. For comparison, figures 1 and 2 also show the corresponding values obtained for periodic boundary conditions. It is clear that the expectation values depend on the boundary conditions for $N=3$, consistent with the behaviour observed for the ground state energy [17]. However, sufficiently away from criticality it appears that these values are no longer dependent on the boundary conditions. Convergence to the exact results (64) and (65) for $N=3$ can be seen in figures 3 and 4, which plot the ground-state expectation values at the centre of the open chain for increasing chain size $L$. Similar convergence is observed for the $N=2$ Ising case for same size chains.

6 Discussion

Motivated by the more general free parafermion $Z(N)$ model (2), we have studied elementary correlations for the quantum Ising model (1) subject to open boundary conditions. The results for the Ising model have been given in Section 3. Our aim has been to generalise these results to the $Z(N)$ model. We have been able to do this in Section 4 for the particular ground-state expectation values obtained by applying the Hellmann-Feynman theorem. The key results are given in equations (64) and (65). We tested these results numerically in Section 5 for $N=3$ and believe them to be correct for general $N$.

As tempting as it may be, we have not referred to these results as correlations. This is because in general non-Hermitian systems are non-unitary, leading to the problem of formulating a meaningful quantum mechanical description with physical correlations. The well known exception is the class of non-Hermitian systems which are ${\cal P}{\cal T}$ symmetric, with a real eigenspectrum for which a definite metric is guaranteed [34]. Here one can also define quasi-Hermitian and pseudo-Hermitian systems. For both types the eigenvalues are real, but no definite conclusions can be made with regard to the existence of a definite metric for the latter.

Of further interest are the Ising model end-to-end correlations (21) and (40). The simple result (21) was obtained long ago by Pfeuty [23]. The corresponding end-to-end correlation $\rho_{x}$ given by (40) does not appear to have been discussed before. Rather than generalising in a relatively simple way for the $Z(N)$ model, like for the order parameters (7) or the results (64) and (65) presented here, such end-to-end ‘correlations’ are seen to diverge with increasing system size at the critical value $\lambda_{c}=1$. In contrast to Hermitian systems, such divergent behaviour should not be surprising for the physics of non-Hermitian systems. An example from optics is the behaviour of the Petermann factor, which depends on the overlap between the left and right eigenvectors of a non-unitary wave operator. This factor takes very large values at resonances which can be explained by degeneracies in the eigenspectrum [35]. Another example relates to optical singularities in non-Hermitian chiral crystals [36]. This divergent behaviour will be discussed in detail elsewhere.

We remain optimistic that it should be possible to derive exact results for various ‘correlations’ of the $Z(N)$ model described by free parafermions, including the results obtained here via the Hellmann-Feynman theorem. Recalling the summary of the free fermion formalism given in Section 2, it is possible to obtain analogous eigenvector components $\psi$ and $\phi$ for the $Z(N)$ model via the corresponding $L\times L$ matrices $C$ and $C^{\dagger}$ [17]. However, there is a conceptual difficulty in relating them back to the spin ‘correlations’ of interest for the $Z(N)$ model. This is related to the problem of establishing (or knowing) an analog of Wick expansion for parafermions. The Wick expansion is the key ingredient in the calculation of correlations in the free fermion XY and Ising chains [24, 25]. Whether or not such a crucial step can be made for free parafermions, allowing the exact calculation of spin ‘correlations’ and dynamics, remains to be seen. There is a glimmer of hope, however. This relates to the observation that the $Z(N)$ model in question satisfies the property of reflection positivity [37]. Reflection positivity has been noted as a perturbative form of Wick rotation.

In concluding, we remark that after all these years it remains a major challenge to calculate the eigenvectors of the $N$-state chiral Potts model. Some development was made in calculating the eigenvectors for the superintegrable case [38, 39, 40, 41], but not sufficiently to allow the explicit calculation of correlation functions, and for example the calculation of the order parameters (7) via this approach. In comparison it would seem that the exact calculation of the eigenvectors of the free parafermion $Z(N)$ model may provide an easier challenge. Indeed, we note that some results have been reported for the eigenvectors of this model [42], but the results are far from transparent. We hope that the progress reported here will inspire further work on this model.

We are grateful to Francisco Alcaraz, Rodney Baxter, Michael Berry, Fabian Essler and Paul Fendley for helpful comments at various stages of this work. This project has been supported by the Australian Research Council Discovery Project DP180101040 and the National Natural Science Foundation of China Grant No. 11575037.