text:
# On the Direct Detectability of the Cosmic Dark Ages: 21-cm Emission from Minihalos

## 1 Introduction

No direct observation of the universe during the period between the recombination epoch at redshift $`z10^3`$ and the reionization epoch at $`z6`$ has yet been reported. While a number of suggestions for the future detection of the reionization epoch, itself, have been made, this period prior to the formation of the first stars and quasars – the cosmic “dark ages” (e.g. Rees 1999) – has been more elusive. Standard Big Bang cosmology in the CDM model predicts that nonlinear baryonic structure first emerges during this period, with virialized halos of dark and baryonic matter which span a range of masses from less than $`10^4M_{}`$ to about $`10^8M_{}`$ which are filled with neutral hydrogen atoms. The atomic density $`n_H`$ and kinetic temperature $`T_K`$ of this gas are high enough that collisions populate the hyperfine levels of the ground state of these atoms in a ratio close to that of their statistical weights (3:1), with a spin temperature $`T_S`$ that greatly exceeds the excitation temperature $`T_{}=0.0681\mathrm{K}`$. Since, as we shall show, $`T_S>T_{\mathrm{CMB}}`$, the temperature of the Cosmic Microwave Background (CMB), as well, for the majority of the halos, these “minihalos” can be a detectable source of redshifted 21-cm line emission. The direct detection of minihalos at such high redshift would be an unprecedented measure of the density fluctuations in the baryons and of the total matter power spectrum at small scales, which will not be probed by other methods yet discussed (e.g. CMB anisotropy).

The possibility of 21-cm line emission or absorption by neutral H at high redshift has been considered before (Hogan & Rees, 1979; Scott & Rees, 1990; Subramanian & Padmanabhan, 1993; Kumar, Padmanabhan, & Subramanian, 1995; Bagla, Nath, & Padmanabhan, 1997; Madau, Meiksen & Rees, 1997; Shaver, et al., 1999; Tozzi et al., 2000). Prior to the release of radiation by nonlinear baryonic structures which condense out of the background universe, the spin temperature of $`\mathrm{H}\mathrm{I}`$ in the diffuse, uncollapsed gas in the intergalactic medium (IGM) is coupled to the CMB, so that $`T_S=T_{\mathrm{CMB}}`$ and neither emission nor absorption in the 21-cm line is possible. Recently, attention has focused on the possibility that radiation from early stars and quasars might decouple $`T_S`$ from $`T_{\mathrm{CMB}}`$ by “Ly$`\alpha `$ pumping” – resonant scattering in the H Ly$`\alpha `$ transition followed by decay of the upper state $`n=2`$ to the ground state $`n=1`$ into one or the other of the hyperfine levels (Madau et al. 1997; Tozzi et al. 2000). This mechanism, it has been suggested, will operate on $`\mathrm{H}\mathrm{I}`$ in the diffuse, uncollapsed IGM during reionization, first to make $`T_S<T_{\mathrm{CMB}}`$, so that the 21-cm transition can be seen in absorption against the CMB, until the same Ly$`\alpha `$ scattering heats the gas shortly thereafter and makes $`T_S>T_{\mathrm{CMB}}`$, thereby causing 21-cm emission in excess of the CMB, before reionization finally destroys the $`\mathrm{H}\mathrm{I}`$. In what follows, however, we show that a substantial fraction of the baryons in the universe may already have condensed out of the diffuse IGM into virialized minihalos, prior to and during reionization. Under these conditions, collisional excitation alone is sufficient to decouple $`T_S`$ from $`T_{\mathrm{CMB}}`$ and cause 21-cm emission in excess of the CMB, thereby providing a signature of the cosmic “dark ages” and of their retreat during reionization.

## 2 21-cm Emission from Individual Minihalos

The 21-cm emission from a single halo depends upon its internal atomic density, temperature, and velocity structure. We model each CDM minihalo here as a nonsingular, truncated isothermal sphere (“TIS”) of dark matter and baryons in virial and hydrostatic equilibrium, in good agreement with the results of gas and N-body simulations from realistic initial conditions (Shapiro, Iliev & Raga, 1999; Iliev & Shapiro, 2001, 2002). This model uniquely specifies the internal structure of each halo (e.g. total and core sizes, $`r_t`$ and $`r_0`$, central total mass density $`\rho _0`$, dark matter velocity dispersion $`\sigma _V=(4\pi \rho _0r_0^2)^{1/2}`$, and gas temperature $`T_K=\mu m_p\sigma _V/k_B`$, where $`\mu `$ is the mean molecular weight), for a given background cosmology as functions of two parameters, the total mass $`M`$ and collapse redshift $`z_{\mathrm{coll}}`$.

The minihalos which contribute significantly to the 21-cm emission span a mass range from $`M_{\mathrm{min}}`$ to $`M_{\mathrm{max}}`$ which varies with redshift. $`M_{\mathrm{min}}`$ is close to the Jeans mass of the uncollapsed IGM prior to reionization, $`M_J=5.7\times 10^3\left(\mathrm{\Omega }_0h^2/0.15\right)^{1/2}\left(\mathrm{\Omega }_bh^2/0.02\right)^{3/5}\times \left[(1+z)/10\right]^{3/2}M_{}`$, while $`M_{\mathrm{max}}=3.95\times 10^7\times (\mathrm{\Omega }_0h^2/0.15)^{1/2}[(1+z)/10]^{3/2}`$ is the mass for which $`T_{\mathrm{vir}}=10^4`$ K according to the TIS model (Iliev & Shapiro 2001) (since halos with $`T_{\mathrm{vir}}10^4\mathrm{K}`$ are largely collisionally ionized). Halos with $`T_{\mathrm{vir}}10^4\mathrm{K}`$ may have radiatively cooled gas inside them which would add to the signal we compute, but such gas is expected to lead to the formation of internal sources of ionizing radiation which will largely offset the effect. Since these additional effects are highly uncertain and are related to the onset of radiative feedback and reionization which we are neglecting in these calculations, we will not consider the role of higher temperature halos further.

