text:
# Simple model for a quantum wire II. Correlations

## Introduction

In a previous paper QWI , hereafter referred to as I, the authors have developed a simple model describing the main features shown by a Quantum Wire, namely: Band Structure when the structure is periodically arranged, Anderson Localization in the random case and easy calculation of Lyapunov exponents, density of electronic states and electrical conductance endowed with an encouraging fractal behaviour. We address the interested reader to reference I for the details of the model and the intrincacies of the computer calculations.

The main drawback of the scheme thereby presented was the lack of a coherent description of the possible correlations which must be present in any realistic model of one dimensional conductance. It is the purpose of the present paper to fill this gap by presenting a surprisingly simple way to incorporate correlations entirely based upon the random nature of the model. This means that in addition to describe correlations with different parameters of the wire (spacing, number and amount of species, … ) the model is flexible enough to be suitable for applying it to different correlation schemes such as the Tight-Binding model, the Anderson model and likely many others correlation prescriptions with sensible physical content.

The paper is organized as follows. In Section 1 we analyze the meaning of the correlations in terms of probabilities of a given atom to be followed by other of a different character. After this analysis we are able to set some limits for the values of these probabilities. The next step is to incorporate those probabilities to the Functional Equation developed in I. How this functional equation is modified by the introduction of correlations in the way just described and how it can be used to calculate the correlation-corrected Lyapunov exponent is the subject of Section 2. In Section 3 some results are shown along with a discussion about the effects of the correlations to close finally with a section of Conclusions.

## I Correlations

Let us consider disordered quantum wires with short-range correlations which will be specified by the probability for a certain atomic species to be followed by another type of atom in the chain sequence. In this model, the system will be characterized by the species concentrations $`\{c_i\}_{i=1,\mathrm{},m}`$ where $`m`$ is the number of different species, and the set $`\{p_{ij}\}_{i,j=1,\mathrm{},m}`$ where $`p_{ij}`$ means the probability for $`i`$ to be followed by $`j`$ (the prob. of finding $`j`$ right after $`i`$). The correlations introduced in this way must be the consequence of the existence of an atomic interaction (or the effect of some physical parameters such as the size of the atoms) which might choose certain spatial sequences of the atoms modifying subtly the otherwise chain’s purely random character. Thus this procedure seems a natural manner to account for correlations which can be present in nature or even those that can be produced inside a manufactured disordered quantum wire in a nonintentional way aside from the fact that of course one can always try to construct a wire exhibiting the desired correlations to recover the theoretical results.

Let us clarify the meaning of $`p_{ij}`$. It is the ratio of the number of -$`ij`$\- clusters and the number of $`i`$ atoms ($`N_i`$). The probability $`p_{ij}`$ is not the probability for $`j`$ to be preceded by $`i`$, because the latter would be obtained dividing by the number of $`j`$ atoms ($`N_j`$) instead of $`N_i`$. In the limit $`N\mathrm{}`$, the system will be determined by $`\{c_i\}`$ and $`\{p_{ij}\}`$ and therefore the density of states will be given by a limiting distribution depending only on those quantities. And that distribution is of course the same whether we move from left to right or from right to left along the chain. So, the correlations we measure have to coincide in both directions. Then it is clear that for the infinite chain $`p_{ij}`$ means the probability for an ’$`i`$’ species atom to be followed or preceded by a ’$`j`$’ species atom.

Once we have understood the meaning of the correlations, let us derive the relations among them. What is the probability of finding at any position of the chain the cluster -iji-? We can write this quantity as,

$`c_ip_{ij}p_{ji}`$ prob. of finding $`i`$ finding $`j`$ after $`i`$
$``$ finding $`i`$ after $`j`$
$`c_jp_{ji}^2`$ prob. of finding $`j`$ finding $`i`$ after $`j`$
$`\text{ finding }i\text{ before }j.`$
Thus the equations for the correlations are,

$`c_ip_{ij}=c_jp_{ji}`$ (1a)
$`{\displaystyle \underset{j=1}{\overset{m}{}}}p_{ij}=1`$ (1b)
$`0p_{ij}1`$ $`i,j=1,\mathrm{},m.`$ (1c)
According to the equations, the correlations matrix $`(p_{ij})`$ is completely known from the above diagonal elements, $`\frac{m(m1)}{2}`$ where $`m`$ is the number of different species. However these elements are not completely independent due to (1c), because when one of the correlations is chosen, the maximum allowed value for some of the rest may be affected. This fact can be clearly seen in two simple examples:

the matrix is completely determined by $`p_{12}`$,

$$p_{11}=1p_{12};p_{21}=\frac{c_1}{c_2}p_{12};p_{22}=1p_{21}$$

but it must be $`p_{21}1p_{12}{\displaystyle \frac{c_2}{c_1}}`$.

Therefore $`p_{12}min\{1,{\displaystyle \frac{c_2}{c_1}}\}`$.

the matrix is determined by $`p_{12},p_{13},p_{23}`$. And we choose their values in that order.

$`p_{12}min\{1,{\displaystyle \frac{c_2}{c_1}}\}`$ which implies

$`p_{13}min\{1p_{12},{\displaystyle \frac{c_3}{c_1}}\}`$.

On the other hand it must be,

$$p_{23}1p_{21}$$

$$p_{32}1p_{31}p_{23}\frac{c_3}{c_2}(1p_{31}),$$

therefore $`p_{23}min\{1{\displaystyle \frac{c_1}{c_2}}p_{12},{\displaystyle \frac{c_3}{c_2}}{\displaystyle \frac{c_1}{c_2}}p_{13}\}`$.

An expression for the general form of the maximum values can be obtained for arbitrary $`m`$. Choosing the above diagonal elements of the correlations matrix by rows,

$$𝐏=\left(\begin{array}{cccccc}p_{11}& \textcolor[rgb]{0.4,0.4,0.4}{p}_{\textcolor[rgb]{0.4,0.4,0.4}{12}}& \textcolor[rgb]{0.4,0.4,0.4}{}& \textcolor[rgb]{0.4,0.4,0.4}{\mathrm{}}& \textcolor[rgb]{0.4,0.4,0.4}{}& \textcolor[rgb]{0.4,0.4,0.4}{p}_{\textcolor[rgb]{0.4,0.4,0.4}{1}\textcolor[rgb]{0.4,0.4,0.4}{m}}\\ p_{21}& p_{22}& \textcolor[rgb]{0.4,0.4,0.4}{p}_{\textcolor[rgb]{0.4,0.4,0.4}{23}}& \textcolor[rgb]{0.4,0.4,0.4}{}& \textcolor[rgb]{0.4,0.4,0.4}{}& \textcolor[rgb]{0.4,0.4,0.4}{p}_{\textcolor[rgb]{0.4,0.4,0.4}{2}\textcolor[rgb]{0.4,0.4,0.4}{m}}\\ \mathrm{}& \mathrm{}& \mathrm{}& \textcolor[rgb]{0.4,0.4,0.4}{\mathrm{}}& & \textcolor[rgb]{0.4,0.4,0.4}{\mathrm{}}\\ \mathrm{}& \mathrm{}& & \mathrm{}& \textcolor[rgb]{0.4,0.4,0.4}{\mathrm{}}& \textcolor[rgb]{0.4,0.4,0.4}{\mathrm{}}\\ p_{(m1)1}& p_{(m1)2}& 2& p_{(m1)(m1)}& \textcolor[rgb]{0.4,0.4,0.4}{p}_{\textcolor[rgb]{0.4,0.4,0.4}{(}\textcolor[rgb]{0.4,0.4,0.4}{m}\textcolor[rgb]{0.4,0.4,0.4}{}\textcolor[rgb]{0.4,0.4,0.4}{1}\textcolor[rgb]{0.4,0.4,0.4}{)}\textcolor[rgb]{0.4,0.4,0.4}{m}}\\ p_{m1}& p_{m2}& 2& p_{m(m1)}& p_{mm}\end{array}\right)$$

