text:
# Role of many-body entanglement in decoherence processes

## Abstract

A pure state decoheres into a mixed state as it entangles with an environment. When an entangled two-mode system is embedded in a thermal environment, however, each mode may not be entangled with its environment by their simple linear interaction. We consider an exactly solvable model to study the dynamics of a total system, which is composed of an entangled two-mode system and a thermal environment, and also an array of infinite beam splitters. It is shown that many-body entanglement of the system and the environment plays a crucial role in the process of disentangling the system.

Decoherence has been studied in the context of quantum-classical correspondence, providing a quantum-to-classical transition of a system Zurek . A single-mode pure state becomes mixed and loses its quantum nature by decoherence in an environment. Although the dynamics of the system has been studied extensively, there has not been a thorough investigation on the quantum correlation between the system and the environment, which is behind the dynamics of the system. The decoherence process can be understood as a process of entanglement between the system and its environment which is composed of a many (normally, infinite) number of independent modes. The increase of the system entropy Zurek2  may be due to the system-environment entanglement.

Most of the studies on decoherence have focused on a single-mode or single-particle system Johnson . This is because if a many-body pure system is initially separable, its decoherence process is a straightforward extension of a single-body system. However, if there is entanglement in the initial pure system, the decoherence mechanism can be of a different nature. For an entangled two-mode pure system, each mode is generically in a mixed state and its passive linear interaction with an environment, which is normally in a mixed state, does not seem to bring about entanglement. What kind of correlation then causes the loss of entanglement initially in the system? In this paper, we answer this question by studying the quantum correlation of a two-mode entangled continuous-variable system with an environment in thermal equilibrium.

A continuous-variable state is defined in an infinite-dimensional Hilbert space and it is convenient to study such a state using its quasiprobability Wigner function Barnett97 , $`W(\stackrel{~}{𝐱})`$, in phase space. For an $`N`$-mode field, the coordinates of phase space are composed of quadrature variables, $`\stackrel{~}{𝐱}=\{q_1,p_1,\mathrm{},q_N,p_N\}`$. Throughout the paper a vector is denoted in bold face and an operator by a hat. A fiber or a free space, through which a light field propagates, is normally considered a thermal environment. The dynamics of the field mode coupled to the thermal environment is, in the Born-Markov approximation, governed by the Fokker-Planck equation Barnett97 ,

$$\frac{W(\stackrel{~}{𝐱})}{\tau }=\frac{\gamma }{2}\underset{i=1}{\overset{2N}{}}\left(\frac{}{\stackrel{~}{x}_i}\stackrel{~}{x}_i+\frac{\stackrel{~}{n}}{4}\frac{^2}{\stackrel{~}{x}_i^2}\right)W(\stackrel{~}{𝐱})$$
(1)
where $`\gamma `$ is the energy decay rate of the system and $`\stackrel{~}{n}=2\overline{n}+1`$ with $`\overline{n}=[\mathrm{exp}(\mathrm{}\omega /k_BT)1]^1`$ is the average number of thermal photons at temperature $`T`$. $`k_B`$ is the Boltzmann constant.

It was shown by one of us that a single-mode Gaussian field interacting with a thermal environment can be modeled by the field passing through an array of infinite beam splitters KimImoto . A beam splitter is a simple passive linear device which keeps the Gaussian nature of an input field. Each beam splitter has two input ports. As the signal field is injected into one input port, it allows a degree of freedom for the other port where noise is injected. The collection of such degrees of freedom forms the environment. We assume a homogeneous thermal environment of temperature $`T`$ with all noise modes having the same physical properties. In Ref. KimImoto , the Fokker-Planck equation (1) for a single-mode field was derived using the beam splitter. The model was used to study the dynamics of entanglement between a single-mode field and its environment Nagaj .

Using the Fokker-Planck equation (1), one may study the dynamics of the system. However, it is hard to know the quantum correlation between the system and the environment as Eq. (1) is obtained by tracing over all environmental variables. In this paper, instead of tracing over all the environmental modes, we keep them to study the dynamics of entanglement of the system and the environment. A two-mode squeezed state, which may be generated by a nondegenerated optical parametric amplifier, is the most renowned and experimentally-relevant entangled state for continuous variables BraunsteinKimble1 . Its degree of entanglement increases as the degree of squeezing $`s`$ increases KimLee02  and it becomes a regularized Einstein-Podolsky-Rosen state when $`s\mathrm{}`$ BraunsteinKimble . In order to simplify the problem, we assume that only the mode $`a_2`$ of the system modes interacts with the environment while the other mode $`a_1`$ is isolated from it.

We consider an exactly solvable model of a two-mode system interacting with a homogeneous thermal environment, which results in the Fokker-Planck equation (1). For the Born-Markov approximation, we employ a time-dependent coupling constant in the model. Let us start with a finite number $`N`$ of environmental modes interacting with the system. The interaction Hamiltonian, in the interaction picture, is

$`\widehat{H}_I(t)={\displaystyle \underset{m=0}{\overset{N1}{}}}i\lambda _N(t)\left(\widehat{a}_2\widehat{b}_m^{}\widehat{a}_2^{}\widehat{b}_m\right),`$ (2)
where $`\widehat{b}_m`$ is the bosonic annihilation operator for environmental mode $`b_m`$ and the coupling constant $`\lambda _N(t)`$ is determined so as to reproduce the Fokker-Planck equation (1).

It is convenient to introduce collective modes $`c_n`$ which are conjugate to $`b_m`$ under the Fourier transformation such that

$`\widehat{c}_n\sqrt{{\displaystyle \frac{2}{N}}}{\displaystyle \underset{m=0}{\overset{N1}{}}}\mathrm{cos}\left({\displaystyle \frac{2\pi }{N}}nm\right)\widehat{b}_m,`$ (3)
where $`\widehat{c}_n`$ is an annihilation operator for a collective mode $`c_n`$. The collective modes are related with the entangling nature of the modes $`b_m`$, for example, the quantum that $`\widehat{c}_n^{}`$ creates from a vacuum is in an entangled state of $`b_m`$ modes.

The collective modes satisfy the boson commutation relation, $`[\widehat{c}_n,\widehat{c}_n^{}^{}]=\delta _{nn^{}}`$ and carry physical properties as bosonic modes. Using the collective mode, a state $`\widehat{\rho }`$ is described by the characteristic function

$`\chi _c(𝐗)=\text{Tr}\widehat{\rho }\mathrm{exp}[i𝐗\widehat{𝐗}^T],`$ (4)
where $`\widehat{𝐗}=(\widehat{Q}_0,\widehat{P}_0,\widehat{Q}_1,\widehat{P}_1,\mathrm{},\widehat{Q}_{N1},\widehat{P}_{N1})`$ with $`\widehat{Q}_n=(\widehat{c}_n+\widehat{c}_n^{})/\sqrt{2}`$ and $`\widehat{P}_n=i(\widehat{c}_n^{}\widehat{c}_n)/\sqrt{2}`$ and $`𝐗=(P_0,Q_0,P_1,Q_1,\mathrm{},P_{N1},Q_{N1})`$. It is straightforward to show that, for a given density operator $`\widehat{\rho }`$, $`\chi _c`$ is the same as the usual characteristic function, $`\chi _b`$, which is obtained in terms of modes $`b_m`$: $`\chi _c(𝐗)=\chi _b(𝐱)`$ where $`𝐱`$ is conjugate to $`𝐗`$ by Fourier transformation (3). The collective modes $`c_n`$ provide a different perspective from the modes $`b_m`$, preserving all physical properties for a given state.

The time-evolution operator $`\widehat{U}_I(\tau )`$ for the interaction Hamiltonian (2) is equivalent to a beam splitter operator with the system mode $`a_2`$ and a collective mode $`c_0`$ as its input ports. That is,

$`\widehat{U}_I(\tau )=\mathrm{exp}[\theta (\tau )(\widehat{a}_2\widehat{c}_0^{}\widehat{a}_2^{}\widehat{c}_0)],`$ (5)
where $`\theta (\tau )=\sqrt{N/2}_0^\tau \lambda _N(t)𝑑t`$ determines the transmittivity, $`t^2(\tau )=\mathrm{cos}^2\theta (\tau )`$. We take the limit, $`N\mathrm{}`$, keeping the transmitted energy finite, $`t^2(\tau )=\mathrm{exp}(\gamma \tau )`$, in accordance with the Fokker-Planck equation. We find an important fact that the interaction of the system with the infinite modes $`b_m`$ of the environment can be reduced into the interaction with the single collective mode $`c_0`$. The properties of collective mode $`c_0`$ changes due to the interaction but each environmental mode $`b_m`$ hardly changes which is reflected in no change of other collective modes $`c_n`$.

A Gaussian field has the characteristic function in the form of $`\chi (𝐱)=\mathrm{exp}(\mathrm{𝐱𝐕𝐱}^T/4)`$ where $`𝐕`$ is the correlation matrix whose elements determine the mean quadrature values of the field: $`V_{ij}=(\widehat{x}_i\widehat{x}_j+\widehat{x}_j\widehat{x}_i)`$. Note that we neglected linear displacement terms in the Gaussian characteristic function as they do not play a crucial role in the entanglement. The entanglement nature of a Gaussian field is thus uniquely represented by its correlation matrix $`𝐕`$. For a two-mode squeezed state, the correlation matrix $`𝐕_s`$ is simply KimLee02

$$𝐕_s=\left(\begin{array}{cc}\mathrm{cosh}(2s)𝟙& \mathrm{sinh}(2s)𝝈_z\\ \mathrm{sinh}(2s)𝝈_z& \mathrm{cosh}(2s)𝟙\end{array}\right)$$
(6)
where $`𝟙`$ is the $`2\times 2`$ unit matrix and $`𝝈_z`$ is the Pauli matrix.

As the system interacts only with the collective mode $`c_0`$ in the homogeneous thermal environment, it suffices to consider the correlation matrix of the two system modes $`a_1`$ and $`a_2`$ and the collective mode $`c_0`$. The collective mode $`c_0`$ is initially in a thermal state with the average number $`\overline{n}`$ of the collective bosons. Thus, the correlation matrix of $`a_1`$, $`a_2`$, and $`c_0`$ before the interaction is given by $`𝐕_0=𝐕_s\stackrel{~}{n}𝟙`$. The evolution operator $`\widehat{U}_I(\tau )`$ is now described by the matrix