The flux per unit frequency, $`𝔉_\nu (dF/d\nu )_{\mathrm{rec}}`$, received at redshift $`z=0`$ at frequency $`\nu _{\mathrm{rec}}`$ from a minihalo at redshift $`z`$ which emits at frequency $`\nu _{\mathrm{em}}=\nu _{\mathrm{rec}}(1+z)`$ is expressed in terms of the brightness temperature $`T_{b,\mathrm{em}}=T_{b,\mathrm{rec}}(1+z)`$ according to

$$𝔉_{\nu _{\mathrm{rec}}}=\frac{2\nu _{\mathrm{rec}}^2}{c^2}k_BT_{b,\mathrm{rec}}(\mathrm{\Delta }\mathrm{\Omega })_{\mathrm{halo}},$$
(1)
where $`(\mathrm{\Delta }\mathrm{\Omega })_{\mathrm{halo}}=\pi r_t^2/D_A^2=\pi (\mathrm{\Delta }\theta _{\mathrm{halo}}/2)^2`$ is the solid angle subtended by the minihalo, and $`D_A`$ is the angular diameter distance. The brightness temperature $`T_{b,\mathrm{em}}`$ is determined by solving the equation of radiative transfer to derive the brightness profile of the minihalo and integrating this profile over the projected surface area, as follows.

The brightness temperature along a line of sight thru a minihalo at projected distance $`r`$ from the center obeys the equation

$$T_b(r)=T_{\mathrm{CMB}}e^{\tau (r)}+_0^{\tau (r)}T_Se^\tau ^{}𝑑\tau ^{},$$
(2)
where quantities are defined in the comoving frame of the minihalo, frequency $`\nu `$ refers here and henceforth to $`\nu _{\mathrm{em}}`$, $`\tau (r)`$ is the total optical depth thru the halo, and the effective absorption coefficient $`\kappa _\nu `$ is given when $`T_{}T_S`$ by:

$$\kappa _\nu =\frac{3c^2A_{10}n_{\mathrm{HI}}}{32\pi \nu ^2}f(\nu )\frac{T_{}}{T_S}$$
(3)
(Field, 1958), where $`A_{10}=2.85\times 10^{15}\mathrm{s}^1`$ is the Einstein $`A`$-coefficient for the 21-cm transition and $`f(\nu )`$ is the normalized line profile. The spin temperature $`T_S`$ is determined by the balance between collisional and radiative excitation and de-excitation by atoms and electrons and by CMB and Ly$`\alpha `$ photons, respectively, according to

$$T_S=\frac{T_{CMB}+y_cT_K+y_\alpha T_\alpha }{1+y_c+y_\alpha },$$
(4)
where $`T_\alpha `$ is the color temperature of the Ly$`\alpha `$ photons, and $`y_\alpha `$ and $`y_c`$ are radiative and collisional excitation efficiencies, respectively (Purcell & Field, 1956; Field, 1958, 1959). The efficiency $`y_c`$ includes contributions from $`\mathrm{H}^0\mathrm{H}^0`$ collisions, $`y_H`$, and from $`\mathrm{e}^{}\mathrm{H}^0`$ collisions, $`y_e`$. Prior to the reionization epoch, Ly$`\alpha `$ pumping is unimportant and collisional excitation alone must compete with excitation by the CMB. This is only possible for gas that is highly nonlinear and sufficiently hot. Such conditions are achieved only inside virialized halos. As shown in Figure 1, the optical depth of an individual halo is not negligible, particularly for smaller-mass halos (due to their lower $`T_S`$). Since $`T_S`$ varies with radial position inside the halo, as a result of its significant central concentration, we must integrate equation (2) numerically. The face-averaged $`T_b`$ of this single halo is given by $`T_b_{\mathrm{halo}}(T_b(r)𝑑A)/A`$, where $`A(M,z)`$ is the geometric cross-section of a halo of mass $`M`$ and collapse redshift $`z`$. The observed flux from an individual halo is then expressed with respect to the CMB by the differential antenna temperature $`\delta T_b[T_b_{\mathrm{halo}}T_{\mathrm{CMB}}(z)]/(1+z)`$.

The line-integrated flux $`F(M,z)`$ received from this minihalo is equal to the flux calculated for $`\nu =\nu _0`$ multiplied by a redshifted effective line-width $`\mathrm{\Delta }\nu _{\mathrm{eff}}(z)`$, defined by $`\mathrm{\Delta }\nu _{\mathrm{eff}}(z)(𝔉𝑑\nu )/𝔉_{\nu _0}`$. For an optically thin minihalo, $`\mathrm{\Delta }\nu _{\mathrm{eff}}`$ reduces to $`\mathrm{\Delta }\nu _{\mathrm{eff}}(z)=[f(\nu _0)(1+z)]^1`$. In that case, for a thermal-Doppler-broadened line profile, $`\mathrm{\Delta }\nu _{\mathrm{eff}}(z)=[(2\pi \mu )^{1/2}\nu _0\sigma _V/c](1+z)^1`$. We have checked that this approximation is adequate even for the optically thicker halos at the small-mass end of the mass function. The differential line-integrated flux $`\delta F(M,z)`$ is given by replacing $`T_{b,\mathrm{rec}}`$ in equation (1) by $`\delta T_b`$ and integrating over frequency as described above.

Our results for individual minihalos are summarized in Figure On the Direct Detectability of the Cosmic Dark Ages: 21-cm Emission from Minihalos. Line profiles of different minihalos along the same line of sight should not typically overlap. The proper mean free path $`\lambda _{\mathrm{mfp}}=n_{\mathrm{halo}}\sigma _{\mathrm{halo}}^1`$ for photons to encounter minihalos in $`\mathrm{\Lambda }`$CDM is 160 kpc at $`z=9`$ (Shapiro 2001), corresponding to a frequency separation, $`\mathrm{\Delta }\nu _{\mathrm{sep}}\nu _0H(z)\lambda _{\mathrm{mfp}}/[c(1+z)]0.1\mathrm{MHz}\mathrm{\Delta }\nu _{\mathrm{eff}}10\mathrm{kHz}`$.