then

$$p_{ij}min\{1\frac{1}{c_i}\underset{k=1}{\overset{i1}{}}c_kp_{ki}\underset{k=i+1}{\overset{j1}{}}p_{ik};\frac{c_j}{c_i}\frac{1}{c_i}\underset{k=1}{\overset{i1}{}}c_kp_{kj}\}$$
(2)
where the limits are given by the previous chosen correlations.

Once the concentrations and $`𝐏`$ are fixed, the characterization of the quantum wire can be given in a compact form using the matrix

$$𝐐=(q_{ij})(c_ip_{ij}),$$
(3)
which is symmetric ($`𝐐=𝐐^t`$) and it determines uniquely the wire and its correlations:

$$\underset{j=1}{\overset{m}{}}q_{ij}=c_i;\underset{i=1}{\overset{m}{}}q_{ij}=c_j;p_{ij}=\frac{q_{ij}}{c_i};$$

$`q_{ij}`$ means the probability of finding the cluster -ij- (-ji-) at any position inside the wire.

This model naturally includes the situation in which the disorder in the wire is completely random111The case in which the wire is obtained “throwing” the atoms in the given concentrations randomly and there are no mechanisms (interactions) which select some privileged clusters or series of atoms., that is just defined by the values $`p_{ij}=c_ji=1,\mathrm{},m`$. Figure 1 shows an example of the correlation space for 2 species.

## II The functional equation

We will now obtain the correlation form of the functional equation, which enable us to calculate the density of states. In this process we are using some of the results and expressions obtained in reference I.

The relations for the phase transmission inside the wire read,

$$𝒯_j^1(\theta )=\mathrm{arctan}\left(2h_j(ϵ)\frac{1}{\mathrm{tan}\theta }\right)$$
(4a)
$$𝒯_j^1(\theta +n\pi )=𝒯_j^1(\theta )+n\pi \theta [0,\pi ).$$
(4b)
where $`h_j(ϵ)=\mathrm{cos}(ϵ)+\frac{(a/a_j)}{ϵ}\mathrm{sin}(ϵ)`$. To carry out the phase average we rearrange the terms depending on the species of the last atom we have passed,

$$\begin{array}{c}\mathrm{\Delta }\theta =c_1\frac{1}{N_1}\underset{j_1}{}\left\{𝒯_j(\theta _j)\theta _j\right\}+\mathrm{}\hfill \\ \hfill \mathrm{}+c_m\frac{1}{N_m}\underset{j_m}{}\left\{𝒯_j(\theta _j)\theta _j\right\}\end{array}$$
(5)
where $`N_i`$ is the number of atoms of species $`i`$ and the sumatory over $`j_i`$ means we are summing over all the positions $`j`$ such that the atom at $`(j1)`$ is of $`i`$-kind:

$$\mathrm{\Delta }\theta =c_1\mathrm{\Delta }\theta _1+\mathrm{}+c_m\mathrm{\Delta }\theta _m=\underset{\gamma }{}c_\gamma \mathrm{\Delta }\theta _\gamma .$$
(6)
From now on the greek indices will refer to species while the latin indices mean the sites of the chain. In this case there exist correlations in the sequence of atoms, that is the probability for one species to be a neighbour of another type can be different for all of them. This implies that each species will give rise to a particular distribution function for the phase.

$$\mathrm{\Delta }\theta _\gamma =\frac{1}{N_\gamma }\underset{j_\gamma }{}_0^\pi \frac{dW_j^\gamma (\theta )}{d\theta }\left\{𝒯_j(\theta )\theta \right\}𝑑\theta $$
(7)
where $`W_j^\gamma (\theta )`$ is the distribution function of the phase at the position $`j`$ generated by the species $`\gamma `$. The equations for the distributions ((B.14) in I) become here,

$$W_j^\gamma (\theta )=W_{j1}^{()}\left(𝒯_\gamma ^1(\theta )\right)W_{j1}^{()}\left(\frac{\pi }{2}\right)+1$$
(8a)
$$W_j^\gamma (\theta +r\pi )=W_j^\gamma (\theta )+r\theta [0,\pi )r$$
(8b)
for all species and positions, where the superscript “$`()`$” means that we do not know the species of the atom preceding the $`\gamma `$ atom, placed at position $`(j1)`$; we only know the probability for a certain species to be there. Let us introduce now the average over all possible sequences with the given concentrations.

$$\mathrm{\Delta }\theta _\gamma =\frac{1}{N_\gamma }\underset{j_\gamma }{}_0^\pi \frac{dW_j^\gamma (\theta )}{d\theta }\underset{\beta }{}p_{\gamma \beta }\left\{𝒯_\beta (\theta )\theta \right\}d\theta $$
(9)
and,

$$W_j^\gamma (\theta )=\underset{\beta }{}p_{\gamma \beta }\left\{W_{j1}^\beta \left(𝒯_\gamma ^1(\theta )\right)W_{j1}^\beta \left(\frac{\pi }{2}\right)\right\}+1.$$
(10)
The next step consists in taking the limit $`N\mathrm{}`$ to approach the limiting distribution.

We define $`𝐖_\gamma (\theta )=\underset{N_\gamma \mathrm{}}{lim}{\displaystyle \frac{1}{N_\gamma }}{\displaystyle \underset{j_\gamma }{}}W_j^\gamma (\theta )`$ which is the solution of:

$$𝐖_\gamma (\theta )=\underset{\beta }{}p_{\gamma \beta }\left\{𝐖_\beta \left(𝒯_\gamma ^1(\theta )\right)𝐖_\beta \left(\frac{\pi }{2}\right)\right\}+1$$
(11a)
$$𝐖_\gamma (\theta +r\pi )=𝐖_\gamma (\theta )+r\theta [0,\pi )r.$$
(11b)
And finally we have to calculate,

$$\mathrm{\Delta }\theta =\underset{\gamma ,\beta }{}c_\gamma p_{\gamma \beta }_0^\pi \frac{d𝐖_\gamma (\theta )}{d\theta }\left\{𝒯_\beta (\theta )\theta \right\}𝑑\theta .$$
(12)
To carry out the integral we shall make use of the existence of a value $`\theta _0`$ such that $`𝒯_\gamma (\theta _0)=\theta _1\gamma `$, and proceeding by parts,