$$𝐔_I=\left(\begin{array}{ccc}& 0& 0\\ 0& t𝟙& r𝟙\\ 0& r𝟙& t𝟙\end{array}\right)$$
(7)
where $`r^2=1t^2`$. Then the correlation matrix for the system and environment after the interaction is obtained as $`𝐕_c=𝐔_I𝐕_0𝐔_I^T`$.

A separability condition was derived by Simon Simon  that a two-mode Gaussian state is separable if and only if the partially momentum-reversed correlation matrix (or equivalently the partially transposed density operator) satisfies the uncertainty principle. The condition was extended to a biseparability condition between a single mode and a group of $`N`$ modes by Werner and Wolf WernerWolf , which reads that a Gaussian field of $`1\times N`$ modes is biseparable if and only if

$$𝚲𝐕𝚲\frac{1}{2}𝝈_y^{(N+1)}0,$$
(8)
where $`𝚲`$ is a partial momentum-reversal matrix, $`𝐕`$ is the correlation matrix of $`1\times N`$ modes and $`𝝈_y`$ is the Pauli matrix. Here, $`𝝈_y^{(N+1)}𝝈_y𝝈_y\mathrm{}𝝈_y`$.

We start with a short discussion on the dynamics of the entanglement for the system. In order to consider quantum statistical properties of the field of modes $`a_1`$ and $`a_2`$, we trace the total density operator over all environmental modes, which is equivalent to considering the correlation matrix $`𝐕_c`$ only for the modes $`a_1`$ and $`a_2`$:

$$𝐕_c(a_1,a_2)=\left(\begin{array}{cc}\mathrm{cosh}(2s)𝟙& t\mathrm{sinh}(2s)𝝈_z\\ t\mathrm{sinh}(2s)𝝈_z& (t^2\mathrm{cosh}2s+r^2\stackrel{~}{n})𝟙\end{array}\right).$$
(9)
It has to be emphasized that this correlation matrix is exactly the same as the solution of the Fokker-Planck equation (1). Using Simon’s criterion in Eq. (8) KimLee02 ; Simon , the field of the modes $`a_1`$ and $`a_2`$ is separable when the transmittivity of the beam splitter is $`t^2\frac{\overline{n}}{1+\overline{n}}`$ for the squeezing parameter $`s0`$. Note that the separability condition does not depend on the initial entanglement of the system as far as there is any entanglement in the initial instance. The separability of the two-mode squeezed state depends only on the temperature of the environment and the overall transmittivity of the beam splitters.

We now study the entanglement of the system and the environment. Here, instead of the entanglement of the system with an individual mode $`b_m`$ of the environment, we are interested in the biseparability of a system mode and the collection of the environmental modes. Let us first consider the entanglement of the modes $`a_1`$ and $`c_0`$. The correlation matrix $`𝐕_c(a_1,c_0)`$ is equivalent to $`𝐕_c(a_1,a_2)`$ in Eq. (9) if $`r`$ and $`t`$ are interchanged. The separability condition is found to be $`t^2\frac{1}{1+\overline{n}}`$, which is again independent from the initial entanglement of the system as far as $`s0`$.

The entanglement of modes $`a_2`$ and $`c_0`$ is easily discussed using a quasiprobability $`P`$ function Barnett97 , the existence of which is a sufficient condition for entanglement LeeJeong . Tracing over mode $`a_1`$ of the two-mode squeezed state, the other mode $`a_2`$ is in a thermal state with $`\stackrel{~}{n}_s2\overline{n}_s+1=\mathrm{cosh}2s`$. It is well-known that thermal states have positive definite $`P`$ functions and the action of a beam splitter transforms only the coordinates of the input $`P`$ functions: $`P_{a_2}(\alpha )P_{c_0}(\beta )\stackrel{bs}{}P_{a_2}(t\alpha r\beta )P_{c_0}(t\beta +r\alpha )`$. Thus the field of modes $`a_2`$ and $`c_0`$ is always separable. A simple decoherence picture, that a system decoheres as it entangles with its environment, does not hold any longer. However this does not mean that there is no entanglement involved in the decoherence process because the other mode of the system has to be taken into account.

Fig. 1 presents the entanglement structure for a two-mode squeezed state interacting with a thermal environment. The solid lines are the boundaries of entanglement of the system mode $`a_1`$ and the collective mode of the environment $`c_0`$ and of the two system modes $`a_1`$ and $`a_2`$. These lines are obtained by the separability condition (8) in the present exactly solvable model. For the comparison, we consider $`N=100`$ beam splitters modeling the interaction with the thermal environment and calculate the biseparability of the $`1\times 100`$-mode field using Giedke et al.’s computational analysis Giedke01 . The computational results of entanglement are denoted by circles and dots. We find that the results are also independent of the squeezing parameter. Fig. 1 shows that the two methods are exactly consistent.

Now, we need to consider what Fig. 1 really tells us. The three modes $`a_1`$, $`a_2`$, and $`c_0`$ compose a tripartite system. Multi-mode entanglement for a continuous-variable system has been studied extensively using beam splitters and single-mode squeezed states Braunstein00 . Giedke et al. classified types of entanglement for a three-mode Gaussian field Giedke01 . When a three-mode field is not biseparable, it is called fully entangled. Here, a three-mode field is biseparable when any grouping of three modes into two are separable. A fully-entangled tripartite qubit system may be entangled in two inequivalent ways Dur00 ; Acin01 , one of which is Greenberg-Horne-Zeilinger (GHZ) entanglement GHZ  and the other is W-entanglement Dur00 . A GHZ-entangled state becomes separable when any one particle is traced over any one particle while a W-entangled state is pairwise entangled for any pair of the three particles. On the other hand, there is another kind of entanglement for a fully-entangled tripartite $`d`$-dimensional system Thap99 . That is, two-way entanglement: One example for three particles labeled as $`a`$, $`b`$ and $`c`$ are pairwise entangled of $`a`$ and $`b`$ and of $`b`$ and $`c`$ but the particles $`a`$ and $`c`$ are separable.

From the biseparability condition in Eq. (8), we find that the tripartite system of $`a_1`$, $`a_2`$, and $`c_0`$ is fully entangled if $`0<t^2<1`$ and $`s0`$. It is notable that the tripartite system being fully entangled is independent of the temperature of the environment and thus over whole region in the Fig. 1. In Fig. 1, it is shown that the tripartite system is two-way entangled, i.e., $`a_1`$-$`a_2`$ and $`a_1`$-$`c_0`$ modes are entangled for $`\overline{n}<1`$ if $`\frac{\overline{n}}{1+\overline{n}}<t^2<\frac{1}{1+\overline{n}}`$. There are two regions where only one pair of the three-mode field is entangled.

We focus on the region of $`\overline{n}>1`$ described by

$$\frac{1}{1+\overline{n}}t_0^2\frac{\overline{n}}{1+\overline{n}},$$
(10)
where there is no pairwise entanglement. Because the tripartite system of $`a_1`$, $`a_2`$ and $`c_0`$ are fully entangled, if there is no pairwise entanglement, the only possibility is GHZ-entanglement among the three modes. Here, we draw the conclusion that GHZ entanglement of an entangled system with an environment is a source to its loss of the initial entanglement. This is clearly seen in Fig. 1 when we consider the dynamics of the entanglement as the interaction time $`\tau `$ increases or equivalently the transmittivity $`t^2`$ decreases. For $`\overline{n}>1`$, the two modes of the system are initially entangled. As they come to the region where they are GHZ-entangled with the environment, the two modes of the system lose their initial entanglement. Finally, the entanglement is transferred to between the system mode $`a_1`$ and the collection of the environmental modes.

When system mode $`a_2`$ is neither GHZ-entangled nor pairwise entangled with $`a_1`$, what happens to $`a_2`$? In order to understand it, let us consider each thermal mode of the environment as a part of a pure system. Then the total system comes to be in a pure state, which is easier to analyze. As pointed out earlier, by tracing over one mode of a two-mode squeezed state, we have a single-mode thermal field with $`\stackrel{~}{n}=\mathrm{cosh}2s_n`$. We may consider a purification of the thermal environmental state such that the thermal mode $`b_m`$ results from tracing a pure two-mode squeezed state $`\widehat{\rho }_{b_mb_m^{}}`$ over its counter mode $`b_m^{}`$. Letting $`c_0^{}`$ be the counter collective mode of $`c_0`$, the modes of $`a_1`$, $`a_2`$, $`c_0`$ and $`c_0^{}`$, that involve in the interaction, compose the pure total system. As the system mode $`a_2`$ interacts with the collective mode $`c_0`$, there should be an analogy between the relations of the modes $`a_1`$, $`a_2`$, and $`c_0`$ and those of the modes $`c_0^{}`$, $`c_0`$, and $`a_2`$. We have already seen that modes $`a_1`$ and $`c_0`$, which have never met, may get entangled by the long interaction. This argument can be applied to the entanglement of the modes $`c_0^{}`$ and $`a_2`$. Here, we obtain an interesting result that the system mode $`a_2`$, which was initially entangled with the system mode $`a_1`$, loses its initial entanglement and becomes entangled not with the interacting environmental modes but with the hidden environmental modes, $`c_0^{}`$. The entanglement of $`c_0`$ and $`c_0^{}`$ can be lost but the entanglement of $`b_m`$ and $`b_m^{}`$ is almost kept unchanged as other collective modes $`c_n`$ and $`c_n^{}`$ are still well-entangled. For the comparison with the present model, we calculate the separability of the system mode $`a_2`$ and the group of the hidden modes $`b_m^{}`$ by employing $`100`$ beam splitters and using Giedke et al.’s computational analysis. This results in the same as the dots in Fig. 1 with $`\overline{n}`$ of the $`x`$-axis replaced by $`(\stackrel{~}{n}_s1)/2`$. On the other hand, if $`\overline{n}`$ of the environment is zero, the environment is in a pure vacuum state and there is no hidden environmental modes. In this case, the system of $`a_1`$ and $`a_2`$ becomes never separable.

In this paper, we have studied the decoherence mechanism by highlighting the entanglement of the continuous-variable system with its environment. We showed that the homogeneous thermal environment can be summarized by a single collective mode with respect to the interaction with the system for the study of entanglement. As the two-mode entangled system is interacting with the thermal environment, the two modes may lose their entanglement by the GHZ-entanglement with the group of the environmental modes. The initial entanglement of the two-mode system transfers to the entanglement of one of the system modes and the group of environmental modes. We also studied the nature of entanglement for the other system mode.