These results predict a “21-cm forest” of minihalo emission lines. At $`z=9`$, for example, there are about 160 minihalo lines per unit redshift along a typical line of sight in an untilted $`\mathrm{\Lambda }`$CDM universe (Shapiro 2001). Detecting the stronger lines would require sub-arcsecond spatial resolution, $`1\mathrm{kHz}`$ frequency resolution, and $`\mathrm{nJy}`$ sensitivity. SKA is expected to have sufficient resolution for such observation, but probably not sufficient sensitivity.

## 3 21-cm Radiation Background from Minihalos

The average differential flux per unit frequency relative to that of the CMB from all the minihalos observed within a given beam of angular size $`\mathrm{\Delta }\theta _{\mathrm{beam}}`$ and frequency bin $`\mathrm{\Delta }\nu _{\mathrm{obs}}`$ is:

$$\overline{\delta 𝔉}_\nu (z)=\frac{\mathrm{\Delta }z(\mathrm{\Delta }\mathrm{\Omega })_{\mathrm{beam}}}{\mathrm{\Delta }\nu _{\mathrm{obs}}}\frac{d^2V(z)}{dzd\mathrm{\Omega }}_{M_{\mathrm{min}}}^{M_{\mathrm{max}}}\delta F\frac{dn}{dM}𝑑M,$$
(5)
where $`d^2V(z)/dzd\mathrm{\Omega }`$ is the comoving volume per unit redshift per unit solid angle, the solid angle $`(\mathrm{\Delta }\mathrm{\Omega })_{\mathrm{beam}}=\pi (\mathrm{\Delta }\theta _{\mathrm{beam}}/2)^2`$, and $`\mathrm{\Delta }\nu _{\mathrm{obs}}/\mathrm{\Delta }z=\nu _0/(1+z)^2`$. We calculate the comoving density of halos at different redshifts using the Press-Schechter (PS) approximation for the halo mass function $`dn/dM`$. If we define the beam-averaged “effective” differential antenna temperature $`\overline{\delta T}_b`$ using $`\overline{\delta 𝔉}_\nu =2\nu ^2k_B\overline{\delta T_b}(\mathrm{\Delta }\mathrm{\Omega })_{\mathrm{beam}}/c^2`$, then

$$\overline{\delta T}_b=\frac{c(1+z)^4}{\nu _0H(z)}_{M_{\mathrm{min}}}^{M_{\mathrm{max}}}\mathrm{\Delta }\nu _{\mathrm{eff}}\delta T_{b,\nu _0}A\frac{dn}{dM}𝑑M.$$
(6)
We consider the currently-favored, flat CDM model with cosmological constant (“$`\mathrm{\Lambda }`$CDM”, $`\mathrm{\Omega }_0=0.3`$, $`\lambda _0=0.7`$, COBE-normalized, $`\mathrm{\Omega }_bh^2=0.02`$, $`h=0.7`$), for three values of the primordial power spectrum index $`n_p=0.9`$, 1, and 1.1, using the primordial power spectrum of Eisenstein & Hu (1999).

Results for $`\overline{\delta 𝔉}_\nu `$ and $`\overline{\delta T}_b`$ are plotted in Figure On the Direct Detectability of the Cosmic Dark Ages: 21-cm Emission from Minihalos. In principle, the variation of $`\overline{\delta T_b}`$ with observed frequency implied by the redshift variations in Figure On the Direct Detectability of the Cosmic Dark Ages: 21-cm Emission from Minihalos should permit a discrimination between the 21-cm emission from minihalos and the CMB and other backgrounds, due to their very different frequency dependences. However, the average differential brightness temperature of this minihalo background is very low and its evolution is fairly smooth, so such measurement may be difficult in practice with currently planned instruments like LOFAR and SKA. The angular fluctuations in this emission, on the other hand, should be much easier to detect, as discussed in the next section.

## 4 Angular Fluctuations in the 21-cm Emission Background

The amplitude of $`q`$-$`\sigma `$ angular fluctuations (i.e. $`q`$ times the rms value) in the differential antenna temperature is given in the linear regime by

$$\frac{\delta T_b^2^{1/2}}{\overline{\delta T_b}}=qb(z)\sigma _p,$$
(7)
where $`\sigma _p`$ is the rms mass fluctuation at redshift $`z`$ in a randomly placed cylinder which corresponds to the observational volume defined by the detector angular beam size, $`\mathrm{\Delta }\theta _{\mathrm{beam}}`$, and frequency bandwidth, $`\mathrm{\Delta }\nu _{\mathrm{obs}}`$, and $`b(z)`$ is the bias factor which accounts for the clustering of rare density peaks relative to the mass. We assume $`b(z)`$ is the flux-weighted average over the mass function of $`b(M,z)=1+(\nu _\mathrm{h}^21)/\delta _c`$, the linear bias factor, where $`\nu _\mathrm{h}=\delta _c/\sigma (M)`$, $`\delta _c`$ is the value of the linearly extrapolated value of overdensity $`\delta \rho /\rho `$ corresponding to the epoch when a top-hat collapse reaches infinite density, and $`\sigma (M)`$ is the the standard deviation of the density contrast filtered on mass scale $`M`$ (e.g. Mo & White 1996). For a cylinder of comoving radius $`R=\mathrm{\Delta }\theta _{\mathrm{beam}}(1+z)D_A(z)/2`$, and length $`L(1+z)cH(z)^1(\mathrm{\Delta }\nu /\nu )_{\mathrm{obs}}`$, we have :

$`\sigma _p^2`$ $`=`$ $`{\displaystyle \frac{8D^2(z)}{\pi ^2R^2L^2}}{\displaystyle _0^{\mathrm{}}}𝑑k{\displaystyle _0^1}𝑑x{\displaystyle \frac{\mathrm{sin}^2(kLx/2)J_1^2[kR(1x^2)^{1/2}]}{x^2(1x^2)}}`$ (8)
$`\times (1+fx^2)^2{\displaystyle \frac{P(k)}{k^2}}`$
(Tozzi et al., 2000) (with several typos in the corresponding expression in that paper corrected here), where $`D(z)\delta _+(0)/\delta _+(z)`$ is the linear growth factor, $`P(k)`$ is the linear power spectrum at $`z=0`$, and the factor $`(1+fx^2)^2`$, with $`f[\mathrm{\Omega }(z)]^{0.6}`$, is the correction to the cylinder length for the departure from Hubble expansion due to peculiar velocities (Kaiser, 1987).