$$\begin{array}{c}\mathrm{\Delta }\theta =(\theta _1\theta _0)\stackrel{\text{}}{\stackrel{}{\underset{\gamma ,\beta }{}c_\gamma p_{\gamma \beta }_{\theta _0}^{\theta _0+\pi }𝐖_\gamma (\theta )\frac{d𝒯_\beta (\theta )}{d\theta }𝑑\theta }}\hfill \\ \hfill +\underset{\gamma }{}c_\gamma _{\theta _0}^{\theta _0+\pi }𝐖_\gamma (\theta )𝑑\theta .\end{array}$$
(13)
To evaluate the middle term we use the average of equations (11a) which is

$$\underset{\gamma }{}c_\gamma 𝐖_\gamma (\theta )+\underset{\beta }{}c_\beta 𝐖_\beta \left(\frac{\pi }{2}\right)1=\underset{\gamma ,\beta }{}c_\gamma p_{\gamma \beta }𝐖_\beta \left(𝒯_\gamma ^1(\theta )\right).$$
(14)
With a simple change of variable, the use of (14) and the equations for the correlations we can write,

$$\begin{array}{c}\text{}=\underset{\gamma ,\beta }{}c_\beta p_{\beta \gamma }_{\theta _1}^{\theta _1+\pi }𝐖_\gamma \left(𝒯_\beta ^1(\theta )\right)𝑑\theta =\hfill \\ \hfill =\underset{\gamma }{}c_\gamma _{\theta _1}^{\theta _1+\pi }𝐖_\gamma (\theta )𝑑\theta +\pi \underset{\beta }{}c_\beta 𝐖_\beta \left(\frac{\pi }{2}\right)\pi ,\end{array}$$
(15)
and going back to (13),

$$\begin{array}{c}\mathrm{\Delta }\theta =(\theta _1\theta _0)+\pi \pi \underset{\beta }{}c_\beta 𝐖_\beta \left(\frac{\pi }{2}\right)+\hfill \\ \hfill +\underset{\gamma }{}c_\gamma \left\{_{\theta _0}^{\theta _0+\pi }𝐖_\gamma (\theta )𝑑\theta _{\theta _1}^{\theta _1+\pi }𝐖_\gamma (\theta )𝑑\theta \right\}.\end{array}$$
(16)
It is not hard to see with the help of (11b) that the last term is $`\theta _0\theta _1`$. Therefore,

$$\frac{\mathrm{\Delta }\theta (ϵ)}{\pi }=1\underset{\gamma }{}c_\gamma 𝐖_\gamma \left(\frac{\pi }{2}\right).$$
(17)
Thus the density of states reads ((B.8) in I),

$$g(ϵ)=\left|\underset{\gamma }{}c_\gamma \frac{d𝐖_\gamma \left(\frac{\pi }{2}\right)}{dϵ}\right|$$
(18)
which can be calculated numerically solving eqs. (11a) for several values of the energy.

### II.1 The Lyapunov exponent

The quantity that provides information about the degree of spatial localization of the electronic states is the Lyapunov exponent, which gives the exponential growth rate of the wave function, and thus its inverse is a measure of the localization length of a state exponentially localized inside the wire. The Lyapunov coefficient can be written for the infinite wire as

$$\mathrm{\Lambda }=\underset{N\mathrm{}}{lim}\frac{1}{N}\underset{j=0}{\overset{N}{}}\mathrm{log}\left|\frac{\mathrm{\Psi }_{j+1}}{\mathrm{\Psi }_j}\right|$$
(19)
where $`\mathrm{\Psi }_j`$ is the amplitude of the wave function at the $`j`$th atomic site. Using a mapping technique originally proposed in Ref.map,  we are able to write $`\mathrm{\Lambda }`$ as a function of the phase (Appendix A),

$$\mathrm{\Lambda }=\frac{1}{2}\underset{N\mathrm{}}{lim}\frac{1}{N}\underset{j=0}{\overset{N}{}}\mathrm{log}F_{\gamma _j}(\theta _j)$$
(20)
where $`F_{\gamma _j}(\theta )=12h_{\gamma _j}(ϵ)\mathrm{sin}(2\theta )+4h_{\gamma _j}^2(ϵ)\mathrm{cos}^2\theta `$ and $`\gamma _j`$ stands for the species of the $`j`$th atom. Due to the periodicity of $`F(\theta )`$ we can use the distribution functions of the phase to calculate the average (20):

$$\begin{array}{cc}\hfill \mathrm{\Lambda }& =\frac{1}{2}\underset{\gamma ,\beta }{}c_\gamma p_{\gamma \beta }_0^\pi \frac{d𝐖_\gamma (\theta )}{d\theta }\mathrm{log}F_\beta (\theta )𝑑\theta =\hfill \\ & =\frac{1}{2}\underset{\gamma }{}c_\gamma \mathrm{log}F_\gamma (\pi )\frac{1}{2}\underset{\gamma ,\beta }{}c_\gamma p_{\gamma \beta }_0^\pi 𝐖_\gamma (\theta )\frac{F_\beta ^{}(\theta )}{F_\beta (\theta )}𝑑\theta .\hfill \end{array}$$
(21)
From this last expression one can obtain numerically the Lyapunov exponent as a function of the energy.

## III Results

Once we have built the mathematical framework lying under the model, one can try to answer the physical questions that obviously arise from this new configuration of the disordered quantum wire: how strong is the effect of this short-range correlations?, does it produce a measurable change in the density of states or in the localization of the electrons? As in Ref. I, we solve the functional equation numerically to obtain the density of states and to calcultate the Lyapunov exponent. It seems really hard to find analytical solutions for equations (11) if they exist at all.

Let us analize in the first place the density of states (DOS) for our quantum wires. In Figure 2 the evolution of the DOS for a binary wire as a function of the correlation $`p_{12}`$ is shown for a certain value of the concentrations. This distribution of states is drastically changed from the initial situation in which the probability to find the cluster -12- is low to the final stage when we impose the atoms of species 1 to appear always isolated. Note however that the concentrations are the same in both cases!. Another example can be seen in Figure 3 for a wire composed of 4 species. As the number of species increases, the number of correlations grows as $`\frac{m(m1)}{2}`$, and the distribution of states can adopt a lot of different shapes: the exploration of the whole correlations and concentrations space for $`m>2`$ can take a long time. From these graphics we conclude that the correlations can unbalance the spectra of the disorder system quite far from the completely random configuration. In fact the correlation can be tuned to open (close) an energy gap or to increase (decrease) the number of available states in a certain energy interval without changing the concentrations of the atomic components. The variation of the DOS with the correlations is of course a function of the concentrations being the chain with homogenized participation the one whose distribution has the strongest dependence on $`𝐏`$. Let us also remark that the fractal behaviour of the density of states in certain energy ranges, reported in reference I, still manifests itself for the different correlations regimes.