###### Acknowledgements.

We thank Prof. V. Bužek for discussions.


format: text

text:
# Competing interactions in the 𝑋𝑌𝑍 model

## Abstract

We study the interplay between a $`XY`$ anisotropy $`\gamma `$, exchange modulations and an external magnetic field along the $`z`$ direction in the $`XYZ`$ chain using bosonization and Lanczos diagonalization techniques. We find an Ising critical line in the space of couplings which occur due to competing relevant perturbations which are present. More general situations are also discussed.

PACS numbers: 75.10.Jm, 75.60.Ej

Introduction. The competition between different relevant perturbations which can render a system critical in a certain domain of the couplings space has been studied by many authors . Within a low energy description of many different 1D lattice models, one often finds the double-frequency sine-Gordon model , that is, a $`U(1)`$ scalar field with two perturbations of different frequencies. In the case where both perturbations are relevant and for a certain frequency ratio, it has been conjectured that an Ising criticality can arise . Such a situation has been found in the study of different lattice systems (see e.g.  and references therein). More recently, the so-called self-dual sine-Gordon model has been studied in .

A simple realization of this situation is found in an exactly solvable case, i.e. in the dimerized $`XY`$ chain  and a qualitative analysis of the phase diagram has been presented for more general cases in a field . The effect of an $`XY`$ anisotropy has also been studied for the Fibonacci $`XY`$ chain in a magnetic field where it was shown that the rather involved staircase structure of the magnetization curve gradually disappears by increasing the anisotropy of the spin exchange interactions .

In the present paper we present a unified picture of the above mentioned effects and describe, by using bosonization as well as numerical techniques, generic situations in which a spin gap opening mechanism  competes with an $`XY`$ anisotropy to render Ising criticality. This issue was studied first in . In  the interplay between a gap-opening perturbation and $`XY`$ anisotropy was analyzed and a qualitative phase diagram was proposed. In this paper we focus on the Hamiltonian

$$H=\underset{n}{}\left((1+\gamma )S_n^xS_{n+1}^x+(1\gamma )S_n^yS_{n+1}^y+\delta (1)^n\left(S_n^xS_{n+1}^x+S_n^yS_{n+1}^y\right)+\mathrm{\Delta }S_n^zS_{n+1}^z\right)h\underset{n}{}S_n^z,$$
(1)
which captures most of the essential aspects referred to above, and provide quantitative results which support the statements made in . A more general situation as that arising in the Fibonacci $`XYZ`$ chain in a magnetic field is discussed in the conclusions.

Free Fermion Results. In preparation for the analysis of more general situations, first we discuss the $`XY`$ or $`\mathrm{\Delta }=0`$ case which already bears some generic features. In this case the model is exactly solvable as it reduces to a bilinear fermionic form. As is well known, after the usual (1-$`d`$) Jordan-Wigner transformation the Hamiltonian can be readily diagonalized by means of a Bogoliubov transformation, leading to

$$H=\underset{k=1}{\overset{N/4}{}}\left\{E_{k,0}(\delta ,\gamma ,h)\left(d_{k,0}^{}d_{k,0}+1\right)+E_{k,1}(\delta ,\gamma ,h)\left(d_{k,1}^{}d_{k,1}+1\right)\right\}$$
(2)
where

$$E_{_k}^{}{}_{(0,1)}{}^{}=\frac{1}{\sqrt{2}}\left[1+\mathrm{cos}(k)+(\delta ^2+\gamma ^2)(1\mathrm{cos}(k))+2h^2\pm f(\delta ,\gamma ,k,h)\right]^{1/2}$$
(3)
and

$$f(\delta ,\gamma ,k,h)=2\left[\delta ^2\gamma ^2(1\mathrm{cos}(k))^2+2h^2\left(1+\mathrm{cos}(k)+\delta ^2(1\mathrm{cos}(k))\right)\right]^{1/2}.$$
(4)
The mean value of the magnetization is thus given by

$$M=\frac{H_{GS}}{h}=\frac{1}{L}\underset{k=1}{\overset{N/4}{}}\left[\frac{E_{k,0}}{h}(\delta ,\gamma ,h)+\frac{E_{k,1}}{h}(\delta ,\gamma ,h)\right]$$
(5)
where GS denotes the ground state, which is the vacuum of Bogoliubov particles’ Fock space. (It should be noticed that the external field does not act as a chemical potential for these particles).

From Eqs. (3) and (4) one easily obtains the critical field values where one of the modes becomes gapless, i.e where the Ising transition occurs

$$h_c=\pm \sqrt{\delta ^2\gamma ^2},$$
(6)
which shows that the Ising transition can only occur for $`\delta >\gamma `$.

The exact magnetization curve shown in Fig. 1 displays all the features we want to discuss (we focus in the region of $`h>0`$). First of all, one notices that there is no actual plateau since the magnetization starts to increase as soon as $`h`$ is turned on, and this is due to the breaking of the $`U(1)`$ symmetry for $`\gamma 0`$. In the region of fields below $`h_c`$ \[ Eq.(6) \], we are in what we call the “pseudo plateau region”, where the slope of the magnetization curve is small and we have dominant density wave correlations. At the critical field one observes the Ising transition where the magnetization should behave as $`MM_c(hh_c)\left(\mathrm{log}|hh_c|1\right)`$. For $`h>h_c`$, the system is in the $`XY`$ phase where the magnetization increases more rapidly. Finally, a second Ising transition is observed before saturation, which occurs at $`h=1`$ independently of the values of $`\delta `$ and $`\gamma `$. For $`\delta <\gamma `$ the first transition does not occur.

Bosonization Analysis. We discuss the results obtained so far within the bosonization approach given below, where we argue that the same picture is valid for arbitrary $`\mathrm{\Delta }`$, provided $`\delta `$ and $`\gamma `$ are suitable renormalized. For small $`\gamma `$ and $`\delta `$ one can study the effects of such perturbations using bosonization, which allows one to include $`\mathrm{\Delta }`$ and $`h`$ exactly through the Bethe Ansatz solution for the Luttinger parameter $`K`$.

The large scale behaviour of the $`XXZ`$ chain can be described by a $`U(1)`$ free boson theory with Hamiltonian

$$H_0=\frac{1}{2}𝑑x\left(vK(_x\stackrel{~}{\varphi })^2+\frac{v}{K}(_x\varphi )^2\right),$$
(7)
which corresponds to the Tomonaga-Luttinger Hamiltonian. The field $`\varphi `$ and its dual $`\stackrel{~}{\varphi }`$ are given by the sum and difference of the light-cone components, respectively. The constant $`K`$ governs the conformal dimensions of the bosonic vertex operators and can be obtained exactly from the Bethe Ansatz solution of the $`XXZ`$ chain (see e.g.  for a detailed summary and references therein). We have $`K=1`$ for the $`SU(2)`$ symmetric case ($`\mathrm{\Delta }=1`$) and is related to the radius $`R`$ of  by $`K^1=2\pi R^2`$. In (7) $`v`$ corresponds to the Fermi velocity of the fundamental excitations of the system. In terms of these fields, the spin operators read

$`S_x^z`$ $``$ $`{\displaystyle \frac{1}{\sqrt{2\pi }}}_x\varphi +a:\mathrm{cos}(2k_Fx+\sqrt{2\pi }\varphi ):+{\displaystyle \frac{M}{2}},`$ (8)
$`S_x^\pm `$ $``$ $`(1)^x:e^{\pm i\sqrt{2\pi }\stackrel{~}{\varphi }}\left(b\mathrm{cos}(2k_Fx+\sqrt{2\pi }\varphi )+c\right):,`$ (9)
where the colons denote normal ordering with respect to the ground state with magnetization $`M`$. The Fermi momentum $`k_F`$ is related to the magnetization of the chain as $`k_F=(1M)\pi /2`$. The effect of an $`XYZ`$ anisotropy and/or the external magnetic field is then to modify the scaling dimensions of the physical fields through $`K`$ and the commensurability properties of the spin operators, as can be seen from (8), (9). The constants $`a`$, $`b`$ and $`c`$ were numerically computed in the case of zero magnetic field  (see also ).

The bosonized Hamiltonian including the perturbations then reads

$$H_{bos}=H_0+\lambda _1𝑑x\mathrm{cos}(\sqrt{2\pi }\varphi )+\lambda _2𝑑x\mathrm{cos}(\sqrt{8\pi }\stackrel{~}{\varphi })$$
(10)
where $`\lambda _1\delta `$, $`\lambda _2\gamma `$.

The scaling dimensions of the perturbations in (10) are $`K/2`$ and $`2/K`$ respectively, which in the $`XX`$ case ($`K=2`$) are both equal to unity. This allows one to understand within this approach the appearance of the Ising transition  since in this self-dual case and for $`h=0`$ one can map the bosonic system into two Majorana fields, and the critical line is given by $`\lambda _1=\pm \lambda _2`$ where one of the masses of the two Majorana fields vanishes rendering criticality. As soon as the external field is turned on the masses of these Majorana fermions change and we then need a bigger value of $`\lambda _1`$ (i.e. $`\delta `$) to find the transition  \[see (6)\].

This can be extended for $`\mathrm{\Delta }>0`$, since its only effect is to modify the scaling dimensions of the two perturbations through $`K(\mathrm{\Delta },h)`$. It is interesting to note, however, that for $`\mathrm{\Delta }>0`$ the cosine of the dual field is less relevant, becoming barely marginal for $`\mathrm{\Delta }=1`$ and $`h=0`$ ($`K(1,0)=1`$), opening the question if one should still expect the Ising transition to occur. One can argue that this is so by analyzing the one-loop renormalization group (RG) equations, where it can be seen that due to the renormalization of $`K`$, which (up to this order in the RG) flows to the $`XX`$ value at large scales for certain initial values of the couplings, the competition between the two cosines can still happen giving rise to the Ising criticality even in this limiting situation. From this qualitative analysis one expects a critical line for $`h=0`$ given by