Illustrative results are plotted for 3-$`\sigma `$ fluctuations as a function of $`\mathrm{\Delta }\theta _{\mathrm{beam}}`$, for $`z=7`$ and 8.5, in Figure On the Direct Detectability of the Cosmic Dark Ages: 21-cm Emission from Minihalos, along with the expected sensitivity limits for the planned LOFAR (300 m filled aperture) and SKA (1 km filled aperture) arrays. We plot in Figure On the Direct Detectability of the Cosmic Dark Ages: 21-cm Emission from Minihalos the predicted spectral variation of these fluctuations vs. redshift $`z`$ for illustrative beam sizes of $`\mathrm{\Delta }\theta _{\mathrm{beam}}=9^{}`$ and $`25^{}`$. These 3-$`\sigma `$ fluctuations should be observable with both LOFAR and SKA with integration times of between 100 and 1000 hours. For a $`25^{}`$ beam, for example, 3-$`\sigma `$ fluctuations can be detected for untilted $`\mathrm{\Lambda }`$CDM by both with a 100 h integration for $`z67.5`$ and a 1000 h integration for $`z11.5`$, while for a $`9^{}`$ beam, SKA can detect them after 100 h for $`z9`$ and after 1000 h for $`z13`$. Results for different values of $`z`$ and $`\mathrm{\Delta }\theta _{\mathrm{beam}}`$ are available upon request.

## Acknowledgments

We are grateful to E. Scannapieco for clarifying the effects of bias and referee P. Tozzi for his thoughtful comments. This work was supported by European Community RTN contract HPRN-CT2000-00126 RG29185 and grants NASA ATP NAG5-10825 and NAG5-10826 and Texas Advanced Research Program 3658-0624-1999.


format: text

text:
# Specular reflection of matter waves from a rough mirror

## Acknowledgments

We acknowledge the assistance of C. Aussibal and thank C. Henkel and J. Thywissen for useful discussions. This work was supported by the European Union under grants IST-1999-11055, and HPRN-CT-2000-00125, and by the DGA grant 99.34.050.


text:
# Scaling of Entanglement close to a Quantum Phase Transitions

## Abstract

In this Letter we discuss the entanglement near a quantum phase transition by analyzing the properties of the concurrence for a class of exactly solvable models in one dimension. We find that entanglement can be classified in the framework of scaling theory. Further, we reveal a profound difference between classical correlations and the non-local quantum correlation, entanglement: the correlation length diverges at the phase transition, whereas entanglement in general remains short ranged.

Classical phase transitions occur when a physical system reaches a state below a critical temperature characterized by a macroscopic order . Quantum phase transitions occur at absolute zero; they are induced by the change of an external parameter or coupling constant , and are driven by fluctuations. Examples include transitions in quantum Hall systems , localization in Si-MOSFETs (metal oxide silicon field-effect transistors; Ref. ) and the superconductor-insulator transition in two-dimensional systems . Both classical and quantum critical points are governed by a diverging correlation length, although quantum systems possess additional correlations that do not have a classical counterpart. This phenomenon, known as entanglement , is the resource that enables quantum computation and communication . The role of the entanglement at a phase transition is not captured by statistical mechanics – a complete classification of the critical many-body state requires the introduction of concepts from quantum information theory . Here we connect the theory of critical phenomena with quantum information by exploring the entangling resources of a system close to its quantum critical point. We demonstrate, for a class of one-dimensional magnetic systems, that entanglement shows scaling behaviour in the vicinity of the transition point.

There are various questions that emerge in the study of this problem. Since the ground state wave-function undergoes qualitative changes at a quantum phase transition, it is important to understand how its genuine quantum aspects evolve throughout the transition. Will entanglement between distant subsystems be extended over macroscopic regions, as correlations are? Will it carry distinct features of the transition itself and show scaling behaviour? Answering these questions is important for a deeper understanding of quantum phase transitions and also from the perspective of quantum information theory. So results that bridge these two areas of research are of great relevance. We study a set of localized spins coupled through exchange interaction and subject to an external magnetic field (we consider only spin-1/2 particles), a model central both to condensed matter and information theory and subject to intense study .

O’Connors and Wootters showed that in the Heisenberg chain the maximization of the entanglement at zero temperature is related to the energy minimization . It is known that Werner states  can be generated in a one dimensional XY model  and that temperature and magnetic field can increase the entanglement of the systems as shown for the Ising and Heisenberg models . Finally we mention the study on the role of the entanglement in the density matrix renormalization group flow and the introduction of entanglement-preserving renormalization schemes . Here we address the problem of the relation between macroscopic order, classical correlations, and quantum correlations. Therefore we analyze the entanglement near the critical point of the XY model in a transverse field. Because of the universality principle – the critical behaviour depends only on the dimension of the system and the symmetry of the order parameter – our results have much broader validity. We find that in the vicinity of a quantum phase transition the entanglement obeys scaling behaviour. On the other hand this analysis provides a clear distinction between the role of entanglement and correlations in quantum systems close to a critical point. (We have been made aware that similar work to that reported here is being performed by T. Osborne and M. Nielsen; T. Osborne, personal communication.)

The system under consideration is a spin-1/2 ferromagnetic chain with an exchange coupling J in a transverse magnetic field of strength $`h`$. The Hamiltonian is