Let us have a look at the localization of the electronic states. In Figures 4 and 5 the behaviour of the Lyapunov exponent for a binary chain changing the correlations and the concentrations is shown. In Figure 3 the different shapes of the degree of localization for different correlations can also be noticed for a 4-species wire. The effective influence of the correlations on the localization of the electrons inside the quantum wires is established. In the study of correlated disorder two main kinds of correlations have been taken into account namely short-range and long-range. The former, mainly represented by the random dimer model RD , has been widely accepted to be the reason for the emergence of resonant energies and a set of states close to the resonant one for which the localization length becomes larger that the system size kp-adame , improving greatly the transport properties bellani . The long-range correlations seem to be able to include in the spectrum mobility edges Moura , that can be entirely controlled by the correlators Izkro , yielding an energy interval of a complete transparency. In our type of wires we have not observed mobility edges nor resonant energies. However an important effect on the localization of the states is shown that seems to act globaly in the whole energy range in contrast to the more restricted effect of others short-range correlated disorder models. In fact for certain energies the Lyapunov exponent can be decreased an 80-90% of its maximun value changing the correlations at the expense of an increasing behaviour for other energies. Whether these correlations might cause $`\mathrm{\Lambda }`$ to go below the inverse of the length of a finite sample for some energy intervals is not certainly clear but it could be a point for future discussions. One can now argue if this influence of the correlations on the localization as well as on the DOS will be a measurable effect on short finite wires that represent a simple practical implementation, or on the contrary is just theoretically achievable only in infinite unreal systems. In order to answer that question and ascertain ourselves of the importance of the results we learn to generate finite sequences showing this kind of correlations.

In Figure 6 we show the Lyapunov exponent of an infinite wire for a certain concentration and correlation calculated using the functional equation, and the one obtained numerically from equation (19) for a correlated disordered chain of 1000 atoms. The convergence is excellent with only one realization of the correlated disorder using one thousand of atoms, and we have also checked it for others couplings of the species and other correlations obtaining a high degree of agreement. Thus it means that the described results are real effects with some practical importance on experimental realizations and on the other hand that the method of the functional equation, although designed for an infinite wire, can reproduce quite well the tendency in finite samples.

## IV Concluding remarks

To summarize, we have extended the model proposed in reference I and completed the formalism introduced there including the possible correlations which must be present in any realistic model of one dimensional conductance. The method described here is the most natural manner to account for correlations that can manifest in nature or even inside a manufactured quantum wire in a nonintentional way.

It has been established how due to the correlated disorder the density of states is essentially modified and the localization of the electronic states changes in a non negligible proportion which might alter some macroscopic properties of these structures.

In addition to the obtained results we are preparing a future report in which we will construct the functional equation formalism in a more coherent way relating it to others methods currently used. There it will be shown how this mathematical framework is flexible enough to be suitable for applying it to different disorder models such as the Tight-Binding model as well as other correlation schemes such as the random dimer.

###### Acknowledgements.

We acknowledge with thanks the support provided by the Research in Science and Technology Agency of the Spanish Goverment (DGICYT) under contract BFM2002-02609.

## Appendix A Lyapunov exponent

It can be shown that the electronic wave function satisfies,

$$2h_{\gamma _j}(ϵ)\mathrm{\Psi }_j=\mathrm{\Psi }_{j+1}+\mathrm{\Psi }_{j1}$$
(22)
where $`\mathrm{\Psi }_j`$ is the amplitude of the wave function at the $`j`$th atomic site of the chain. This equation can be visualized as a two dimensional mapping in polar coordinates with the identification $`\mathrm{\Psi }_{j+1}=r_{j+1}\mathrm{cos}\theta _{j+1}`$ and $`\mathrm{\Psi }_j=r_{j+1}\mathrm{sin}\theta _{j+1}`$,

$$\left(\begin{array}{c}r_{j+1}\mathrm{cos}\theta _{j+1}\\ r_{j+1}\mathrm{sin}\theta _{j+1}\end{array}\right)=\left(\begin{array}{cc}2h_{\gamma _j}(ϵ)& 1\\ 1& 0\end{array}\right)\left(\begin{array}{c}r_j\mathrm{cos}\theta _j\\ r_j\mathrm{sin}\theta _j\end{array}\right)$$
(23)
which yields the known transmission relationship for the phase and the following for the moduli $`r`$,

$$\frac{r_{j+1}^2}{r_j^2}F_{\gamma _j}(\theta _j)=12h_{\gamma _j}(ϵ)\mathrm{sin}(2\theta _j)+4h_{\gamma _j}^2(ϵ)\mathrm{cos}^2\theta _j.$$
(24)
If now we introduce the polar coordinates into the definition of the Lyapunov exponent we obtain:

$$\begin{array}{c}\mathrm{\Lambda }=\underset{N\mathrm{}}{lim}\frac{1}{N}\underset{j=0}{\overset{N}{}}\mathrm{log}\left|\frac{\mathrm{\Psi }_{j+1}}{\mathrm{\Psi }_j}\right|\hfill \\ \hfill =\underset{N\mathrm{}}{lim}\frac{1}{N}\underset{j=0}{\overset{N}{}}\mathrm{log}\left(\frac{r_{j+1}}{r_j}\right)+\underset{N\mathrm{}}{lim}\frac{1}{N}\underset{j=0}{\overset{N}{}}\mathrm{log}\left|\frac{\mathrm{cos}\theta _{j+1}}{\mathrm{cos}\theta _j}\right|\\ \hfill =\underset{N\mathrm{}}{lim}\frac{1}{N}\underset{j=0}{\overset{N}{}}\mathrm{log}\left(\frac{r_{j+1}}{r_j}\right)\end{array}$$
(25)
because according to the meaning of $`\mathrm{\Lambda }`$ the cosinus term must be zero. And using (24) we can calculate the Lyapunov coefficient as a function of the phase, as:

$$\mathrm{\Lambda }=\frac{1}{2}\underset{N\mathrm{}}{lim}\frac{1}{N}\underset{j=0}{\overset{N}{}}\mathrm{log}F_{\gamma _j}(\theta _j).$$
(26)
which is our expression (20).


format: text

text:
# MONTE CARLO STUDY OF FOUR-SPINON DYNAMIC STRUCTURE FUNCTION IN ANTIFERROMAGNETIC HEISENBERG MODEL

## I Introduction

The spin $`s=\frac{1}{2}`$ Heisenberg quantum spin chain describes the magnetic properties of quasi-one-dimensional antiferromagnetic compounds like KCuF3 hirakawa-kurogi . The spin dynamics is experimentally investigated using inelastic neutron scattering Exp . From a theoretical standpoint, the quantity of interest is the dynamic structure function (DSF) $`S`$ of two local spin operators. This is because the magnetic scattering cross section per magnetic site is directly proportional to it squires-lovesey .

The Heisenberg model has been studied quite intensively PSD ; Baxter , and because of the presence of the $`\mathrm{U}_q(\widehat{\mathrm{sl}}_2)`$ symmetry JimMi , a number of exact results are available. Static properties need only the Yang-Baxter relation Baxter , whereas dynamic correlation functions require the additional notion of vertex operators and exploit bosonization methods CorFun . One thus obtains compact expressions for form factors JimMi .