$$\delta \gamma ^\nu $$
(11)
with $`\nu =\frac{K_{eff}(1K_{eff}/4)}{(K_{eff}1)}`$ and $`K_{eff}`$ is the large scale value of $`K`$. As mentioned above, for $`0\mathrm{\Delta }<1`$, the one loop RG equations indicate that $`K_{eff}=2`$, therefore implying that $`\nu 1`$. A thorough RG analysis should be performed nevertheless, since this large scale value of $`K_{eff}`$ is attained in a strong coupling regime, where the analysis ceases to be valid, but this is out of the scope of the present paper. We provide numerical evidence for the validity of eq. (11) with $`\nu =1`$ below.

Numerical Results. In what follows, we examine this latter conjecture by studying numerically the extreme $`SU(2)`$ case ($`\mathrm{\Delta }=1`$). On the other hand, this enables an independent test of the bosonization scenario within non-perturbative regimes.

First, we computed the gap of Hamiltonian (1) at $`M=0`$ by means of Lanczos diagonalizations  of finite chain lengths $`L`$ with periodic boundary conditions. If $`L/4`$ is an integer it turns out that the ground state is even, i.e. $`\mathrm{exp}\left(i\pi _n\sigma _n^+\sigma _n^{}\right)1`$, otherwise the gap spectrum has to be computed within the odd subspace of $`H`$ . In Fig. 2 we display the results so obtained for $`14L22`$ (see  ). As expected, finite size effects become noticeable on approaching the critical regime encompassed within the pronounced gap minima, whereas convergence towards the thermodynamic limit becomes typically logarithmic . Interestingly, the locations of these minima result however fairly independent of the system size, thus facilitating the estimation of the critical line. For $`L12`$ however, the spectrum gap increases monotonically with $`|\gamma |`$. Therefore, to avoid misleading results we dismiss the data of those small sizes. In turn, this impedes the usage of most of the recursive extrapolation algorithms encountered in the literature , which show their full strength only if the can be iterated several times on a large sequence of sizes. Thus, we content ourselves with a standard (logarithmic) gap extrapolation of the form $`\mathrm{g}\mathrm{gap}(L)+A/L^B`$. Though semi-quantitative \[i.e. $`A=A(\gamma ,\delta )`$, $`B=B(\gamma ,\delta )`$ \], the results shown in the inset of Fig. 2 nevertheless indicate a gap closing at the minima of the finite size data.

In studying numerically the massless regime obtained from such minima however, it should be borne in mind that finite size corrections to the gap of the homogeneous Heisenberg chain ($`\gamma =\delta =0`$), vary slowly as $`\mathrm{ln}(\mathrm{ln}L)/\mathrm{ln}^2L`$, thus affecting the results over a wide range of sizes . Therefore, we are confronted to restricting considerations to regimes where either $`|\gamma |`$ or $`|\delta |`$ are not too small. Fig. 3 exhibits our estimation of the critical line down to $`\delta =0.15`$ above which however, a wide linear regime shows up, in agreement with the bosonization approach. Similar results were observed for $`0<\mathrm{\Delta }<1`$ as well.

Secondly, we turn to the GS magnetization curves of Hamiltonian (1). These were calculated numerically as $`\frac{1}{L}_jGS|\sigma _j^z|GS`$, where in general the parity $`\mathrm{exp}\left(i\pi _n\sigma _n^+\sigma _n^{}\right)`$ of $`|GS`$ comes out to be both field and size dependent. In spite of the existence of a whole massive regime away the critical line, the magnetization grows linearly for small applied fields $`h`$ (a feature which holds also for $`\mathrm{\Delta }=0`$), rather than displaying a plateau at $`M=0`$. This is shown by the inset of Fig. 4(a). Since this latter regime is dominated by a finite gap width (i.e. short range correlations), finite size effects become negligible there. This renders our low field results most reliable, until the first parity change in $`|GS`$ occurs (signaled by an abrupt increase of the magnetization), probably associated to the emergence of gapless modes such as those referred to in the free fermion case. When selecting the Hamiltonian parameters in a massless (zero field) point, the magnetization still grows linearly with $`h`$ though with a substantially bigger amplitude, as is shown by the inset of Fig. 4(b). This results in the removal of all the pseudo-plateaux observed in Fig. 4(a) with a massive $`\gamma \delta `$ point. Also, notice that near the brink of saturation in either case the susceptibility tends to diverge, alike the $`XY`$ situation.

Finally, we address to the predicted Ising behavior of the magnetization curves near the critical fields $`h_c`$ discussed in the previous sections. Upon estimating the former with the fields $`h_c(L)`$ involving the first parity jump in the GS of a finite chain, we obtained a fair logarithmic regime actually applying over more than two decades in $`(hh_c)\left(\mathrm{log}|hh_c|1\right)`$. This is displayed in Fig. 5. The value of the upper line slope results almost independent of the system size and is indicative of the logarithmic regimes entailed both by the free fermion and bosonization approaches.

Conclusions. To summarize, we have studied how the tendencies towards the formation of massive spin excitations (through dimerized $`\delta `$-exchanges), and towards the $`XY`$-ordering (via pairing of $`\gamma `$-interactions), compete with each other. In the bosonization picture, these tendencies manifest themselves in the existence of two (competing) relevant interactions and bring about the Ising transition, which in fact is connected to the same one driven by an external field. Due to the breaking of the $`U(1)`$ symmetry, the latter does not couple to any conserved quantity and therefore the magnetization process gets essentially modified with respect to the $`\gamma =0`$ case, i.e. it ever increases in spite of the presence of massive regimes.

Both the bosonization and the numerical analyses support the picture that the physical mechanism rendering the Ising transition is quite general and valid both in the strong and weak coupling limits. Following the ideas given in  one can conjecture that a similar picture for each plateaux in the XYZ Fibonacci chain holds. The mechanism in that case is the same, where the relevant operator coming from dimerization is replaced for each plateaux by the operator commensurate at the corresponding frequency of the Fibonacci Fourier spectrum. From this one concludes that a similar picture as that found in the XY case  is also valid for generic XYZ Fibonacci chains.

The numerical analysis lent further support to the bosonization results within several non-perturbative situations. We found evidence of a massless line along with the expected Ising like behaviour near critical fields. In turn, this suggests that most basic features of the fully interacting system can be captured by the free fermion picture discussed above. However, the issue as to whether or not all of the GS parity changes induced by the field \[Fig. 4(a) \] correspond to Ising transitions in the thermodynamic limit, remains quite open.

Acknowledgments. It is a pleasure to acknowledge fruitful discussions with A. Honecker and A.A. Nersesyan. The authors acknowledge partial financial support of CONICET and Fundación Antorchas.

$`\mathrm{\#}`$ On leave from Universidad Nacional de La Plata and Universidad Nacional de Lomas de Zamora.


text:
# Criticality and finite size effects in a simple realistic model of stock market

\[

## Abstract

We discuss a simple model based on the Minority Game which reproduces the main stylized facts of anomalous fluctuations in finance. We present the analytic solution of the model in the thermodynamic limit and show that stylized facts arise only close to a line of critical points with non-trivial properties. By a simple argument, we show that, in Minority Games, the emergence of critical fluctuations close to the phase transition is governed by the interplay between the signal to noise ratio and the system size. These results provide a clear and consistent picture of financial markets as critical systems.

\]

Understanding the origin of the anomalous collective fluctuations arising in stock markets poses novel and fascinating challenges in statistical physics. Stock market prices are characterized by anomalous collective fluctuations – known as stylized facts – which are strongly reminiscent of critical phenomena: Prices do not follow a simple random walk process, but rather price increments are fat tailed distributed and their absolute value exhibits long range auto-correlations, called volatility clustering.

The connection with critical phenomena is natural, because financial markets are indeed complex systems of many interacting degrees of freedom — the traders. By means of agent based modeling, it has been realized  that stylized facts are due to the way in which the trading activity of agents interacting in a market “dresses” the fluctuations arising from economic activity – the so-called fundamentals. Ref.  has shown that very simple models based on the Minority Game  can reproduce a quite realistic and rich behavior. Their simplicity makes an analytical approach to these models possible, using tools of statistical physics. Whether Minority Games describe interacting traders is a matter of debate . At any rate, the emergence of anomalous fluctuations in these models, besides providing a scenario for the behavior of real markets, poses questions in statistical physics which deserve interest of their own.

In this Letter, we first introduce the simplest possible Grand Canonical Minority Game (GCMG) which reproduces the main stylized facts, i.e. fat tails and volatility clustering. Then we present the analytic solution of this model in the relevant thermodynamic limit. It shows that the behavior of GCMG, in this limit, exhibits Gaussian fluctuation for all parameter values but on a line of critical points which marks a discontinuous phase transition. For finite size systems, numerical simulations reveal that stylized facts emerge close to the transition line, but they abruptly disappear as the system size increases. Remarkably, the vanishing of stylized facts when the system’s size increases also occurs in a variety of models of financial markets . We present a theory of finite size effects which is fully confirmed by numerical simulations. This allows us to conclude that anomalous fluctuations are properties of the critical point in GCMG. The phase transition is quite unique as it mixes features which are typical of first order phase transitions – as discontinuities and phase coexistence – and of second order phase transitions – such as the divergence of correlation volumes and finite size effects.

In the market described by the Minority Game , agents $`i=1,\mathrm{},N`$ submit a bid $`b_i(t)`$ to the market in every period $`t=1,2,\mathrm{}`$. Agents whose bid has the opposite sign of the total bid $`A(t)=_ib_i(t)`$, win whereas the others lose. The bids of agents depend on the value $`\mu (t)`$ of a public information variable, which is drawn uniformly from the integers $`1,\mathrm{},P`$. In other words, agents have trading strategies which prescribes to agent $`i`$ a bid $`a_i^\mu \pm 1`$ for each information $`\mu `$. Each agent is assigned one such strategy, randomly chosen from the set of $`2^P`$ possible strategies of this type. Agents are adaptive and may decide to refrain from playing if their strategy is not good enough . More precisely, the bids of agents take the form $`b_i(t)=\varphi _i(t)a_i^{\mu (t)}`$ where $`\varphi _i(t)=1`$ or $`0`$ according to whether agent $`i`$ trades or not. In order to assess the performance of their strategy, agents assign scores $`U_i(t)`$ which they update by

$$U_i(t+1)=U_i(t)a_i^{\mu (t)}A(t)ϵ_i.$$
(1)
where