$$H=\text{}\frac{J}{2}\underset{i=1}{\overset{N}{}}(1\text{}\gamma )\sigma _i^x\sigma _{i+1}^x+(1+\gamma )\sigma _i^y\sigma _{i+1}^y\text{}h\underset{i=1}{\overset{N}{}}\sigma _i^z$$
(1)
where $`\sigma ^a`$ are the Pauli matrices ($`a=x,y,z`$) and $`N`$ is the number of sites. We assume periodic boundary conditions. It is convenient, for later purposes, to define a dimension-less coupling constant $`\lambda =J/2h`$. For $`\gamma =1`$ Eq. (1) reduces to the Ising model whereas for $`\gamma =0`$ it is the XY model. For all the interval $`0<\gamma 1`$ the models belong to the Ising universality class and for $`N=\mathrm{}`$ they undergo a quantum phase transition at the critical value $`\lambda _c=1`$. The magnetization $`\sigma ^x`$ is different from zero for $`\lambda >1`$ and it vanishes at the transition. On the contrary the magnetization along the $`z`$-direction $`\sigma ^z`$ is different from zero for any value of $`\lambda `$. At the phase transition the correlation length $`\xi `$ diverges as $`\xi |\lambda \lambda _c|^\nu `$ with $`\nu =1`$ (Refs.  and ). We confine our interest to the entanglement between two spins, on position $`i`$ and $`j`$, in the chain. All the information needed is contained in the reduced density matrix $`\rho (i,j)`$ obtained from the ground state wave-function after all the spins except those at positions $`i`$ and $`j`$ have been traced out. The resulting $`\rho (i,j)`$ represents a mixed state of a bipartite system; a good deal of work has been devoted to quantify the entanglement in this case. We use the concurrence between sites $`i`$ and $`j`$, related to the “entanglement of formation”, and defined as

$$C(i,j)=\mathrm{max}\{0,r_1(i,j)\text{}r_2(i,j)\text{}r_3(i,j)\text{}r_4(i,j)\}$$
(2)
In Eq.(2) $`r_\alpha (i,j)`$ are the square roots of the eigenvalues of the product matrix $`R=\rho (i,j)\stackrel{~}{\rho }(i,j)`$ in descending order; the spin flipped matrix $`\stackrel{~}{\rho }`$ is defined as $`\stackrel{~}{\rho }\sigma ^y\sigma ^y\rho ^{}\sigma ^y\sigma ^y`$. In the previous definition, the eigenstates of $`\sigma ^z\{||\}`$ should be used. Translation invariance implies that $`C(i,j)=C(|ij|)`$. The concurrence will be evaluated as a function of the relative position $`|ij|`$ between the spins and the distance $`|\lambda \lambda _c|`$ from the critical point. The structure of the reduced density matrix is obtained by exploiting symmetries of the model (see caption of Fig. 1). The non-zero entries of $`\rho `$ can then be related to the various correlation functions, and the concurrence of the ground state is evaluated exactly starting from the results in Refs. , , and .

First we look at the Ising Model ($`\gamma =1`$). The first question we consider is the range of the entanglement $`\xi _E`$, that is, the maximum distance between two spins at which the concurrence is different from zero. The result is surprising: even at the critical point, where spin-spin correlations extend over a long range (the correlation length is diverging for an infinite system), the concurrence vanishes unless the two sites are at most next-nearest neighbors. The truly non-local quantum part of the two-point correlations is nonetheless very short-ranged.

In order to quantify the change of the many-body wave-function when the system crosses the critical point, we look at the derivatives of the concurrence as a function of $`\lambda `$. In this case we need to consider only the nearest neighbour and next nearest neighbour concurrence. We first discuss the behaviour of the nearest neighbor concurrence.

The results for systems of different size (including the thermodynamic limit) are presented in Fig.1. For the infinite chain $`_\lambda C(1)`$ diverges on approaching the critical value as

$$_\lambda C(1)=\frac{8}{3\pi ^2}\mathrm{ln}|\lambda \lambda _c|+const.$$
(3)
Equation (3) quantifies non-local correlations in the critical region. One aspect of this system, particularly relevant for quantum information is the study of the precursors of the critical behaviour in finite samples. This study is known as finite size scaling . In Fig.1 the derivative of $`C(1)`$ respect to $`\lambda `$ is considered for different system sizes. As expected, there is no divergence for finite $`N`$, but there are clear anomalies. The position of the minimum $`\lambda _m`$ scales as $`\lambda _m\lambda _c+N^{1.87}`$ and its value diverges logarithmically with increasing system size as

$$_\lambda C(1)|_{\lambda _m}=0.2702\mathrm{ln}N+const.$$
(4)
According to the scaling ansatz  in the case of logarithmic singularities, the ratio between the two prefactors of the logarithm in Eqs.(3) and (4) is the exponent which governs the divergence of the correlation length $`\xi `$. In this case ($`8/3\pi ^20.2702`$) it follows that $`\nu =1`$, as it is known from the solution of the Ising model . By proper scaling  and taking into account the distance of the minimum of $`C(1)`$ from the critical point, it is possible to make all the data for different $`N`$ collapse onto a single curve (Fig.2).

This figure contains the data for the lattice size ranging from $`N=41`$ up to $`N=2701`$. These results show that all the key ingredients of the finite size scaling are present in the concurrence. We note that finite size scaling scaling is fulfilled over a very broad range of values of $`N`$ which are of interest in several protocols in quantum information.

A similar analysis can be carried on for the next-nearest neighbour concurrence $`C(2)`$. Since $`_\lambda C(2)|_{\lambda _c}=0`$, the logarithmic singularity here appears first in the second derivative with respect to $`\lambda `$ (see Fig.3 legend). In the thermodynamic limit $`_\lambda ^2C(2)=0.108\mathrm{ln}|\lambda \lambda _c|+const`$. In this case also, the data collapse and finite size scaling (Fig.3) agree with the expected scaling behaviour, $`\nu =1`$. This completes the analysis of entanglement for the one-dimensional Ising model.