Regarding the dynamic structure function, the focus has so far been on $`S_2`$, the two-spinon contribution to the total $`S`$. First there has been the Anderson (semi-classical) spin-wave theory And , an approach based on an expansion in powers of $`1/s`$ and hence, exact only in the classical limit $`s=\mathrm{}`$. It can describe with some satisfaction compounds with higher spins TMMC , but fails in the quantum limit $`s=\frac{1}{2}`$. For this latter system, the Müller ansatz has been proposed, built mainly from finite-chain calculations and symmetry considerations Muller . It is an approximate expression for $`S_2`$ that accounts satisfactorily for many aspects of the phenomenology Exp . More recently, an exact expression for $`S_2`$ has been obtained, making extensive use of the $`\mathrm{U}_q(\widehat{\mathrm{sl}}_2)`$ symmetry BCK , and a comparison with the Müller ansatz shows that it gives a better account of the data KMB ; BKM .

Beyond the two-spinon DSF $`S_2`$ are of course all the other $`S_{2p>2}`$ contributions. The first one to look at is the four-spinon dynamic structure function $`S_4`$, and the purpose of the present article is to describe its behavior. An exact expression for $`S_4`$ has been derived in ABS , see also Boug . We intend to use that expression to describe the behavior of $`S_4`$ as a function of the neutron energy transfer $`\omega `$ and neutron momentum transfer $`k`$.

We must mention that a preliminary investigation into the behavior of $`S_4`$ has already been initiated in ABSS . There we have used quadratures to perform the integrations involved in the expression of $`S_4`$. But because of slow convergence of the algorithms we wrote, we could describe the behavior of $`S_4`$ only as a function of $`k`$, and only for relatively small values of $`\omega `$. In the present work, the integrations are performed using Monte Carlo methods, and this makes it possible not only to study $`S_4`$ as a function of $`k`$ for a wider range of values of $`\omega `$, but obtain a description of $`S_4`$ as a function of $`\omega `$ for a wide range of values of $`k`$ as well. It is however important to note that, for the same values of $`\omega `$, the behavior of $`S_4`$ as a function of $`k`$ we obtain here by Monte Carlo methods is similar and consistent with the one we obtained in ABSS  using quadratures.

Also, we systematically carry a comparison of each result we obtain regarding $`S_4`$ to a corresponding one regarding $`S_2`$. The reason is that there is good familiarity in the literature with the two-spinon DSF and so, such a comparison allows faster acquaintance with the four-spinon contribution. To have the discussion as clear as possible, we scale both $`S_4`$ and $`S_2`$ to one. All results concerning $`S_2`$ are already known BCK ; KMB ; only those concerning $`S_4`$ are new.

This article is organized as follows. After these introductory remarks, we describe in the next section the Heisenberg model and give the definition of the dynamic structure function. We write its decomposition in $`2p`$-spinon contributions and give the expressions of $`S_2`$ and $`S_4`$. In section three, we first determine the four-spinon continuum, the region in the $`(k,\omega )`$ plane outside which $`S_4`$ is identically zero. We will see that it extends beyond the spin-wave des Cloizeaux and Pearson boundaries. To the best of our knowledge, this result is a first direct and explicit exact theoretical confirmation that the total $`S`$ for the infinite chain tails outside the spin-wave continuum. Next in this section is a description of the behavior of $`S_4`$ as a function of $`\omega `$ for fixed values of $`k`$ followed by a comparison with the corresponding behavior of $`S_2`$. It is seen that there is consistency with the four-spinon continuum and that the overall shape (not the details) of $`S_4`$ is similar to the one of $`S_2`$, each in its own respective continuum. Last is carried a description of the behavior of $`S_4`$ as a function of $`k`$ for fixed values of $`\omega `$ and a comparison with the corresponding one of $`S_2`$. Here too similarity between the two overall shapes is found as well as consistency with the four-spinon continuum. Section four includes concluding remarks and indicates few directions in which this work can be carried forward.

This paper is a continuation of the work ABS . We use the same notation except for a slight modification in (18) where here we introduce the function $`h`$ instead of the function $`f`$ with the relation $`f\mathrm{exp}\left(h\right)`$.

## II Four-spinon dynamic structure function

The antiferromagnetic spin-$`\frac{1}{2}`$ $`XXX`$ Heisenberg chain is defined as the isotropic limit of the $`XXZ`$ anisotropic Heisenberg Hamiltonian:

$$H=\frac{1}{2}\underset{n=\mathrm{}}{\overset{\mathrm{}}{}}\left(\sigma _n^x\sigma _{n+1}^x+\sigma _n^y\sigma _{n+1}^y+\mathrm{\Delta }\sigma _n^z\sigma _{n+1}^z\right).$$
(1)
$`\mathrm{\Delta }=(q+q^1)/2`$ is the anisotropy parameter and the isotropic antiferromagnetic limit is obtained as $`\mathrm{\Delta }1^{}`$, or equivalently $`q1^{}`$. Here $`\sigma _n^{x,y,z}`$ are the usual Pauli matrices acting at the site $`n`$ of the chain. The exact diagonalization of this Hamiltonian is performed directly in the thermodynamic limit. This is necessary if we want to exploit the $`\mathrm{U}_q(\widehat{\mathrm{sl}}_2)`$ quantum group symmetry present in the model JimMi . One consequence is the appearance of two vacuum states $`|0_i`$, $`i=0,1`$ due to two different boundary conditions on the infinite chain. The Hilbert space $``$ consists of $`n`$-spinon energy eigenstates $`|\xi _1,\mathrm{},\xi _n_{ϵ_1,\mathrm{},ϵ_n;i}`$ such that:

$$H|\xi _1,\mathrm{},\xi _n_{ϵ_1,\mathrm{},ϵ_n;i}=\underset{j=1}{\overset{n}{}}e(\xi _j)|\xi _1,\mathrm{},\xi _n_{ϵ_1,\mathrm{},ϵ_n;i},$$
(2)
where $`e(\xi _j)`$ is the energy of spinon $`j`$ and $`\xi _j`$ is a spectral parameter living on the unit circle. In the above relation, $`ϵ_j=\pm 1`$. The translation operator $`T`$ which shifts the spin chain by one site acts on the energy eigenstates in the following manner:

$$T|\xi _1,\mathrm{},\xi _n_{ϵ_1,\mathrm{},ϵ_n;i}=\underset{i=1}{\overset{n}{}}\tau (\xi _i)|\xi _1,\mathrm{},\xi _n_{ϵ_1,\mathrm{},ϵ_n;1i},$$
(3)
where $`\tau (\xi _j)=e^{ip(\xi _j)}`$ and $`p(\xi _j)`$ is the lattice momentum of spinon $`j`$. The exact expressions of the spinon energy and lattice momentum in terms of the spectral parameter are known in the literature Fadeev-Takhtajan ; JimMi ; ABS . We are interested in their $`XXX`$ limit and it is given below in eq (12). The completeness relation in $``$ reads:

$$𝐈=\underset{i=0,1}{}\underset{n0}{}\underset{\{ϵ_j=\pm 1\}_{j=1,n}}{}\frac{1}{n!}\underset{j=1}{\overset{n}{}}\frac{d\xi _j}{2\pi i\xi _j}|\xi _1,\mathrm{},\xi _n_{ϵ_1,\mathrm{},ϵ_n;ii;ϵ_1,\mathrm{},ϵ_n}\xi _1,\mathrm{},\xi _n|.$$
(4)
The two-point dynamic structure function is defined as the Fourier transform of the zero-temperature vacuum-to-vacuum two-point function. The transverse DSF is therefore defined by:

$$S^{i,+}(\omega ,k)=_{\mathrm{}}^{\mathrm{}}𝑑t\underset{m}{}e_i^{i(\omega t+km)}0|\sigma _m^+(t)\sigma _0^{}(0)|0_i,$$
(5)
where $`\omega `$ and $`k`$ are the neutron energy and momentum transfer respectively, and $`\sigma ^\pm `$ denotes $`(\sigma ^x\pm i\sigma ^y)/2`$. The DSF satisfies the following relations:

$$S^{i,+}(\omega ,k)=S^{i,+}(\omega ,k)=S^{i,+}(\omega ,k+2\pi ),$$
(6)
expressing reflection symmetry and periodicity. Inserting the completeness relation (4) and using the Heisenberg relation:

$$\sigma _m^{x,y,z}(t)=\mathrm{exp}\left(iHt\right)T^m\sigma _0^{x,y,z}(0)T^m\mathrm{exp}\left(iHt\right),$$
(7)
we can write the transverse DSF as the sum of $`n`$-spinon contributions:

$$S^{i,+}(\omega ,k)=\underset{n\mathrm{even}}{}S_n^{i,+}(\omega ,k),$$
(8)
where the $`n`$-spinon DSF $`S_n`$ is given by:

$`S_n^{i,+}(\omega ,k)`$ $`=`$ $`{\displaystyle \frac{2\pi }{n!}}{\displaystyle \underset{m}{}}{\displaystyle \underset{ϵ_1,\mathrm{},ϵ_n}{}}{\displaystyle \underset{j=1}{\overset{n}{}}\frac{d\xi _j}{2\pi i\xi _j}e^{im(k+_{j=1}^np_j)}\delta \left(\omega \underset{j=1}{\overset{n}{}}e_j\right)}`$ (9)
$`\times X_{ϵ_n,\mathrm{},ϵ_1}^{i+m}(\xi _n,\mathrm{},\xi _1)X_{ϵ_1,\mathrm{},ϵ_n}^{1i}(q\xi _1,\mathrm{},q\xi _n),`$
a relation in which $`X^i`$ denotes the form factor:

$$X_{ϵ_1,\mathrm{},ϵ_n}^i(\xi _1,\mathrm{},\xi _n)_i0|\sigma _0^+\left(0\right)|\xi _1,\mathrm{},\xi _n_{ϵ_1,\mathrm{},ϵ_n;i}.$$
(10)
In relation (9), $`i+m`$ is to be read modulo 2. Note that each $`S_n`$ satisfies the symmetry relations (6).

An exact expression for the form factor $`X^i`$ is known smirnov ; JimMi . To arrive at the result, one has to exploit extensively the infinite-dimensional representation of $`\mathrm{U}_q(\widehat{\mathrm{sl}}_2)`$ and bosonize the relevant vertex operators in order to be able to manipulate in a systematic way traces of these operators which ultimately yield the correlation functions. Using this form factor, it is possible to give an exact expression for the $`n`$-spinon DSF in the anisotropic case ABS  and determine its isotropic limit Boug , obtained via the replacement JimMi ; ABS :

$$\xi =ie^{2i\epsilon \rho };q=e^\epsilon ,\epsilon 0^+,$$
(11)
where $`\rho `$ becomes the spectral parameter suited for this limit. The expressions of the energy $`e`$ and momentum $`p`$ in terms of $`\rho `$ then read:

$$e(\rho )=\frac{\pi }{\mathrm{cosh}(2\pi \rho )}=\pi \mathrm{sin}p;\mathrm{cot}p=\mathrm{sinh}(2\pi \rho );\pi p0.$$
(12)
It turns out that the transverse two-spinon DSF $`S_2`$ does not involve a contour integration, see (9). Its exact expression has been derived in BCK . It reads:

$$S_2^+(\omega ,k\pi )=\frac{1}{4}\frac{e^{I(\rho )}}{\sqrt{\omega _{2u}^2\omega ^2}}\mathrm{\Theta }(\omega \omega _{2l})\mathrm{\Theta }(\omega _{2u}\omega ),$$
(13)
where $`\mathrm{\Theta }`$ is the Heaviside step function and the function $`I(\rho )`$ is given by:

$$I(\rho )=_0^+\mathrm{}\frac{dt}{t}\frac{\mathrm{cosh}(2t)\mathrm{cos}(4\rho t)1}{\mathrm{sinh}(2t)\mathrm{cosh}(t)}e^t.$$
(14)
$`\omega _{2u(l)}`$ is the upper (lower) bound of the two-spinon excitation energies called the des Cloizeaux and Pearson (dCP) Muller ; BCK  upper (lower) bound or limit. They read:

$$\omega _{2u}=2\pi \mathrm{sin}\frac{k}{2};\omega _{2l}=\pi |\mathrm{sin}k|.$$
(15)
The quantity $`\rho `$ is related to $`\omega `$ and $`k`$ by the relation:

$$\mathrm{cosh}\pi \rho =\sqrt{\frac{\omega _{2u}^2\omega _{2l}^2}{\omega ^2\omega _{2l}^2}},$$
(16)
which is obtained using eq (12) and the energy-momentum conservation laws. The properties of $`S_2`$ have been discussed in KMB ; BKM  where a comparison with the Müller ansatz Muller  is carried.

The four-spinon DSF $`S_4`$ involves only one contour integration and its expression is given in ABS . For $`0k\pi `$ it reads:

$$S_4^+(\omega ,k\pi )=C_4_\pi ^0𝑑p_3_\pi ^0𝑑p_4F(\rho _1,\mathrm{},\rho _4).$$
(17)
For other values of $`k`$, it extends by symmetry using (6). $`C_4`$ is a numerical constant irrelevant for the present work since we will scale $`S_4`$ to unity, and the integrand $`F`$ is given by:

$$F(\rho _1,\mathrm{},\rho _4)=\underset{(p_1,p_2)}{}\frac{\mathrm{exp}[h(\rho _1,\mathrm{},\rho _4)]_{\mathrm{}=1}^4|g_{\mathrm{}}(\rho _1,\mathrm{},\rho _4)|^2}{\sqrt{W_u^2W^2}}.$$
(18)
The different quantities involved in this expression are as follows.