$$A(t)=\underset{i=1}{\overset{N}{}}b_i(t)=\underset{i=1}{\overset{N}{}}\varphi _i(t)a_i^{\mu (t)}.$$
(2)
So if $`a_i^{\mu (t)}A(t)`$ is large enough, i.e., larger than $`ϵ_i`$, the score $`U_i`$ increases. The larger $`U_i`$, the more likely it is that the agent trades ($`\varphi _i=1`$). Here we suppose that

$$\mathrm{Prob}\{\varphi _i(t)=1\}=\frac{1}{1+e^{\mathrm{\Gamma }U_i(t)}}$$
(3)
where $`\mathrm{\Gamma }>0`$ is a constant. A good strategy prescribes bids $`a_i^\mu `$ which tend to coincide with those $`b(t)=\text{sign}A(t)`$ of the minority of agents. The connection with markets goes along the lines of Refs. , which show that $`A(t)`$ is proportional to the difference of price logarithms; here, we take $`\mathrm{log}p(t+1)=\mathrm{log}p(t)+A(t)`$.

The threshold $`ϵ_i`$ in Eq. (1) models the incentives of agents for trading in the market. Some investors may have incentives to trade because they need the market for exchanging goods or assets. This corresponds to $`ϵ_i<0`$. On the contrary, speculators who only trade for profiting of price fluctuations typically have $`ϵ_i>0`$. Of course there may be prudent investors with $`ϵ_i>0`$ or risk-lover speculators with $`ϵ_i<0`$ and a whole range of other type of traders. Here we focus, for simplicity, on the case

$`ϵ_i`$ $`=`$ $`ϵ\text{for }iN_s`$
$`ϵ_i`$ $`=`$ $`\mathrm{}\text{for }N_s<iN`$
The $`N_p=NN_s`$ agents who have $`ϵ_i=\mathrm{}`$ — we call them producers after Refs.  — trade no matter what, whereas the remaining $`N_s`$ — the speculators — trade only if the cumulated performance of their active strategy increases more rapidly than $`ϵt`$.

If the conditional time average $`A|\mu `$ of $`A(t)`$ given $`\mu (t)=\mu `$ is non-zero, then the knowledge of $`\mu (t)`$ allows a statistical prediction of the sign of $`A(t)`$. A measure of predictability is hence given by

$`H_0=\overline{A^2}={\displaystyle \frac{1}{P}}{\displaystyle \underset{\mu =1}{\overset{P}{}}}A|\mu ^2`$

where we introduced the notation $`\overline{(\mathrm{})}`$ for averages over $`\mu `$ ($`\mathrm{}`$ denotes averages on the stationary state). When $`H_0=0`$ the market is unpredictable or informationally efficient. Volatility is instead defined as $`\sigma ^2=\overline{A^2}`$ and it measures market’s fluctuations. A further quantity of interest is the number $`N_{\mathrm{act}}(t)=_i\varphi _i(t)`$ of active speculators in the market.

Exact results can be obtained in the thermodynamic limit, which is defined as the limit $`N_s,N_p,P\mathrm{}`$, keeping constant the reduced number of speculators and producers $`n_s=N_s/P`$, $`n_p=N_p/P`$. In this limit, both $`\sigma ^2`$ and $`H_0`$ diverge with the system size, since $`A(t)\sqrt{N}`$. Hence we shall consider the rescaled quantities $`H_0/P`$ or $`\sigma ^2/P`$. A detailed account of the calculation will be given elsewhere . Here we just discuss the main step and the results. Following Ref. , we derive an Ito stochastic differential equations for the strategy scores $`y_i(\tau )=U_i(t)`$ in the rescaled continuous time $`\tau =t/N`$

$$\frac{\text{d}y_i}{\text{d}\tau }=\overline{a_iA_y}ϵ+\eta _i.$$
(4)
Here $`\eta _i`$ is a zero average Gaussian noise term with

$$\eta _i(\tau )\eta _j(\tau ^{})=\frac{1}{N}\overline{a_ia_jA^2_y}\delta (\tau \tau ^{}).$$
(5)
In Eqs. (4,5) averages $`\mathrm{}_y`$ are taken on the distribution of $`\varphi _i(t)`$ in Eq. (3), which depends on $`y_i(\tau )`$ in a non-linear way: $`\mathrm{Prob}\{\varphi _i(t)=1\}=1/[1+e^{\mathrm{\Gamma }y_i(\tau )}]`$. Hence Eq. (4) is a quite complex system of non-linear equations with a noise strength proportional to the time dependent volatility $`\overline{A^2_y}`$. This feedback will be responsible for the emergence of volatility build-ups.

Following Refs.  we find that the fraction $`\varphi _i`$ of times that agent $`i`$ plays his active strategy in the stationary state is the solution of the minimization of the function

$$H_ϵ=\frac{1}{P}\underset{\mu =1}{\overset{P}{}}\left[\underset{i=1}{\overset{N}{}}\varphi _ia_i^\mu +\underset{i=N_s+1}{\overset{N_s+N_p}{}}a_i^\mu \right]^2+2ϵ\underset{i}{}\varphi _i$$
(6)
with respect to $`\varphi _i`$. Note that for $`ϵ=0`$ this function reduces to the predictability $`H_0`$. For $`ϵ0`$, the solution to this problem, and hence the stationary state, is unique. An exact statistical mechanics description of the solution $`\{\varphi _i\}`$ can be carried out with the replica method, because the replica symmetric ansatz is exact. Furthermore the solution to the Fokker-Planck equation corresponding to Eq. (4) can be well approximated by a factorized ansatz for $`ϵ0`$. This means that the off-diagonal correlations vanish \[$`(\varphi _i\varphi _i)(\varphi _j\varphi _j)=0`$ for $`ij`$\] and, as a consequence, the volatility turns out to be given by $`\sigma ^2=\overline{A^2}=H_0+_{i=1}^{N_s}\varphi _i(1\varphi _i)`$. The solution $`\{\varphi _i\}`$ of the minimization of $`H_ϵ`$ provides a complete description of the model in the limit $`N\mathrm{}`$ for $`ϵ>0`$. In particular the behavior is independent of $`\mathrm{\Gamma }`$.

Fig. 1 shows that all these conclusions are perfectly supported by numerical simulations: With a fixed number $`n_p`$ of producers, as the number $`n_s`$ of speculators increases, the market becomes more and more unpredictable, i.e. $`H_0`$ decreases. At the same time also the volatility $`\sigma ^2`$ decreases. In a market with few speculators ($`n_s<1`$ in Fig. 1), most of the fluctuations in $`A(t)`$ are due to the random choice of $`\mu (t)`$ (i.e. $`\sigma ^2H_0`$) and the number $`n_{\mathrm{act}}`$ of active speculators grows approximately linearly with $`n_s`$.

When $`n_s`$ increases further, the market reaches a point where it is barely predictable. Then, for $`ϵ>0`$ the number of active traders decreases and finally converges to a constant. This means that the market becomes highly selective: Only a negligible fraction of speculators trade ($`\varphi _i(t)=1`$) whereas the majority is inactive ($`\varphi _i(t)=0`$). The volatility $`\sigma ^2`$ also remains constant in this limit.

For $`ϵ<0`$ we see a markedly different behavior: The number of active speculators continues growing with $`n_s`$ even if the market is unpredictable $`H_00`$. The volatility $`\sigma ^2`$ has a minimum and then it increases with $`n_s`$ in a way which depends on $`\mathrm{\Gamma }`$. In other words, $`ϵ=0`$ for $`n_sn_s^{}(n_p)`$ ($`=4.15\mathrm{}`$ for $`n_p=1`$) is the locus of a first order phase transition across which $`N_{\mathrm{act}}`$ and $`\sigma ^2`$ exhibit a discontinuity. This same picture applies to a wider range of GCMG models such as that of Ref. .

Numerical simulations reproduce anomalous fluctuations similar to those of real financial markets close to the phase transition line. As shown in Fig. 2, the distribution of $`A(t)`$ is Gaussian for small enough $`n_s`$, and has fatter and fatter tails as $`n_s`$ increases; the same behavior is seen for decreasing $`ϵ`$. In particular the distribution of $`A(t)`$ shows a power law behavior $`P(|A|>x)x^\beta `$ with an exponent which we estimated as $`\beta 2.8,1.4`$ for $`n_s=20,200`$ respectively and $`ϵ=0.01`$. Note that a realistic value $`\beta 3`$  is obtained for $`n_s=20`$.

This is inconsistent, at first sight, with the theoretical results discussed previously for $`N\mathrm{}`$. Indeed, if the distribution of $`\varphi _i`$ factorizes, $`A(t)`$ is the sum of $`N_s`$ independent contributions and it satisfies the Central Limit Theorem. This implies that for $`ϵ0`$ the variable $`A(t)/\sqrt{N}`$ converges in distribution to a Gaussian variable with zero average and variance $`\sigma ^2/N`$ in the limit $`N\mathrm{}`$. There are no anomalous fluctuations and no stylized facts. Fig. 3 indeed shows that the anomalous fluctuations of Fig. 2 are finite size effects which disappear abruptly as the system size increases (or if $`\mathrm{\Gamma }`$ is small).

In order to understand these finite size effects, we note that volatility clustering arises because the noise strength in Eqs. (4,5) is proportional to the time dependent volatility $`\overline{A^2_y}`$. The noise term is a source of correlated fluctuations because $`\overline{a_ia_jA^2_y}/N1/\sqrt{N}`$ is small but non zero, for $`ij`$. It is reasonable to assume that the dynamics will sustain collective correlated fluctuations in the $`y_i`$ only if the correlated noise is larger than the signal $`\overline{a_iA_y}ϵ`$ which agents receive form the deterministic part of Eq. (4). Time dependent volatility fluctuations would be dissipated by the deterministic dynamics otherwise. A quantitative translation of this insight goes as follows: The noise correlation term is of order $`\overline{a_ia_jA^2_y}/N\sigma ^2/P^{3/2}`$ for $`ij`$. This should be compared to the square of the deterministic term of Eq. (4) $`[\overline{a_iA_y}+ϵ]^2\left[\sqrt{H_0/P}+ϵ\right]^2`$. Rearranging terms, we find that volatility clustering sets in when