A cornerstone in the theory of critical phenomena is the concept of universality – that is, the critical properties depend only on the dimensionality of the system and the broken symmetry in the ordered phase. Universality in the critical properties of entanglement was verified considering the properties of the family of models defined in Eq.(1) with $`\gamma =1/2`$. The range of entanglement $`\xi _E`$ is not universal. The maximum possible distance between entangled pairs increases and tends to infinity as $`\gamma `$ tends to zero. From the asymptotic behaviour of the reduced density matrix we find that $`\xi _E`$ goes as $`\gamma ^1`$. This however has no dramatic consequences; the total concurrence $`_nC(n)`$ stored in the chain is a weakly increasing function of $`\gamma `$ (for $`0<\gamma 1`$, $`0<_nC(n)<0.2`$). More interesting is the critical behaviour of the concurrence. To be specific we consider $`C(1)`$ in the case $`\gamma =0.5`$, shown in Fig.4. As it was obtained for the Ising model, scaling is fulfilled with the critical exponent $`\nu =1`$ in agreement with the universality hypothesis.

The analysis of the ’resource’ entanglement for a condensed matter system close to a quantum critical point allows to characterize both quantitatively and qualitatively the change in the wave-function of the ground state on passing the phase transition. A notable feature which emerges is that though the entanglement itself is not an indicator of the phase transition, an intimate connection exists between entanglement, scaling and universality. In a way this analysis allows us to discern what is genuinely quantum in a zero-temperature phase transition. The results presented here might be tested by measuring different correlation functions, for example with neutron scattering, and extracting from these the entanglement properties of the ground state close to the critical point .

We finally discuss the perspectives of this work towards quantum computation. First, system sizes considered here ($`10^3`$) could be those of a realistic quantum computer. Second, the scaling behaviour found could be a powerful tool to evaluate (and hence to use) entanglement in systems having different number of qubits. In particular, close to the critical point, the entanglement depends strongly on the field – so it could be tuned, realizing an ’entanglement switch’. Last, long range correlations, typical of the critical region, might be of great importance in stabilizing the system against errors due to imperfections.

Acknowledgments. The authors would like to thank G.M. Palma, F. Plastina, and J. Siewert for helpful discussions. This work was supported by the European Community (IST-SQUIBIT), and by INFM-PRA-SSQI.


text:
# On DSS Method for a Fast Identification of the Static and Dynamic Responses of ServovalvesThis research was supported by the Russian Found for Basic research, INTAS gants (No.97-804, No.00-221) and by DIETZ automation GmbH.

## 1 Introduction

In 1950 Bill Moog invented a very reliable device to control flow in hydraulic systems by electric signals. Since that time electrohydraulic servovalves are key elements of various automatic control systems. Now it is impossible to imagine rockets, aircrafts, satellites, certain machine tools and control systems without servovalves. Such automatic systems require high performance of the servovalve. That is why research and development of such control systems as well as their maintenance require to test the dynamic and static parameters of servovalves. Unfortunately up to now, test equipment for servovalves has numerous of lacks. First of all, the test equipment needs powerful hydraulic pump stations, which require water cooling, noise insolation and frequent maintenance.111The power of pump stations reaches 100 kW and more. Moreover, this equipment requires large space and tubing for the distribution of the hydraulic fluid. The installation of such test equipment takes time and moving it is difficult. We must notice, that at the present time the complete test equipment costs more than $150,000 and a test of a servovalve takes about 15 minutes. The cost may be increased many times if to take into account that the different servovalves use different kinds of hydraulic fluids (yellow oil, red oil, kerosene, skydrol, water and so on). Notice that noise from the pump stations creates an additional problem for hydraulic test equipment. The point is that the presence of noise superimposes restrictions on the construction of buildings and the laboratories where the test stand is located. Today, we know how to make a quiet, compact, mobile and cheap hydraulic test stand.222The estimated expenses for the equipment are less $25,000.

In the present paper an approach for testing dynamic and static characteristics of servovalves and principles of construction of modern test equipment are described. The basic idea of the new approach is to divide the test process into many test sub-processes, each of which takes a very short interval of time. This idea implies the possibility to use the energy of a hydraulic accumulator to produce a short subtest and, therefore, to make obsolete powerful pump stations. Using this idea, we have developed a series of mathematical methods, wrote computer programs and made computer simulations for various servovalves. The experiments showed a possibility to identify the static characteristic of a servovalve for 1 to 2 seconds of the test time. To identify the dynamic characteristics we use a number of sub-tests each of which requires just some milliseconds of the test time. So, using a 4kW motor to fill a three-liter accumulator, the complete test of a servovalve takes less than minute. It is necessary to note that the new approach is based on the construction of a dynamic model of the servovalve which is under test and this fact turns out to be extremely useful. The knowledge of the dynamical system opens unlimited possibilities for designers and developers of control systems. They can use Mathematica, MathCad, MatLab, Simulink, LabView and other software for fast simulations on a computer as well as for complex mathematical analysis.

## 2 Statement of the problem

A servovalve is a complex dynamical system and the complete description of all its parameters, obviously, is not possible. However, most of engineers and developers are interested in some basic dynamic and static characteristics. To describe dynamics of servovalves, at the present time, people usually consider amplitude-phase frequency responses. In this case, obviously, it is assumed implicitly that dynamics of a servovalve can be described well by a linear differential equation or by a linear transfer function.333In a non-linear case the amplitude-phase response depends on amplitude and form of the input signal and, therefore, does not make big sense. In practice, to describe the dynamic parameters, engineers usually use the following system of differential equations of the second order (see for instance , , ).

$$A\frac{d^2x}{dt^2}+B\frac{dx}{dt}+Cx=u(t),\frac{dx}{dt}(0)=x(0)=0.$$
(1)
Here, $`u(t)`$ is an input electrical signal, $`x(t)`$ is the output flow of a liquid, $`A`$, $`B`$ and $`C`$ are some fixed parameters of the servovalve.

The static characteristic $`x=F(u)`$ of the servovalve is close to the linear function $`x=u/C`$, however, in practice, there is the necessity to know the static behaviour much better. In this case, by definition, to construct the static response, it is necessary to know the output flow $`x`$ at any fixed control $`u`$ . Here, of course, it is supposed, that at fixed control the servovalve is a structurally stable dynamic system with a unique stable equilibrium point. The such direct construction of the static characteristic, obviously, requires time because it is necessary to check the valve at all the points of the control range and each of the points requires a delay for the transition time. In practice, to find the statics approximately, people usually construct the Lissajous figure $`(u(t),x(t))`$, where $`u(t)`$ is a very low-frequency harmonic control signal (usually less $`0.01`$Hz). But in any case it takes time.