$`W`$ $`=`$ $`\omega +\pi (\mathrm{sin}p_3+\mathrm{sin}p_4);`$
$`W_u`$ $`=`$ $`2\pi \left|\mathrm{sin}\left(K/2\right)\right|,K=k+p_3+p_4;`$
$`\mathrm{cot}p_j`$ $`=`$ $`\mathrm{sinh}(2\pi \rho _j),\pi p_j0.`$ (19)
The function $`h`$ is given by:

$$h(\rho _1,\mathrm{},\rho _4)=\underset{1j<j^{^{}}4}{}I(\rho _{jj^{}}),$$
(20)
where $`\rho _{jj^{}}=\rho _j\rho _j^{}`$ and the function $`g_{\mathrm{}}`$ reads:

$`g_{\mathrm{}}`$ $`=`$ $`(1)^{\mathrm{}+1}(2\pi )^4{\displaystyle \underset{j=1}{\overset{4}{}}}\mathrm{cosh}(2\pi \rho _j)`$ (21)
$`\times {\displaystyle \underset{m=\mathrm{\Theta }(j\mathrm{})}{\overset{\mathrm{}}{}}}{\displaystyle \frac{_i\mathrm{}(m\frac{1}{2}\mathrm{\Theta }(\mathrm{}i)+i\rho _{ji})}{_{ij}\pi ^1\mathrm{sinh}(\pi \rho _{ji})}}{\displaystyle \underset{i=1}{\overset{4}{}}}{\displaystyle \frac{\mathrm{\Gamma }(m\frac{1}{2}+i\rho _{ji})}{\mathrm{\Gamma }(m+1+i\rho _{ji})}},`$
where $`\mathrm{\Theta }`$ is the Heaviside step function. In (18), the sum $`_{(p_1,p_2)}`$ is over the two pairs $`(p_1,p_2)`$ and $`(p_2,p_1)`$ solutions of the energy-momentum conservation laws:

$$W=\pi (\mathrm{sin}p_1+\mathrm{sin}p_2);K=p_1p_2.$$
(22)
They read:

$$(p_1,p_2)=(K/2+\mathrm{arccos}\left(W/\left[2\pi \mathrm{sin}\left(K/2\right)\right]\right),K/2\mathrm{arccos}\left(W/\left[2\pi \mathrm{sin}\left(K/2\right)\right]\right)).$$
(23)
Note that the solution in (23) is allowed as long as $`W_lWW_u`$ where $`W_u`$ is given in (19) and:

$$W_l=\pi |\mathrm{sin}K|.$$
(24)
The (analytic) behavior of the function $`F`$ in (18) is discussed in ABS . It is shown that the series $`g_{\mathrm{}}`$ is convergent. It is also shown that $`g_{\mathrm{}}`$ stays finite when two $`\rho _i`$’s or more get equal. Since the function $`\mathrm{exp}\left(h\right)`$ goes to zero in these regions KMB , the integrand $`F`$ of $`S_4`$ is regular there. Furthermore, it is shown that $`F`$ is exponentially convergent when one of the $`\rho _i`$’s goes to infinity, which means the two integrals over $`p_3`$ and $`p_4`$ in (17) do not yield infinities. All these analytic results help secure safe numerical manipulations.

## III Behavior of exact four-spinon DSF

From now on, we restrict ourselves to the interval $`0k\pi `$. All forthcoming results can be carried to the other intervals of $`k`$ using the symmetry relations (6). Also, we scale both $`S_4`$ and $`S_2`$ to appropriate units in order to display conveniently their respective behaviors.

### III.1 Four-spinon continuum

The first feature we discuss is the ‘four-spinon continuum’, by analogy with the two-spinon (or the spin-wave) continuum. It is the extent of the region in the $`(k,\omega )`$-plane outside which $`S_4`$ is identically zero. Remember that from (13), $`S_2`$ is confined to the region $`\omega _{2l}(k)\omega \omega _{2u}(k)`$, where $`\omega _{2l,u}(k)`$ are the dCP boundaries given in (15). From the condition $`W_lWW_u`$ mentioned after (23), we deduce that in order for $`S_4`$ to be nonzero identically, we must have $`\omega _{4l}(k)\omega \omega _{4u}(k)`$, where:

$`\omega _{4l}(k)`$ $`=`$ $`3\pi \mathrm{sin}(k/3)\mathrm{for}0k\pi /2;`$
$`\omega _{4l}(k)`$ $`=`$ $`3\pi \mathrm{sin}(k/3+2\pi /3)\mathrm{for}\pi /2k\pi ;`$
$`\omega _{4u}(k)`$ $`=`$ $`4\pi \mathrm{cos}(k/4)\mathrm{for}0k\pi .`$ (25)
We see that $`\omega _{4l}(k)`$ and $`\omega _{4u}(k)`$ are a sort of four-spinon dCP boundaries for $`S_4`$. The two and four-spinon continua are drawn in figure 1 below. We immediately notice that the four-spinon continuum is not restricted to the region between the two-spinon dCP branches, which means, a fortiori, that the full $`S`$ is also not confined to the spin-wave continuum. This is a direct and explicit theoretical confirmation of the ‘tail’ of the dynamic structure function observed outside the spin-wave continuum in finite-chain numerical calculations Muller  and the phenomenology Exp . Due to the six and higher spinon contributions, it is actually legitimate to expect the full DSF to tail even further outside the four-spinon continuum, with arguably much smaller values. For example, in the interval $`0k/\pi <0.51741`$, there is a narrow region between $`\omega _{2u}`$ and $`\omega _{4l}`$ inside which both $`S_2`$ and $`S_4`$ are identically zero whereas the total $`S`$ may have (very small) nonzero values, something that could eventually be checked in finite-chain calculations. However, it is not possible to estimate exactly the continua corresponding to the $`S_{n>4}`$ without manageable explicit formulae for these. But what is already clear from our study is the fact that indeed, the spin-wave continuum is not restrictive to the total dynamic structure function.

There are further general features we can also read from figure 1 without having recourse to detailed calculations. For example, we see that inside the interval $`0k/\pi <0.51741`$, the four-spinon continuum lies entirely above the two-spinon continuum. This means that for this interval, $`S_2`$ may be accepted as a good approximation for the total $`S`$ between the spin-wave boundaries and $`S_4`$ between $`\omega _{4l}`$ and $`\omega _{4u}`$. For $`0.51741k/\pi 1`$ however and as $`k`$ increases, there is increasing overlap between the two continua so that inside the spin-wave continuum, we may expect $`S`$ to divert in perhaps a non-negligible way from $`S_2`$.

### III.2 Behavior as a function of the energy transfer

Next we describe the behavior of $`S_4`$ as a function of $`\omega `$ for fixed values of $`k`$. Figure 2 displays the shapes of $`S_4`$ for respectively $`k/\pi =1/4,1/2,3/4`$ and $`1`$. As an illustration, take for example the case $`k/\pi =1/4`$. We see on figure 2 that $`S_4`$ is zero until $`\omega /\pi 0.68`$, which corresponds to the beginning of the four-spinon continuum $`\omega _{4l}(\pi /4)/\pi =3\mathrm{sin}(\pi /12)=0.77646`$. We obtain a slightly smaller value when reading directly from figure 2 because of fitting. Note that it is important to check each time consistency with the four-spinon continuum. This is because the lower and upper branches $`\omega _{4l,u}(k)`$ are not imposed explicitly in the integration algorithm of $`S_4`$ unlike the case of $`S_2`$ where the corresponding expression (13) incorporates explicitly the two-spinon dCP boundaries $`\omega _{2l,u}(k)`$. For $`S_4`$, only the conditions $`W_lWW_u`$ noted just before (24) are imposed.