$$\frac{H_0}{\sigma ^2}+2ϵ\sqrt{\frac{H_0}{P}}\frac{P}{\sigma ^2}+ϵ^2\frac{P}{\sigma ^2}\frac{K}{\sqrt{P}}$$
(7)
where $`K`$ is a constant. This prediction is remarkably well confirmed by Fig. 4: In the lower panel we plot the two sides of Eq. (7) as a function of $`n_s`$, for different system sizes. The upper panel shows that the volatility $`\sigma ^2/N`$ starts deviating from the analytic result exactly at the crossing point $`n_s^c(P)`$ where Eq. (7) holds true. Furthermore the inset shows that the region $`n_s>n_s^c(P)`$ is described by a novel type of scaling limit. Indeed the curves of Fig. 4 collapse one on top of the other when plotted against $`n_s/n_s^c(P)`$.

The non-linearity of the response of agents is crucial for the onset of volatility time dependence. If $`\mathrm{\Gamma }`$ is small the response becomes smooth and anomalous fluctuations disappear (see Fig. 3). This picture is not affected by the introduction of a finite memory in the learning process of agents, as e.g. in Ref. . In particular the exponents of Fig. 2 do not depend on the memory.

The fact that, in finite systems, stylized facts arise only close to the phase transition is reminiscent of finite size scaling in the theory of critical phenomena: In $`d`$-dimensional Ising model, for example, at temperature $`T=T_c+\epsilon `$ critical fluctuations (e.g. in the magnetization) occur as long as the system size $`N`$ is smaller than the correlation volume $`\epsilon ^{d\nu }`$. But for $`N\epsilon ^{d\nu }`$ the system shows the normal fluctuations of a paramagnet.

Eq. (7) and $`H_0/Pϵ^2`$ imply that the same occurs in the GCMG with $`d\nu =4`$. In other words, the critical window shrinks as $`N^{1/4}`$ when $`N\mathrm{}`$. However, because of the long range nature of the interaction, anomalous fluctuations either concern the whole system or do not affect it at all, as clearly shown in Fig. 3. In the critical region the Gaussian phase coexists probabilistically with a phase characterized by anomalous fluctuations. This and the discontinuous nature of the transition at $`ϵ=0`$, are usually typical of first order phase transitions.

The picture of a phase transition controlled by the signal to noise ratio appears to be universal for Minority Games. Finite size effects close to the phase transition of the standard MG  are indeed explained by the same generic argument: When the signal to noise ratio $`H_0/\sigma ^2`$ is of order $`1/\sqrt{P}`$ self-sustained collective fluctuations arise.

Volatility clustering in real markets is known to be due to wild fluctuations in the volume of trades . Volume is the number of active traders $`N_{\mathrm{act}}+N_p`$ in the GCMG. Hence wild volume fluctuations require correlated collective fluctuations in the behavior of agents which only arise close to criticality. This suggests that real markets operate close to a phase transition. Numerical simulations suggest that exponents vary continuously on the line of critical points. This raises the question of why real markets self-organize close to the critical surface with $`\alpha 3`$.

We conclude that the GCMG exhibits a quite peculiar type of phase transition which mixes properties of continuous and discontinuous transitions. Finite size effects clearly relate the occurrence of stylized facts to the analytic nature of the phase transition. The extension of renormalization group approaches to this system promises to be a quite interesting challenge.


text:
# Effects of electron correlation on the photocurrent in quantum dot infrared photodetectors∗

## Abstract

The effect of electron correlation on the photocurrent of self-assembled InAs/InGaAs quantum dot infrared photo-detector (QDIPs) is studied. It is found that Coulomb interaction and level mixing in the many-body open system lead to double peaks associated with the intra-band transitions involving two lowest levels of the quantum dot. Furthermore, the photocurrent is a nonlinear function of the steady-state carrier density and it displays a plateau due to Coulomb blockade.

This work was supported by a subcontract from the University of Southern California under the MURI program, AFSOR, Contract No. F49620-98-1-0474.

Recently, many efforts have been devoted to the understanding of transport properties in quantum dot photo-detectors (QDIPs). The advantage of the QDIP over quantum well photo-detectors (QWIPs) is that light can be directly coupled to the electrons in the normal incidence geometry due to the effect of QD confinement in directions perpendicular to the growth axis and the dark current is smaller for the same detection wave length considered. Other significant features that are unique to QDs include the Coulomb blockade effect and phonon bottleneck.

Due to the localized nature of electrons in QDs, it is essential to take into account the effects of Coulomb blockade in the analysis of photoresponse of QDIPs, which in general can be ignored in QWIPs. For the non-equilibrium system considered here, it is convenient to use the Keldysh Green function to calculate the transport and optical properties while including the electron correlation. This technique has been used extensively in the study of nonlinear transport properties of quantum systems.

We solve the Anderson Hamiltonian of a two level system coupled with leads in the presence of an electromagnetic radiation with frequency $`\omega `$. We adopt the approach given by Jouho et al and extend it to the present case with asymmetric tunneling rates. We find that the time averaged tunneling current is given by

$$J(t)_j=(\mathrm{\Gamma }_j^R\mathrm{\Gamma }_j^L)\frac{eN_j(t)}{2}e\frac{dϵ}{\pi }[\mathrm{\Gamma }_j^Lf_L\mathrm{\Gamma }_j^Rf_R]ImG_j^r(ϵ,\omega )$$
(1)
where $`N_j(t)`$ is the electron occupation number at QD, $`f_L=f(ϵ\mu _L)`$ and $`f_R=f(ϵ\mu _R)`$ are the Fermi distribution function of the left lead and right lead, respectively. The chemical potential difference between these two leads is related to the applied bias ($`V_a`$) via $`\mu _L\mu _R=eV_a`$. $`\mathrm{\Gamma }_j^L`$ and $`\mathrm{\Gamma }_j^R`$ denote the tunneling rates from the QD to the left and right leads, respectively, for electrons in level $`j`$. $`G_j^r(ϵ,\omega )`$ is the retarded Green function for an electron in level j of the QD. The first term in Eq. (1) for the excited state $`(j=2)`$ provides the photo-induced tunneling current which exists only when $`\mathrm{\Gamma }_2^R\mathrm{\Gamma }_2^L`$, a condition that can occur in a system with asymmetric potential.

Since the incident radiation considered in QDIP application is usually very weak, we can ignore the renormalization of the retarded Green’s function due to electron-photon interaction . Thus, we have (within the Hartree-Fock approximation)

$`G_1^r(ϵ)={\displaystyle \frac{1N_1}{ϵE_1+i\mathrm{\Gamma }_1/2}}+{\displaystyle \frac{N_1}{ϵE_1U_{11}+i\mathrm{\Gamma }_1/2}},`$

$`G_2^r(ϵ)={\displaystyle \frac{1N_1}{ϵE_2+i\mathrm{\Gamma }_2/2}}+{\displaystyle \frac{N_1}{ϵE_2U_{12}+i\mathrm{\Gamma }_2/2}},`$

where $`U_{11}`$ denotes the Coulomb interaction between two electrons in level 1 and $`U_{12}`$ denotes that for one electron in level 1 and the other in level 2. To calculate the steady-state electron occupation number for the QD, we solve the semiconductor Bloch equations for the two-level system coupled to leads and the electromagnetic radiation. We obtain (for $`N_2N_1`$)

$`\mathrm{\Gamma }_1N_1=Im𝒳{\displaystyle \frac{dϵ}{\pi }[\mathrm{\Gamma }_1^Lf_L(ϵ)+\mathrm{\Gamma }_1^Rf_R(ϵ)]ImG_1^r(ϵ)},`$

$`\mathrm{\Gamma }_2N_2=Im𝒳{\displaystyle \frac{dϵ}{\pi }[\mathrm{\Gamma }_2^Lf_L(ϵ)+\mathrm{\Gamma }_2^Rf_R(ϵ)]ImG_2^r(ϵ)},`$

and

$`𝒳(\omega )=2\lambda ^2(N_2N_1)\{{\displaystyle \frac{1N_1}{\omega _r\omega +i\mathrm{\Gamma }/2}}+{\displaystyle \frac{N_1}{\omega _r+U_{12}U_{11}\omega +i\mathrm{\Gamma }/2}}\},`$

where $`\omega _rE_2E_1`$ is the resonant frequency, $`\mathrm{\Gamma }\mathrm{\Gamma }_1+\mathrm{\Gamma }_2`$, and $`\lambda `$ is the momentum matrix element for the inter-level optical transition.

Note that $`𝒳(\omega )`$ consists of two poles corresponding to optical transitions from initial states at $`E_1`$ and ($`E_1+U_{11}`$) to final states at $`E_2`$ and $`(E_2+U_{12})`$, respectively. The two transitions have different relative strengths which depend on the averaged steady -state occupation number in the ground state, $`N_1`$. The first term corresponds to the inter-level transition of a single electron in the QD \[which occurs with a relative probability $`(1N_1)`$\], while the second term corresponds to the inter-level transition of a second electron in the QD (which occurs with a relative probability $`N_1`$) under the influence of the first electron, which remains in the ground state at all times. In the latter case, the Coulomb repulsion between the two electrons give rise to an energy shift $`U_{12}U_{11}`$. The fact that we have a fractional occupancy $`N_1`$ in a single quantum dot is attributed to the level mixing effect (coupling of the QD level to the continuum states in the leads) in the many-body open system.

The theory is applied to a realistic self-assembled quantum dot (SAQD) device. We consider an InAs/InGaAs SAQD system with conical shape. The SAQD is embedded in a slab of InGaAs with finite width, $`W`$. The slab is then placed in contact with heavily doped InGaAs to form an n-i-n structure for infrared detection. Within the effective-mass model, the QD electron is described by the equation

$`[{\displaystyle \frac{1}{2m^{}(\rho ,z)}}+V(\rho ,z)eFz]\psi (\rho ,\varphi ,z)=E\psi (\rho ,\varphi ,z),`$

$`m^{}(\rho ,z)`$ is the position-dependent effective mass, which takes on values of $`m_G^{}=0.067m_e`$ (for GaAs) and $`m_I^{}=0.024m_e`$ (for InAs). The potential $`V(\rho ,z)`$ is equal to $`0`$ in the InGaAs barrier region and $`V_0`$ inside the InAs QD region. The potential in the depletion layers (which separate the slab from the leads) are modelled by an electrostatic potential $`V_d(z)`$