In the present work, to take into account the nonlinearity of the static characteristic, we shall consider the following quasilinear dynamic model of a servovalve.

$$A\frac{d^2x}{dt^2}+B\frac{dx}{dt}+f(x)=u(t),\frac{dx}{dt}(0)=x(0)=0.$$
(2)
This system differs from the system (1) by the presence of a non-linear function $`f(x)`$. In general, such a nonlinearity can lead to the appearance of a very complex chaotic behavior.444In 1918 Duffing found a chaos in a similar system with a cubic nonlinearity. See, for instance, . However, a real servovalve is a stable dynamical system and its static response is close to a linear function. Therefore we shall require that the function $`f(x)`$ is a monotone one and is close to the function $`Cx`$. Moreover we shall assume, that the roots of the characteristic equation $`A\lambda ^2+B\lambda +C=0`$ of the corresponded linear system have negative real parts. These conditions guarantee a stability of the system (2) at any constant control $`u(t)=const`$ and, therefore, the static characteristic $`x=F(u)`$ of the servovalve is the non-linear function $`x=f^1(u)`$. The experiments made by the authors of the present work have shown that the system (2) describes well both static and dynamic characteristics of servovalves. Therefore the main problem to test a servovalve is reduced to the identification of the parameters $`A`$, $`B`$ and the function $`f(x)`$ of the system (2). To solve this problem it is necessary to suggest a form of the test signal $`u(t)`$ and give an algorithm for a numerical analysis of $`x(t)`$ which identifies $`A`$, $`B`$ and $`f(x)`$. Notice that we are interested in such methods which give a good solution of the problem for a short enough time of testing. Moreover, the methods should be structurally stable to the presence of noise and measurement errors.

## 3 Amplitude-phase frequency responses of a quasilinear system

Let a solution $`x(t)`$ of the system (2) correspond to the input harmonic signal $`u(t)=\mathrm{sin}(\omega t)`$. After some transition time $`\mathrm{\Delta }t`$ the output signal $`x(t)`$ will be a periodic function with the same frequency $`\omega `$. In this case, obviously, $`x(t)`$ can be represented by the following Fourier series expansion:

$$x(t)=\underset{k=0}{\overset{\mathrm{}}{}}R_k(\omega )\mathrm{sin}(k\omega t+\phi _k(\omega )).$$

For any $`k`$, the amplitude $`R_k(\omega )`$ and the initial phase $`\phi _k(\omega )`$ of the corresponding harmonic are calculated by the formulas:

$$\begin{array}{c}R_k(\omega )=\left|K_k(i\omega )\right|,\phi _k(\omega )=\mathrm{arg}\left(K_k(i\omega )\right),\hfill \\ \\ K_k(i\omega )=\frac{\omega }{2\pi }\underset{\mathrm{\Delta }t}{\overset{\mathrm{\Delta }t+2\pi /\omega }{}}x(t)e^{ik\omega t}𝑑t.\hfill \end{array}$$
(3)
Let us note that the amplitude $`R_1(\omega )`$ and the phase lag $`\phi _1(\omega )`$ of the first harmonic are known as the amplitude-phase characteristics of the dynamic system. Namely these functions represent the most of interest at the research and the development of automatical control systems.

There are many different methods to construct the frequency responses $`R_1(\omega ),`$ $`\phi _1(\omega )`$. The most simple method is to use the formulas (3) to calculate the functions directly. In this case, to construct the frequency response, it is necessary to make measurements for a larger series of frequencies $`\omega =\omega _1,\mathrm{},\omega _N,`$ (usually $`N>20`$) and, then, to find the corresponding interpolation curves $`R_1(\omega ),`$ $`\phi _1(\omega )`$.

The other method is to identify the parameters $`A`$, $`B`$ and function $`f(x)`$ of the quasilinear system (2). In this case the test function $`u(t)`$ can be chosen from a wide class of functions (for example we can use the step signal $`u(t)=1(t)`$). In order to identify $`A`$, $`B`$ and $`f(x)`$ the following optimisation problem, for instance, can be solved:

$$\underset{A,B,f}{\mathrm{min}}\underset{k=1}{\overset{N}{}}\left(x(t_k)x_{A,B,f}(t_k)\right)^2.$$
(4)
Here $`x(t_k)`$ are the results of measurement (samples) of the output signal and $`x_{A,B,f}(t_k)`$ is a solution of the quasilinear system (2) with parameters $`A`$, $`B`$, $`f`$ for the control function $`u(t)`$ at the moments $`t_k(\tau _1,\tau _2)`$. The minimization is performed with respect to $`A`$, $`B`$ and $`f`$, where $`f`$ is a function of a finite-parametric family of functions $`f_\mu `$. So, the minimization is performed with respect to parameters $`A`$, $`B`$ and $`\mu `$. The knowledge of parameters $`A`$, $`B`$ and $`f`$, obviously, allows to construct the frequency responses. This method allows to identify the parameters of the system for an arbitrarily small time interval of measurement $`(\tau _1,\tau _2)`$. However, the method requires a high precision in measurement of the output signal $`x(t_k)`$.

In our case to find the frequency responses $`R_1(\omega ),`$ $`\phi _1(\omega )`$ of servovalves it is expedient to use a method for construction of dynamic characteristics which takes into account the proximity of the system (2) to the linear one (1). In this case, it is natural to suppose, that $`R_1(\omega ),`$ $`\phi _1(\omega )`$ are very similar to the frequency responses of the linear system (1) and our problem is reduced to search the optimum parameters $`A^{}`$, $`B^{}`$ and $`C^{}`$. Of course, $`A^{}`$, $`B^{}`$ and $`C^{}`$ can be found as a solution of (4) but we suggest an other way. It is well known, that the frequency responses of the linear system (1) can be computed by the following formulas:

$$\begin{array}{c}R(\omega )=\left|K(i\omega )\right|,\phi (\omega )=\mathrm{arg}(K(i\omega )),\hfill \\ \\ K(s)=1\mathrm{}\left(As^2+Bs+C\right).\hfill \end{array}$$
(5)
Here $`K(s)`$ is the transfer function of the system (1). In order to find parameters $`A^{}`$, $`B^{}`$ and $`C^{}`$ we suggest to solve the following optimisation problem:

$$\underset{A,B,C}{\mathrm{min}}\underset{k=1}{\overset{N}{}}\left|K(i\omega _k)K_1(i\omega _k)\right|^2,$$
(6)
where $`\omega _1,\mathrm{},\omega _N`$ is a series of test frequencies and $`K_1(i\omega _1),\mathrm{},K_1(i\omega _N)`$ are values of the transfer function which were computed by (3). In contrast to the above mentioned methods the last one does not require high precision in measurement of the output signal $`x(t_k)`$ and gives good results fast enough because the method works for small series of experiments ($`N2`$).

## 4 The DSS method

In order to identify the static response we use the DSS method. The basic idea of this method is the following. Let $`x(t)`$ be an output signal of the dynamic system (2) for a test signal $`u(t)`$. Then, obviously, the following relation is fulfilled:

$$f(x(t))=u(t)A\frac{d^2x}{dt^2}(t)B\frac{dx}{dt}(t).$$
(7)
That means that the vector-function

$$(x(t),f(x(t)))=(x(t),u(t)A\frac{d^2x}{dt^2}(t)B\frac{dx}{dt}(t))^2$$
(8)
is a parametric form of the static characteristic $`u=f(x)`$ (or $`x=f^1(u)`$). We can consider that the term

$$A\frac{d^2x}{dt^2}(t)B\frac{dx}{dt}(t)$$
(9)
in (8) is a dynamical correction of the Lissajous figure $`(x(t),u(t))`$. Namely this term allows to construct the static characteristic of the system for a short time, i.e. without a delay for the relaxation of transition processes. However, to evaluate the dynamic correction (9), it is necessary to know parameters $`A`$ and $`B`$ of the system (2). Below we describe a method of successive approximations to find $`A`$, $`B`$ and $`f(x)`$. First of all, we solve the problem (6) to find parameters $`A^{}`$, $`B^{}`$ and $`C^{}`$. The found parameters $`A^{}`$, $`B^{}`$ and the linear function $`C^{}x`$ is the first approximation $`A_0`$, $`B_0`$ and $`f_0(x)`$ of the corresponding parameters $`A`$, $`B`$ and function $`f(x)`$ of the system (2). To construct the next approximations $`A_{n+1}`$, $`B_{n+1}`$ and $`f_{n+1}(x)`$ ($`n=0,1,2,\mathrm{}`$) we use the following inductive rule. Let assume that the parameters $`A_n`$, $`B_n`$ and function $`f_n`$ are known and $`x_n(t)`$ is a solution of

$$A_n\frac{d^2x_n}{dt^2}+B_n\frac{dx_n}{dt}+f_n(x)=u(t),\frac{dx}{dt}(0)=x(0)=0.$$
(10)
Then:

1. Substituting the parameters $`A_n`$ and $`B_n`$ into the formula (7) instead of $`A`$ and $`B`$, we find the next approximation $`f_{n+1}(x)`$ of the function $`f(x)`$.
2. The next generation $`A_{n+1}`$ and $`B_{n+1}`$ of the parameters $`A`$ and $`B`$ can be found as a solution of the following optimisation:

$$\underset{A_{n=1},B_{n+1}}{\mathrm{min}}\underset{k=1}{\overset{N}{}}\left(x(t_k)x_{A_{n+1},B_{n+1},f_{n+1}}(t_k)\right)^2.$$
(11)
Here $`x(t_k)`$ are the results of measurement (samples) of output signal and $`x_{A_{n+1},B_{n+1},f_{n+1}}(t_k)`$ is a solution of the quasilinear system (2) with parameters $`A_{n+1}`$, $`B_{n+1}`$, $`f_{n+1}`$ for the control function $`u(t)`$ at the moments $`t_k(\tau _1,\tau _2)`$. Notice, the minimization is performed with respect to $`A_{n+1}`$ and $`B_{n+1}`$ only.

The described steps should be repeated until the acceptable approximation will be found. One more significance of our method is that, there are no big needs to compute the derivatives $`{\displaystyle \frac{dx}{dt}}(t_k)`$ and $`{\displaystyle \frac{d^2x}{dt^2}}(t_k)`$ numerically or to use for that special electronics. The point is that we can use the derivatives of the approximated solution $`x_n(t)x_{A_n,B_n,f_n}(t)`$ which was obtained on the $`n^{th}`$ step of the iteration process. In this case the $`n^{th}`$ approximation of the static characteristic is given by the following parametric form:

$$(x(t),f(x(t)))=(x(t),u(t)A_n\frac{d^2x_n}{dt^2}(t)B_n\frac{dx_n}{dt}(t)).$$
(12)
By (10) the formula (12) takes the following simple form:

$$(x(t),f(x(t)))=(x(t),f_n(x_n(t))).$$

So, to construct the static characteristic we do not need to compute the derivatives of the output signal $`x(t)`$. That fact is extremely important. The point is that the measured signal $`x(t)`$ contains noise and the measuring errors which lead to big problems at the differentiation.

In conclusion we must note that we did not investigate all questions about convergence of our DSS process however the computer simulations demonstrated fantastic results. Moreover, the acceptable results can be obtained even on the first step of the DSS process if to use a low frequency input signal ($`12`$ Hz).

Notice also that the DSS method can be very simply generalised onto a quasilinear systems of high order. For example, we tested our approach for the following system of differential equations of the third order:

$$A\frac{d^3x}{dt^3}+B\frac{d^2x}{dt^2}+C\frac{dx}{dt}+f(x)=u(t),\frac{d^2x}{dt^2}(0)=\frac{dx}{dt}(0)=x(0)=0.$$