Starting from $`\omega /4\pi 0.17`$, $`S_4`$ jumps sharply from zero to a maximum at $`\omega /4\pi 0.21`$ (read from figure 2). Then it decreases with two apparent local minima until it becomes negligible at roughly $`\omega /4\pi 0.6`$. The upper branch of the four-spinon continuum at $`k=\pi /4`$ is $`\omega _{4u}(\pi /4)/4\pi =\mathrm{cos}(\pi /16)=0.9808`$. Hence, $`S_4`$ is practically negligible (but not identically zero) for values of $`\omega /4\pi `$ between about $`0.6`$ and $`0.9808`$. The description of $`S_4`$ as a function of $`\omega `$ for the other values of $`k`$ can be carried along the same lines and each time, consistency with the four-spinon continuum is checked. The shapes are all similar to one another: steep increase from zero to a maximum value, followed by a ‘wiggled’ slower decrease to zero.

Furthermore, it is interesting to notice that for all values of $`k`$, the overall shape of $`S_4`$, not the detail, is roughly similar to the one of $`S_2`$, represented as a function of $`\omega `$ for the same values of $`k`$ in figure 3. For instance, for $`k/\pi =1/4`$, $`S_2`$ starts very sharply at about $`\omega /2\pi 0.35`$, which corresponds to the start of the two-spinon continuum at $`\omega _{2l}(\pi /4)/2\pi =0.35355`$. It reaches a maximum before decreasing more slowly towards zero. The only difference worth mentioning is the decrease of $`S_4`$ after its absolute maximum which presents local minima and maxima whereas the decrease of $`S_2`$ is always ‘smooth’. This may originate from the ‘richer’ structure of the expression of $`S_4`$ with respect to that of $`S_2`$.

### III.3 Behavior as a function of the momentum transfer

Last we describe the behavior of $`S_4`$ as a function of $`k`$ for fixed values of $`\omega `$. In figure 4 are plotted the graphs of $`S_4`$ in terms of $`k`$ for $`\omega /\pi =1/4`$, $`1/2`$, $`3/4`$, $`1`$ respectively. Let us describe for example the case $`\omega /\pi =1/4`$. We see on figure 4 that $`S_4`$ is zero until we reach the value $`k/\pi 0.88`$. On the other hand, the four-spinon continuum lies outside the interval $`0.07967k/\pi 0.92033`$. In the region $`0k/\pi \mathrm{0.0796\hspace{0.17em}7}`$, figure 4 shows no discernible finite values for $`S_4`$, only a very thin ‘trace’ that would be more visible with a better resolution. This means that $`S_4`$ is negligible for those small values of $`k`$. For larger values of $`\omega `$ though, $`S_4`$ picks up clear finite values in the interval $`0k3\mathrm{arcsin}(\omega /3\pi )`$; see the other graphs on figure 4. Those values get larger as $`\omega `$ increases.

Returning to the case $`\omega /\pi =1/4`$, we see that $`S_4`$ starts from zero at $`k/\pi 0.88`$ (slightly smaller than the exact value $`0.92033`$ because of fitting), rises sharply to a maximum and then decreases. The behavior of $`S_4`$ for the other values of $`\omega `$ is also consistent with the four-spinon continuum. Take for example the case $`\omega /\pi =1/2`$. The four-spinon continuum indicates that $`S_4`$ is identically zero for $`0.1599k/\pi 0.84008`$. We see indeed small values for $`S_4`$ from $`k/\pi =0`$ up to a little before 0.2, and $`S_4`$ rising again from zero a little after $`k/\pi =0.8`$ to a local maximum, then to an absolute maximum before decreasing.

Here too $`S_2`$ has a similar overall behavior within its own (two-spinon) continuum. Figure 5 displays the behavior of $`S_2`$ with respect to $`k`$ for the same fixed values of $`\omega `$. But as $`\omega `$ increases, we notice the richer structure of $`S_4`$ with respect to the corresponding one for $`S_2`$. This is due to the more involved expression of $`S_4`$. In fact, for larger values of $`\omega `$, the integrand $`F`$ in (17) is nonvanishing in larger and larger areas in the $`(p_3,p_4)`$ plane. For illustration, compare the behavior of $`F`$ shown in figure 6 for $`(k,\omega )=(0.5\pi ,2\pi )`$ with the one shown in figure 7 for $`(k,\omega )=(0.5\pi ,3\pi )`$. In this regard, we recall that the quadrature-based algorithms written in ABSS  were based on this observation and took advantage of the fact that for fairly small values of $`\omega `$, the integrand $`F`$ was negligible in large areas of the $`(p_3,p_4)`$ plane. This is the reason why those algorithms could not be carried to larger values of $`\omega `$, or be able to describe efficiently a behavior as a function of the neutron energy $`\omega `$ itself. As already mentioned, we have used here Monte Carlo techniques and we can do better, but consistency with ABSS  is realized, i.e., we have the same behavior of $`S_4`$ as a function of $`k`$ for the same small values of $`\omega `$.

## IV Conclusion

In this work, we have described the behavior of the exact four-spinon dynamic structure function $`S_4`$ in the antiferromagnetic isotropic Heisenberg quantum spin chain at zero temperature as a function of the neutron energy $`\omega `$ and momentum transfer $`k`$. We have also determined the four spinon continuum, the region outside which $`S_4`$ is identically zero. The discussion was carried in the form of a comparison with the corresponding behavior of the exact two-spinon dynamic structure function $`S_2`$, already known in the literature. Figures 8 and 9 summarize these two behaviors where $`S_4`$ (8) and $`S_2`$ (9) are drawn in the $`(k,\omega )`$ plane. Recall that we have scaled them down to one.

There are four directions in which one may wish to carry forward with this work. The first is the anisotropic case. The model is exactly soluble and we do have generic expressions for $`S_n`$ in the form of contour integrals in the spectral parameters complex planes Boug . The difficulty here is that the integrands involve much more complicated functions which are already present in $`S_2`$, and one should expect intricate complexities in this more general case.

The second direction in which one may want to push forward is the situation where there is an external magnetic field. There are finite-chain calculations in this regard, Muller  and more recently the works KM-KBM . But one has to remember that the model in not exactly solvable in this case. So one may want to try small perturbations around the zero-field limit solution. The third direction is the finite-temperature case. Here too there are finite-chain results and it is interesting to see those effects on $`S_2`$ and $`S_4`$. The fourth direction is to look into the situation of a spin-one chain. The model is still exactly solvable and exploiting the quantum group symmetry, compact expressions for the form factors are available idzumi ; bougourzi-weston ; bougourzi-any-spin .