$`V_d(z)=\{\begin{array}{cc}\frac{V_1}{D}(z+\frac{W}{2})\hfill & \text{for }(D+\frac{W}{2})<z<\frac{W}{2}\hfill \\ \frac{V_1}{D}(z\frac{W}{2})\hfill & \text{for }\frac{W}{2}<z<D+\frac{W}{2}.\hfill \end{array}`$

For the purpose of constructing the approximate wave functions, we place the system in a large cylindrical confining box with length $`L`$ and radius $`R`$ ($`R`$ must be much larger than the radius of the cone, $`r_c`$). We adopt $`R=400`$Å, $`D=350`$ Å, $`V_1=0.205eV`$, and $`W=300`$Åfor all calculations. We solve the eigen-functions of the effective-mass Hamiltonian via the Ritz variational method. The wave functions are expanded in a set of basis functions which are chosen to be products of Bessel functions and sine waves

$`\psi _{nlm}(\rho ,\varphi ,z)=J_l(\alpha _n\rho )e^{il\varphi }\mathrm{sin}[k_m(z+{\displaystyle \frac{L}{2}})],`$

where $`k_m=m\pi /L`$, m=1,2,3… Throughout the paper, the origin of z is set at the middle of the confining box. $`J_l`$ is the Bessel function of order $`l`$ ($`l=0,1,2,\mathrm{},`$ etc.) and $`\alpha _nR`$ is the $`n`$-th zero of $`J_l(x)`$. 40 sine functions multiplied by 15 Bessel functions for each angular function ($`l=0`$ or 1) are used to diagonalize the Hamiltonian.

Fig. 1 shows the energy levels of the confined states with $`l=0`$ (solid line) and 1 (dotted line) as functions of height $`h`$ of the QD with base radius fixed at $`R_0=70`$ Å. The other material parameters used here are: wetting layer thickness $`d=3`$ Å, the conduction-band offset $`V_0=0.4`$ eV (this includes the effect of hydrostatic strain due to the lattice mismatch between InAs and In0.2Ga0.8As), and length of the confining box $`L`$ = 600 Å. At least two bound states for each angular function ($`l=0,`$ or 1) are found. For infra-red detector application, we are seeking an intra-band transition (between the ground and first excited state) at an energy around 0.125 eV, which occurs at $`h=50`$ Å, for $`R_0=70`$ Å. The tunneling rates can be calculated numerically via the stabilization method as described in Ref. 3.

We only consider the zero temperature and low bias case, where the chemical potential at the left leads ($`\mu _L`$) is lower than $`E_2`$, so that the average population in the exited state is small. Fig. 2 shows the photocurrent as a function of frequency for various applied voltages: solid line ($`V_a=0.11V`$), dotted line ($`V_a=0.12V`$), and dashed line ($`V_a=0.13V`$). The parameters used to obtain Fig. 2 are $`E_1=139meV`$, $`E_2=14meV`$, $`U_{11}=10.4meV`$, and $`U_{21}=7.2meV`$, which are all calculated based on the effective-mass model. The Fermi level in the source and drain region is assumed to $`E_F=15meV`$. The broadening of the energy level $`E_1`$ including all tunneling processes (dominated by the acoustic-phonon assisted tunneling in this case) is assumed to be $`\mathrm{\Gamma }_1=0.01meV`$. The precise value of $`\mathrm{\Gamma }_1`$ is not important, since photocurrent is not sensitive to $`\mathrm{\Gamma }_1`$. For the excited state, the broadening parameter is given by $`\mathrm{\Gamma }_2=\mathrm{\Gamma }_2^{}+\mathrm{\Gamma }_2^R`$. $`\mathrm{\Gamma }_2^{}`$ is mainly due to radiative and non-radiative recombination from interacting with phonons and defects. The actual value depends on the sample quality and temperature. Here, we assume $`\mathrm{\Gamma }_2^{}=1meV`$. The other contribution due to the direct tunneling is calculated via the stabilization method as described in Ref. . The values are found to be $`\mathrm{\Gamma }_2^R=0.439meV`$, $`0.545meV`$, and $`0.651meV`$ for $`V_a=0.11,0.12`$, and $`0.13V`$, respectively. The spectrum of photocurrent consists of two peaks centered at frequencies $`\omega =E_2E_1`$, and $`\omega =E_2E_1+U_{12}U_{11}`$. The relative strength of these peaks are determined by the average occupation number in the ground state $`(N_1)`$, which is bias-dependent. As shown in the figure, the Coulomb interaction leads to a double-peak photocurrent spectrum with energy separation related to the intra-level and inter-level Coulomb energies ($`U_{11}`$ and $`U_{12}`$).

For the QDIP characteristics, the photocurrent versus applied bias ($`JV`$ curve) is also of interest\[7-9\]. Fig. 3 shows the photocurrent as a function of bias for frequencies at $`\omega _r`$ (dotted line) and $`\omega _r+U_{12}U_{11}`$ (solid line). Using Eq. (2) for the polarization, we can readily understand the behavior of photocurrent. The behavior of the photocurrent is determined by the prefactors $`(1N_1)N_1`$ and $`N_1N_1`$ for the two poles at $`\omega =\omega _r`$ and $`\omega =\omega _r+U_{12}U_{11}`$, respectively. At very low bias, $`N_1`$ is small; thus, the magnitude of the solid line is much weaker than that of the dotted line. As the applied bias increases, the solid line displays a plateau due to the effect of Coulomb blockade on $`N_1`$. When the applied bias overcomes the charging effect, $`N_10.5`$ and the solid line becomes almost identical to the dotted line.

In this study, we have used tunneling carriers as the source for photocurrent, in contrast to the captured carriers typically used in QWIPs. Due to the phonon bottleneck effect, it is predicted that the capture rate of the electron by the QD will be low. This could reduce the performance of QDIPs which use captured carriers as a photocurrent source. Using tunneling carriers as the photocurrent source will not have this drawback, and in this case we find that the effect of electron correlation leads to a double-peak spectrum and the photocurrent is a highly nonlinear function of the carrier density in the ground state $`N_1`$. Both of these effects must be taken into account in the analysis of photoresponse of QDIPs.

Figure Captions

Fig. 1. Energies of the bound states of a conical InAs/$`In_{0.2}Ga_{0.8}`$As QD as functions of height $`h`$ of the QD with base radius fixed at $`R_0=70`$ Å. Solid lines: ($`l=0`$). Dotted lines: ($`l=1`$).

Fig. 2. Photocurrent as a function of frequency for different applied voltages: $`V_a=0.11V`$ (solid line), $`V_a=0.12V`$ (dotted line) and $`V_a=0.13V`$ (dashed line).

Fig. 3. Photocurrent as a function of bias for incident frequencies at $`\omega =E_2E_1`$ (dotted line) and $`\omega =E_2E_1+U_{12}U_{11}`$ (solid line).


text:
# Electrodynamics of a Coulomb Glass in n-type Silicon

## Abstract

Optical measurements of the real and imaginary frequency dependent conductivity of uncompensated n-type silicon are reported. The experiments are done in the quantum limit, $`\mathrm{}\omega >k_BT`$, across a broad doping range on the insulating side of the Metal-Insulator transition (MIT). The observed low energy linear frequency dependence shows characteristics consistent with theories of a Coulomb glass, but discrepancies exist in the relative magnitudes of the real and imaginary components. At higher energies we observe a crossover to a quadratic frequency dependence that is sharper than expected over the entire dopant range. The concentration dependence gives evidence that the Coulomb interaction energy is the relevant energy scale that determines this crossover.

Doped semiconductors close to the MIT are systems ideally suited for the study of strong electron-electron interactions. Below a critical doping the material is an insulator, i.e. at zero temperature there is no DC conductivity. The charge carriers are localized in the Anderson sense, and screening is reduced from the metallic regime.

Anderson coined the term Fermi glass to describe a non-interacting disordered insulating system whose universal properties, independent of system specific details, are determined by Fermi statistics alone . Mott’s original treatment of such a system , which has a finite density of states at the Fermi level, but where localization makes an insulator, did not address the issue of electron-electron interactions as pointed out by Pollak , and then Efros and Sklovskii (ES) . A Fermi glass that includes interactions between localized electrons has been termed a Coulomb glass. Among other things, such a system is typified by a depletion in the single particle density of states around the Fermi level which was termed by ES the Coulomb gap, $`\mathrm{\Delta }`$ .

By taking into account the mean Coulomb interaction between two sites forming a resonant pair $`U(r_\omega )=e^2/\epsilon _1r_\omega `$, where $`r_\omega =\xi [ln(2I_0/\mathrm{}\omega )]`$ is the most probable hop distance between pairs, $`I_0`$ is the hopping attempt rate and $`\epsilon _1`$ is the dielectric constant, ES derived the real part of the ac conductivity of a Coulomb glass to be:

$$\sigma _1=\beta e^2g_0^2\xi ^5\omega [ln(2I_0/\mathrm{}\omega )]^4[\mathrm{}\omega +U(r_w)].$$
(1)
Here $`\beta `$ is a constant of order one, $`g_0`$ is the non-interacting density of states and $`\xi `$ is the localization length. This formula takes on a different frequency dependence in two limits. When the photon energy, $`\mathrm{}\omega >U(r_\omega )`$, one recovers the same quadratic frequency dependence that Mott derived for a non-interacting Fermi glass. Here the Coulomb glass is indistinguishable from the Fermi glass in so far as the high frequency limit of the conductivity is concerned. In the opposite limit, $`\mathrm{}\omega <U(r_\omega )`$, the conductivity of a Coulomb glass will show an approximately linear dependence on frequency, plus logarithmic corrections. The imaginary component of the complex conductivity, $`\sigma _2`$, as predicted by Efros, should be identical to $`\sigma _1`$ up to a logarithmic factor. We should note that Eq. (1) was derived for the case where $`\mathrm{}\omega >\mathrm{\Delta }`$, the Coulomb gap width. However a linear dependence (albeit with additional logarithmic corrections) and an eventual crossover to Mott’s non-interacting quadratic law is still expected even for the case where $`\mathrm{}\omega <\mathrm{\Delta }`$.

Despite its fundamental importance, there have been very few AC conductivity studies that have been done in the quantum limit $`\mathrm{}\omega >k_BT`$, yet at low enough frequencies to probe the relevant energy ranges of Eq. (1). Only very recently, have measurements been attempted that address these issues. M. Lee et al. found that for concentrations close to the MIT the expected linear to quadratic crossover occurs, but is much sharper than predicted . Other recent work was done on amorphous NbSi, where the general systematics of a linear in frequency Coulomb glass response was found, but a crossover to the non-interacting quadratic dependence was not reached .

In this Letter we report the first measurements of the frequency dependent real and imaginary conductivity on the insulating side of the MIT for a doped, crystalline semiconductor. We find a concentration and frequency dependence of the components consistent with predictions of a Coulomb glass, but a discrepancy arises in the predicted ratio of their relative magnitudes. At higher frequencies a crossover to quadratic Mott-like behavior is observed and due to having measured across a broad dopant range, the concentration dependence of the crossover energy is shown to be consistent with the Coulomb interaction energy $`U`$ and not the Coulomb gap width $`\mathrm{\Delta }`$.

Nominally uncompensated n-type silicon samples were obtained from Recticon Enterprises Inc. A boule of silicon was grown using the Czochralski method and sliced into 1 mm thick slabs. Room temperature resistivity was measured using an ADE 6035 resistivity gauge, and the concentration determined using the Thurber scale. A number of samples were etched with a $`4\%HF+96\%HNO_3`$ solution; this resulted in no difference in the results of frequency dependent conductivity measurements. The Si:P samples discussed in this Letter span a range from 39% to 69%, stated as $`\frac{x}{x_c}`$, a percentage ratio of the sample’s dopant concentration to the critical concentration at the MIT.

At 35 and 60 GHz the real part of the conductivity was evaluated from the measured loss of highly sensitive resonant cavities via the perturbation method. The technique and analysis is well established. In the millimeter spectral range, 80 GHz to 1 THz, backward wave oscillators were employed as coherent sources in a transmission configuration. Fabry-Perot like resonances in the transmission were analyzed uniquely determining both components of the complex conductivity. Cavity measurements were performed down to 1.2 K and in the millimeter spectral range down to 2.8 K. With these base temperatures the quantum limit $`\mathrm{}\omega >k_BT`$ ($`1K20GHz`$) of the system was being investigated.

The real part of the frequency dependent conductivity, $`\sigma _1`$, at $`T=2.8K`$ for two samples is shown in Fig. 1. This data, representative of all the samples, shows an approximately linear dependence at low frequencies and then a sharp crossover to an approximately quadratic behavior at higher frequencies. This is the qualitatively expected behavior from Eq. (1). However, as seen by the overlayed fits, Eq. (1) provides only a rough guide. The solid lines are linear and quadratic fits to the low frequency and high frequency data respectively. Individually these functions fit the data well, however the full crossover function does not. The dotted line is a fit using the same method as Ref. , namely forcing the linear portion to pass through the origin as well as the low frequency data and leaving the pre-factor of the quadratic term as a free variable. The crossover between linear and quadratic portions is much more abrupt than the ES function predicts and is observed over our entire doping range as was observed previously in an analogous system, Si:B, for samples closer to the MIT .

Having shown that the frequency dependence is qualitatively consistent with a Coulomb glass, we focus on the low frequency regime. Both the real part and the magnitude of the imaginary part, $`|\sigma _2|`$. of the complex conductivity are plotted as a function of frequency in Fig. 2. The dotted lines are linear fits and show the relative increase in the magnitude of both components as the MIT is approached.

A complex conductivity obeying a power law can be expressed in a simple Kramers-Kronig compatible form, $`\sigma _1(\omega )+i\sigma _2(\omega )=A(i\omega )^\alpha `$, in order to have an independent means for determining the exponent. To determine the power $`\alpha `$ one can utilize the ratio of $`|\sigma _2|`$ versus $`\sigma _1`$ (with the frequency as a variable). The power $`\alpha `$ is given by,

$$\alpha =\frac{2}{\pi }tan^1\left(\frac{|\sigma _2|}{\sigma _1}\right).$$
(2)
The uppermost panel in Figure 3 shows the ratio of the imaginary to the real part of the conductivity for Si:P. NbSi is an amorphous insulating glass which is characterized by a vastly higher density of states (DOS) at $`E_F`$ . Disorder should give a somewhat similar functional dependence for the AC conductivity, but with very different energy scales. Data from it is included as another example of AC conductivity in a disordered system and for comparison purposes, as Si:P has a very small DOS in its impurity band. First we note that this ratio for Si:P remains large and approximately constant across our range of dopant concentrations. From theory, one expects $`|\sigma _1|`$ to be approximately equal to $`|\sigma _2|`$ to within a factor of 2-5 (with a reasonable estimate for $`I_0`$) as predicted by Efros . Applied to Si:P, Eq. (1) in the $`\mathrm{}\omega <U(r_\omega )`$ limit correctly predicts a linear correspondence between $`\sigma _1`$ and $`|\sigma _2|`$, but the theory incorrectly predicts the measured proportionality by at least a factor of thirty. The proportionality is closer for NbSi but has a dependence on the doping concentration. This may be related to entering the quantum critical (QC) regime as discussed below. Here we have used the susceptibility due to the interacting electrons themselves $`4\pi \chi =\epsilon _1\epsilon _{Si}`$, namely the full measured dielectric constant, $`\epsilon _1`$, minus the static background dielectric constant of the host silicon, $`\epsilon _{Si}`$, to determine the magnitude of the imaginary component of the conductivity shown in Figure 3. The magnitude of the ratio becomes even larger and the discrepancy greater if the full dielectric constant $`\epsilon _1`$ is considered instead.

The middle panel in Figure 3 shows the power $`\alpha `$ as determined by Eq. (2). The values for Si:P are approximately equal to, but slightly less than one, consistent with Fig. 2. This indicates that the prefactor of the real and imaginary components of the complex conductivity have the same concentration dependence. When approaching the MIT, the frequency dependence is expected to cross over to the QC behavior, i.e. $`\sigma _1\omega ^{1/2}`$ in NbSi, when the localization length $`\xi `$ becomes comparable to the characteristic frequency dependent length scale, e.g. the dephasing length, $`\mathrm{}_\omega `$. The crossover is not a phase transition and need not be sharp, therefore looking at a fixed window of frequencies, a broad, smooth crossover from $`\omega \omega ^{1/2}`$ would show an averaged power of the frequency dependence similar to that measured for NbSi shown in the middle panel of Figure 3. The fact that we see an $`\alpha 1`$ across our whole doping range in Si:P, but an $`\alpha `$ that approaches 0.5 in NbSi indicates that the critical regime in Si:P is much narrower. Simple dimensional arguments give a result similar to the non-interacting case that the crossover should be inversely proportional to the dopant DOS. The vastly smaller (a factor of $`10^3`$) dopant density in Si:P relative to NbSi is consistent with a narrower QC regime in Si:P as compared to NbSi.

Because our data spans a large range of concentrations, the doping dependence of the crossover energy scale from linear to quadratic can be analyzed to see whether its dependence is consistent with the functional form of other energy scales, e.g. the Coulomb interaction energy $`U`$ or the Coulomb gap width $`\mathrm{\Delta }`$. Fig. 4 shows data of our experimentally determined crossover energy scale plotted versus the commonly used functional form, $`1x/x_c`$. M. Lee et al. performed AC conductivity measurements on a pair of Si:B samples , and in said measurements, a sharp crossover analogous to our own observations was seen. Our data extends the observed range of the sharp crossover deep into the insulating regime as M. Lee et al.’s original measurements were closer to the critical concentration.

Recall that the Coulomb interaction energy between two sites forming a resonant pair is $`U(r_\omega )=e^2/\epsilon _1r_\omega `$, which is dependent on concentration via the dielectric constant and the localization length dependent most probable hop distance. The Coulomb gap width as determined by ES  is $`\mathrm{\Delta }=e^3g_0^{1/2}/\epsilon _1^{3/2}`$. Although M. Lee et al. postulated that a smaller Coulomb gap governed the crossover, it is not unreasonable to expect that this gap scales with the ES single particle one as both are presumably caused by the long range Coulomb interaction. The magnitude of the gap width is dependent on concentration through the dielectric constant and through the DOS term. Thus it is clear that a separate measurement of the full dielectric constant, $`\epsilon _1`$ as a function of concentration, as we have done, need be performed in order to analyze the crossover energy scale. Although there is no theoretical expectation for $`\omega _c`$ or $`\mathrm{\Delta }`$ to scale as a power law over the whole doping range, we may still parameterize these quantities over our doping range as power laws. We find that the crossover frequency $`\omega _c`$ is proportional to $`(1x/x_c)^{1.65}`$, and the full dielectric constant $`\epsilon _1(1x/x_c)^{0.4}`$. Using these values, it seems highly unlikely that the $`\mathrm{}\omega _c=\mathrm{\Delta }`$, as the required DOS’s concentration dependence for this to be valid seems improbable.

We do find strong support that the crossover energy is determined by the Coulomb interaction energy though. By setting the measured crossover energy scale equal to the Coulomb interaction energy we are able to determine both the magnitude of the localization length and the exponent with which it diverges as a function of concentration. Using an appropriate pre-factor for the overlap integral , $`I_0=10^{13}s^1`$ in the most probable hop distance term, $`r_\omega `$, we find a localization length dependence as plotted in the bottom panel of Fig. 3. We find that the localization length exponent is close to unity, the value originally predicted by McMillan in his scaling theory of the MIT , and the magnitude of the localization length is reasonable. Both these results strongly point towards the Coulomb interaction energy as being the energy scale at which the observed frequency dependent crossover from ES to Mott-like hopping conduction occurs.

In summary, we have observed behavior consistent with a Coulomb glass across our entire range of doping concentrations in Si:P for the frequency and concentration dependence of the real and imaginary components of the conductivity. It is consistent with theoretical predictions except in predicting the relative magnitudes of $`\sigma _1`$ and $`\sigma _2`$ and the sharper than expected crossover from linear to quadratic frequency dependence. From the expected form of $`U(x)`$ and the reasonable behavior of the extracted localization length, it seems likely that the crossover is governed by the Coulomb interaction strength of a resonant pair and not the Coulomb gap.

We would also like to thank Steve Kivelson as well as Boris Shklovskii for useful conversations. This research was supported by the National Science Foundation grant DMR-0102405.