text:
# Mixing and relaxation dynamics of the Henon map at the edge of chaos

## 1 Introduction

Nonlinear one-dimensional dissipative systems at their critical points (like chaos threshold, bifurcation points, etc) have become one of the vivid areas of research since 1997, after pioneering work of Tsallis et al. . One of the relevant reason for this, in our opinion, was that they were able to constitute one of the first methods of inferring the entropic index $`q`$ of nonextensive entropy

$$S_q=\frac{1_{i=1}^Wp_i^q}{q1},$$
(1)
and consequently determining the proper thermostatistics which the system under consideration obeys, from the microscopic dynamics of the system. This method depends on the measure of the divergence of initially nearby orbits. More precisely, the effect of any uncertainty on initial conditions exhibits, for periodic and chaotic orbits, an exponential temporal evolution as $`\xi (t)=lim_{\mathrm{\Delta }x(0)0}\frac{\mathrm{\Delta }x(t)}{\mathrm{\Delta }x(0)}\mathrm{exp}(\lambda t)`$, where $`\lambda `$ is the standard Lyapunov exponent, $`\mathrm{\Delta }x(0)`$ and $`\mathrm{\Delta }x(t)`$ are the uncertainties at times 0 and $`t`$ ($`\lambda <0`$ for periodic cases, and $`\lambda >0`$ for chaotic ones); whereas at critical points (e.g., chaos threshold) this feature is related to a power-law type such as

$$\xi (t)=\left[1+(1q_{mix})\lambda _{q_{mix}}t\right]^{1/(1q_{mix})}$$
(2)
where $`\lambda _{q_{mix}}`$ is a generalized version of the Lyapunov exponent (here, $`mix`$ stands for mixing and the aim of this notation will be transparent soon). This equation corresponds to power-law growth of the upper bounds of a complex time dependence of $`\xi (t)`$ and these upper bounds allow us to infer the proper value of $`q_{mix}`$ for the system under consideration. Let us call this Method I of finding the proper $`q_{mix}`$ of a given dynamical system.

The second method (Method II) is based on the geometrical aspects of the critical attractor of the dynamical system at the chaos threshold. Scaling arguments have shown that the appropriate $`q_{mix}`$ is related to the multifractal structure of the critical dynamical attractor by

$$\frac{1}{1q_{mix}}=\frac{1}{\alpha _{min}}\frac{1}{\alpha _{max}}(q_{mix}<1)$$
(3)
where $`\alpha _{min}`$ and $`\alpha _{max}`$ are the end points of the multifractal singularity spectrum of the critical attractor and characterize the most concentrated and most rarefied regions respectively on the attractor. This fascinating relation relates power-law sensitivity or mixing (left-hand side) with purely geometrical quantities (right-hand side) and constitutes another independent method for inferring the proper $`q_{mix}`$ of the dynamical system under study.

Finally, a third method (Method III) has been proposed for obtaining $`q_{mix}`$, dealing with the entropy increase rates . This method is the one where nonextensive entropy $`S_q`$ enters directly and consequently it provides also a justification for the use of notation $`q`$ within the two previous methods. The procedure for Method III is the following: We first partition the phase space into $`W`$ cells, then we choose one of these cells (randomly or not) and locate $`N`$ initial conditions all inside this cell. At any time $`t`$, we have $`N_i(t)`$ points in each cell (naturally $`_{i=1}^WN_i(t)=N`$) and thus we can define a set of probabilities through $`p_i(t)=N_i(t)/N`$ ($`i`$), which enable us to calculate the entropy $`S_q(t)`$. Using this procedure, one can define a generalized version of the Kolmogorov-Sinai entropy as $`K_qlim_t\mathrm{}lim_W\mathrm{}lim_N\mathrm{}S_q(t)/t`$, and the proper $`q_{mix}`$ value of a given system is the one for which the entropy production $`K_q`$ is finite. If $`q<q_{mix}`$ ($`q>q_{mix}`$), $`K_q`$ diverges (vanishes).

These three independent methods (from now on referred to as mixing-based approach) have already been tested and verified with numerical calculations for a number of one-dimensional dissipative systems such as the logistic map , $`z`$-logistic map family , the circle map , $`z`$-circular map family , asymmetric logistic map family  and the single site map . Quite remarkably, it is found that all three methods give one and the same value for $`q_{mix}`$ of a given system. It is worth mentioning here that all these numerical findings have been supported very recently by an analytical work .

We can now turn our attention to the second approach, which will be referred to as relaxation-based approach. This approach has been developed by one of us and his collaborators  and it is based on investigating numerically the critical temporal evolution of the volume of the phase space occupied by an ensemble of initial conditions spread over the entire phase space. In that sense, this approach differs from the mixing-based approach, where a set of initial conditions spread in the vicinity of the inflection point is used. In the beginning of time, all cells are occupied and $`S_q`$ is bounded by the equiprobability as $`S_q=(W^{1q}1)/(1q)`$. As time evolves, the number of occupied cells (say $`W_{occ}`$) starts contracting as $`\left[1+(q_{rel}1)t/\tau _q\right]^{1/(1q_{rel})}`$, here $`q_{rel}>1`$, $`rel`$ stands for relaxation and $`\tau _q>0`$ is a characteristic relaxation time. In , $`z`$-logistic map family has been analyzed and $`q_{rel}(t\mathrm{})`$ values for $`z`$ parameter has been obtained.

The existence of two different classes of $`q`$ values coming from these two different approaches (namely, mixing-based and relaxation-based approaches) is also present in the results recently obtained for fully developed turbulence . Therefore, all these findings naturally force us to ask an intriguing question: are these classes connected to each other or not ? The affirmative answer to this question came very recently by Borges at al. . Their procedure goes as follows: Beginning with $`N`$ points located inside a single cell (similar to Method III explained before), the dynamics of the system (for $`t\mathrm{}`$) leads it to a state with a stationary value of $`S_{q_{mix}}`$ — this is the attractor of the system for the particular conditions of $`W`$, $`N`$ and the initial cell chosen. The path that the system evolves towards its attractor changes as the initial cell is changed. Some of the initial cells spread the points very fast, increasing the entropy to values much higher than that of the attractor (there is an overshooting in the time evolution of $`S_{q_{mix}}`$). Consequently, for these specific cases, the attractor is achieved from above. The procedure consists in observing the rate in which the system goes to its attractor, for those initial cells which yield large overshoots. Specifically, the time evolution of $`\mathrm{\Delta }S_{q_{mix}}S_{q_{mix}}(t)S_{q_{mix}}(\mathrm{})`$ is followed. The procedure is repeated for increasing values of $`W`$ (and correspondingly increasing values of $`N`$). A power-law regime appears and becomes more pronounced for increasing $`W`$. The slope (in a log-log plot) of the range of the $`\mathrm{\Delta }S_{q_{mix}}`$ curve that follows a power-law is identified with $`1/(q_{rel}(W)1)`$ (with $`q_{rel}>1`$). It is worth to stress that $`q_{rel}=q_{rel}(W)`$, i.e., the rate in which the system reaches its final state depends on the coarse graining ($`W`$) adopted. The connection between the mixing-based and the relaxation-based approaches is established by the scaling law

$$q_{rel}(\mathrm{})q_{rel}(W)W^{|q_{mix}|}.$$
(4)
In this work, our aim is to provide further support to this connection by analyzing, for the first time, a two-dimensional dissipative system within this context. To accomplish this task, we will focus on the Henon map

$`\begin{array}{cc}x_{t+1}=& \\ y_{t+1}=& \end{array}\begin{array}{cc}1ax_t^2+y_t& \\ bx_t& \end{array}`$ (9)
where $`a`$ and $`b`$ are map parameters. This map reduces to the standard logistic map when $`b=0`$, whereas it becomes conservative when $`b=1`$; in between these two cases, it is a two-dimensional dissipative map. Specifically, we focus on small values of the $`b`$ parameter such as $`b=0.001`$, $`0.01`$, $`0.1`$, etc. This map has been studied within the mixing-based approach using the three different methods and it is found that the $`q_{mix}`$ value for all values of $`b`$ parameter coincides with that of the logistic case, as expected . Here, we first investigate the Henon map within the relaxation-based approach and then establish the connection between these two different approaches following the same lines of .

## 2 Numerical Results

As introduced in the previous Section, for the relaxation approach, we follow the dynamical evolution of an ensemble of initial conditions uniformly distributed over the phase space. In practice, a partition of the phase space on $`W`$ cells of equal size is performed and a set of $`N`$ initial copies of system is followed whose initial conditions are uniformly distributed over the entire phase space. Our results for the temporal evolution of the volume occupied by the ensemble are given in Fig. 1 for various $`b`$ values. As it is evident, after the power-law contraction of the volume occupied by the ensemble sets up, the slope \[$`\mu =1/(1q_{rel})`$\] is independent from the $`b`$ values and coincides with that of the logistic map (namely, $`\mu =0.71`$ hence $`q_{rel}=2.41`$).

Now we are at the position to establish the connection between this result and the one obtained in . Fig. 2 shows the temporal evolution of $`S_{q_{mix}}`$ ($`q_{mix}=0.2445`$) for different values of $`W`$, at the edge of chaos for $`b=0.1`$. (It would be recommended that $`NW`$, but we adopted here $`N=W`$, due to restricted computer capabilities. Of course, we compared this case with $`N=10W`$, and we found that the differences are negligible.) It becomes evident that the limits $`lim_W\mathrm{}`$ and $`lim_t\mathrm{}`$ are non-commutable. If the former is taken first, the system will never relax to a stationary state (equilibrium is but a particular case), which leads to the conclusion that reaching equilibrium depends on the coarse graining adopted, and ultimately, depends on the observer! Equilibrium is only possible with incomplete knowledge ($`W<\mathrm{}`$). To arrive to the scaling law stated by Eq. (4), we must take first the temporal limit, and then $`W\mathrm{}`$. When we consider the downhill side of these curves and take the slopes of $`\mathrm{\Delta }S_q`$ versus time (as explained earlier), then we end up with Fig. 3, where the verification of the scaling (4) is evident. On the other hand, the extrapolation value of $`q_{rel}(W\mathrm{})`$ is computed as $`2`$, which is different from the logistic case ($`2.4`$). This discrepancy may be originated, most probably, from the fact that $`b=0.1`$ value is too large for analyzing, as it is also observed from the deterioration seen in Fig. 1 of  (Method I) for increasing values of the $`b`$ parameter. A careful investigation of the effect of increasing $`b`$ values to $`q_{rel}(W\mathrm{})`$ would be a subject of an extensive computational effort.

## 3 Conclusions

After the analysis of two-dimensional dissipative Henon map at the chaos threshold using three distinct methods (Methods I, II and III) , it became evident that the scenario, which we called here mixing-based approach, is valid not only for one-dimensional dissipative map families, but also for two-dimensional cases. At this point, what missing in the literature were (i) to study the relaxation-based approach for a two-dimensional dissipative system making use of the same procedure applied so far to one-dimensional cases , and (ii) to verify the connection (given by the scaling relation (4)) between the mixing-based and relaxation-based approaches for a two-dimensional dissipative system. This work constitutes the first step along this line and provides some further elements that reinforce the idea of weak mixing.

## Acknowledgments

This work has been supported by the Turkish Academy of Sciences, in the framework of the Young Scientist Award Program (UT/TUBA-GEBIP/2001-2-20).

Figure Captions

Figure 1 \- Time evolution of the volume occupied by the ensemble for various values of $`b`$ parameter. Thick line corresponds to a straight line with $`\mu =0.71`$.

Figure 2 \- Time evolution of $`S_{q_{mix}}`$, with $`q_{mix}=0.2445`$ (log-log scale). The Henon map is at its chaos threshold for $`b=0.1`$, that is, $`a_c=1.26359565\mathrm{}`$. Points correspond to $`W=700\times 700`$, $`1000\times 1000`$, $`2000\times 2000`$, $`4000\times 4000`$ and $`5000\times 5000`$.

Figure 3 \- Scaling law for the Henon map (Eq. (4)), at the chaos threshold with $`b=0.1`$. Correlation coefficient is $`0.99`$. To improve the correlation, it becomes necessary to use greater values of $`W`$ and $`N`$, that means better computer machinery. Dashed line indicates extrapolated value of $`q_{rel}(W\mathrm{})1.987`$ (that is different from the value of the logistic ($`b=0`$) map: $`q_{rel}(W\mathrm{})2.4`$).


format: text

text:
# Localization of shocks in driven diffusive systems without particle number conservation

## I Introduction

One-dimensional driven diffusive systems proved to be a rewarding research topic in the past years Habil . They were shown to exhibit boundary induced phase transitions Krug , spontaneous symmetry breaking Evans ; Mukamel95  and phase separation Mukamel1 ; Mukamel2 . Recently, attention turned to the case of systems without particle conservation in the bulk. In Ref. Kolomeisky2,  the effect of a single detachment site in the bulk of an asymmetric simple exclusion process (ASEP) was studied. In Refs.Challet ; Parmeggiani  the interplay of the simplest one-dimensional driven model, the totally asymmetric exclusion process (TASEP) with local absorption/desorption kinetics of single particles acting at all sites, termed ’Langmuir kinetics’ (LK) was considered. These models were inspired by the dynamics of motor proteinsmotorprot4  which move along cytoskeletal filaments in a certain prefered direction while detachment and attachment can also occur between the cytoplasm and the filament, and, in a very different setting, by dynamics of limit orders in a stock exchange market. Being an equilibrium process, LK is well understood, while the combined process of TASEP and LK showed the new feature of a localized shock in the density profile of the stationary state Parmeggiani .

The TASEP is defined on a one-dimensional lattice of size $`L`$. Each site can either be empty or occupied by one particle. In the bulk particles can hop from site $`i`$ to site $`i+1`$ with unit rate, provided the target site is empty. At site 1 particles can enter the lattice from a reservoir with density $`\rho _{}`$ provided the site is empty. They can leave the system from site $`L`$ into a reservoir of density $`\rho _+`$ with rate $`1\rho _+`$. Thus in the interior of the lattice the particle number is a conserved quantity. The phase diagram and steady states of the TASEP as a function of the boundary rates are known exactly Ligg75 ; Domany ; ASEP . Furthermore a theory of boundary induced phase transitions exists, which explains the phase diagram quantitatively in terms of the dynamics of shocks Kolomeisky . In the stationary state these shocks exist as an upward density shock along the coexistence line between the high and low density phases, i.e., they connect a region with low density to the left of the shock position with a high density region to its right. The shock performs symmetric random walk between the boundaries of the system.

One may equip the system with the additional feature of local particle creation at empty sites with rate $`\omega _a`$ and annihilation with rate $`\omega _d`$ (see Fig. 1)Parmeggiani ; Challet . In the thermodynamic limit $`L\mathrm{}`$ there are three regimes to be distinguished: If $`\omega _a`$ and $`\omega _d`$ are of an order larger than $`1/L`$ the steady state of the system will be that of Langmuir kinetics, i.e., there will be a uniform density of $`K=\omega _a/(\omega _a+\omega _d)`$ in the system. In case of $`\omega _a`$ and $`\omega _d`$ being of smaller order than $`1/L`$, the local kinetics is negligible and the system will behave as the TASEP. The case of the local rates being of order $`1/L`$ is the most interesting one and will be investigated further on. Writing

$$\omega _a=\mathrm{\Omega }_a/L,\omega _d=\mathrm{\Omega }_d/L$$
(1)
the phase diagram can be formulated in terms of $`\mathrm{\Omega }_a`$, $`\mathrm{\Omega }_d`$, $`\rho _{}`$ and $`\rho _+`$. In Ref. Parmeggiani,  it was shown that for $`\mathrm{\Omega }_a`$ and $`\mathrm{\Omega }_d`$ fixed, the phase diagram as a function of $`\rho _{}`$ and $`\rho _+`$ does not only exhibit the low-density and high-density phases known from the TASEP, but also a high-low coexistence phase. In this phase the shock does not move in the system but its position is a function of the rates $`\rho _{}`$ and $`\rho _+`$ (see Fig. 2).

Parmeggiani et.al. presented not only Monte Carlo simulations but derived also a mean field equation for the density profile which was shown to coincide with the simulation profiles. We argue here that the mean field approximation can not be used in general. The coincidence with the Monte Carlo in Ref. Parmeggiani,  is due to lack of correlations in true steady state of the TASEP. We claim that the stationary density profile can be derived in general using a hydrodynamic equation taking correlations into account (in case of the TASEP this equation is equal to that obtained with a mean field approach). For the Katz-Lebowitz-Spohn (KLS) model, which is a generic model of interacting driven diffusive systems KLS ; Schmittmann  we show that this hydrodynamic equation correctly describes the density profiles on a quantitative level, while a mean field approach would fail to reproduce even basic qualitative features of the system, e.g., phase separation into three distinct density regimes.

## II Hydrodynamic equation

In the following we are interested in the $`L\mathrm{}`$ limit which we perform by tuning the lattice spacing $`a=1/L0`$ and rescaling of time $`t=t_{\text{lattice}}/L`$ (Eulerian scaling) to get the continuous (hydrodynamic) limit of the model. In this framework $`\mathrm{\Omega }_{a,d}`$ are the attachment/detachment rates per unit length. We claim that the hydrodynamic equation describing the time dependence of the local density $`\rho (x)`$ for a general driven diffusive system with Langmuir kinetics takes the form

$$\frac{}{t}\rho +\frac{}{x}j(\rho )=(\rho ),$$
(2)
where $`j(\rho )`$ is the exact current in a driven diffusive system with homogeneous density $`\rho `$ without LK and $`(\rho )`$ the source term describing the Langmuir kinetics. Here, we consider only that choice of $`(\rho )`$ which corresponds to the process depicted on Fig. 1:

$$(\rho )=\mathrm{\Omega }_a(1\rho (x,t))\mathrm{\Omega }_d\rho (x,t)$$
(3)
Other choices of $`(\rho )`$, which might e.g. describe the local annihilation of particle pairs are to be discussed in a forthcoming publication coagulation .

As is usually done in the rigorous derivation of the hydrodynamic limit of conservative systems Rezakhanlou91 , our nonconservative Eq.(2) implicitly assumes that the system is locally stationary because the exact form of the stationary flux is used. We argue that this assumption is justified since the nonconservative part of the dynamics of the system at macroscopic scale is so slow that locally the system reaches stationarity with respect to the conservative part of the dynamics. Any finite perturbation caused by the nonconservative dynamics would travel a macroscopic distance and hence dissipate before interacting with another perturbation. Hence the hydrodynamic description (after time rescaling $`tϵt`$) is adequate for describing the full dynamics. For physical insight in the formation of shocks one needs other tools which are discussed below.

Rewriting equation (2) by using that $`_t\rho (x,t)=0`$ in the stationary state and $`_xj=j/\rho \rho /x`$ yields for the stationary density profile $`\rho (x)`$:

$$v_c(\rho )\frac{\rho (x)}{x}=(\rho ).$$
(4)
Here, $`v_c=j/\rho `$ is the collective velocity, i.e., the drift velocity of a center of mass of a local density perturbation on a homogeneous stationary background with the density $`\rho `$ (for system with the Langmuir kinetics switched off) Habil ; Kolomeisky . The stationary density profile has to satisfy (4) as well as the boundary conditions $`\rho (0)=\rho _{}`$ and $`\rho (1)=\rho _+`$. As equation (4) is of the first order there will be in general no smooth solution fitting both boundary conditions. In the original lattice model this discrepancy is resolved by appearance of shocks and/or boundary layers. To regularize the problem, one can add to (2) and correspondingly to (4) a vanishing viscosity term

$$v_c(\rho )\frac{\rho (x)}{x}=(\rho )+\nu \frac{^2\rho (x)}{x^2},$$
(5)
where $`\nu >0`$ is of order of $`1/L`$. This term makes the hydrodynamic equation second order and ensures a smooth solution fitting both boundary conditions. The shock has then a width of order $`1/L`$ (see Ref. Parmeggiani, ), i.e., in the thermodynamic limit the rescaled solution becomes discontinuous. We claim that equation (5) gives the same result in the $`L\mathrm{}`$ limit as the Monte Carlo, therefore it can be used as a tool to compute the stationary density profile. The main difference between (5) and the MC is that the former does not take fluctuations into account which leads to a shock width of order $`1/L`$ while in a MC after averaging it is of the order of $`1/\sqrt{L}`$ due to the fluctuation of the shock position.

The stationary density profile for a given $`j(\rho )`$ and parameters $`\mathrm{\Omega }_a`$, $`\mathrm{\Omega }_d`$, $`\rho _{}`$ and $`\rho _+`$ can be derived from the flow-field of the differential equation (4) by using the rules, formulated and explained below:

1. In the interior of the lattice the stationary density profile either follows a line of the flow field of the differential equation (4) or makes a jump. Jumps can only occur between densities yielding the same current, i.e., the current is continuous in the interior of the lattice.
2. Let $`\rho _\pm ^{}`$ be defined as limiting left and right densities with the boundary layers cut away:

$$\rho _{}^{}=\underset{x+0}{lim}\rho (x),\rho _+^{}=\underset{x10}{lim}\rho (x),$$

where $`\rho (x)`$ is the stationary profile in the hydrodynamic limit. The boundary layer at $`x=0`$ ( i.e., if $`\rho _{}\rho _{}^{}`$) has to satisfy the following condition:

$`\text{if }\rho _{}<\rho _{}^{}\text{ then }j(\rho )>j(\rho _{}^{})\text{ for any }\rho (\rho _{},\rho _{}^{})`$ (6)
$`\text{if }\rho _{}>\rho _{}^{}\text{ then }j(\rho )<j(\rho _{}^{})\text{ for any }\rho (\rho _{}^{},\rho _{})`$ (7)
The condition for the stability of the boundary layer at $`x=1`$ (if there is) is similar:

$`\text{if }\rho _+^{}<\rho _+\text{ then }j(\rho _+^{})<j(\rho )\text{ for any }\rho (\rho _+^{},\rho _+)`$ (8)
$`\text{if }\rho _+^{}>\rho _+\text{ then }j(\rho _+^{})>j(\rho )\text{ for any }\rho (\rho _+,\rho _+^{})`$ (9)
3. Shocks between a density $`\rho _l`$ to the left of the shock and $`\rho _r`$ to the right of the shock are stable only if they are stable in the absence of Langmuir kineticsHabil ; Hager .

Remarks:

* Although LK does not conserve locally the number of particles, Eq. (2) with the vanishing viscosity added (5) can be rewritten formally in the form

$$\frac{\rho (x,t)}{t}+\frac{}{x}\stackrel{~}{j}(x,t)=0,\stackrel{~}{j}(x,t)=j(\rho )_A^x(\rho )𝑑x\nu \frac{\rho }{x}(t)$$
(10)
where $`(t)`$ is some time-dependent function. Suppose that there is a shock at the position $`X_0`$ connecting the densities $`\rho _l`$ and $`\rho _r`$. The mass transfer across the shock is

$$\frac{}{t}_{X_00}^{X_0+0}\rho (x,t)𝑑x=\stackrel{~}{j}(X_0+0,t)\stackrel{~}{j}(X_00,t)=j(\rho _r)j(\rho _l),$$
(11)
since the Langmuir term and the viscosity term change only infinitesimally across the shock. In the stationary state, the RHS of (11) vanishes which explains the rule (A).
* The rule (B) is due to the fact that in the boundary layer of the vanishing length $`\delta l0`$, the LK term in (10) can be neglected. Consequently, for the stationary current at the boundaries we have $`\stackrel{~}{j}(x)=j(\rho (x))\nu \frac{\rho }{x}=J`$, which yields the known maximization/minimization principleHabil ; Popkov , and is equivalent to rule (B). Indeed at the left boundary $`J=j(\rho _{}^{})`$ (see (6) for notations), and if, e.g., $`\rho _{}<\rho _{}^{}`$, then $`\frac{\rho }{x}>0`$. Consequently, we obtain $`j(\rho _{})=J+\nu \frac{\rho }{x}>J`$, which is exactly (6). Analogously one obtains (7)-(9).
* The rule (C) is explained by the marginal role the Langmuir kinetics plays locally in space and in time. The first, LK is very slow locally for large $`L`$ (see (1)), and the second, it acts “orthogonally” on the particle distribution, not affecting directly the particle motion. Hence, the local perturbations will still spread with the velocity corresponding to the local density level $`\rho `$, thus rendering the same stability conditions for a shock as for the diffusive system without LK.

Condition (C) is easy to check geometrically through the current-density relation: an upward (downward) shock is stable if the straight line connecting the points $`(\rho _l,j(\rho _l))`$ and $`(\rho _r,j(\rho _r))`$ stays below (above) the $`j(\rho )`$ curveHager ; Popkov . Because of criterion (A) these lines are always horizontal in this case which gives zero mean velocity (but not localization) for the shock in absence of Langmuir kinetics.
* In the cases we have considered (ASEP, KLS model), the rules (A)-(C) define an unique stable solution (see an Appendix) and we believe that this is true also in general case, i.e., for arbitrary $`j(\rho )`$ dependence and given choice (3) of Langmuir kinetics.

In the following we apply the general theory to specific models.

## III Revisiting the ASEP with Langmuir kinetics

Using the differential equation (4) and the rules given above we reconsider the ASEP with Langmuir kinetics Challet ; Parmeggiani . Here, the current-density relation is given by $`j(\rho )=\rho (1\rho )`$, which yields $`v_c(\rho )=12\rho `$. Thus equation (4) becomes

$$\left(12\rho (x)\right)_x\rho (x)=\mathrm{\Omega }_a(\mathrm{\Omega }_a+\mathrm{\Omega }_d)\rho (x),$$
(12)
which is identical with the mean field equation in Ref. Parmeggiani,  in the thermodynamic limit. We would like to stress that this coincidence is caused by the fact that the mean field current-density relation for the TASEP is exact. As is demonstrated below, equation (4) also holds when this is not the case, as e.g. for the one-dimensional KLS model.

Due to rule (A) as stated above (continuity of the current in the interior of the lattice) shocks in the interior can only occur in the case that $`\rho _l=1\rho _r`$, as $`j(\rho )`$ is symmetric to $`\rho =1/2`$. Rule (C) (stability of the shock) furthermore requires that $`\rho _r>\rho _l`$. These observations coincide with the findings of Parmeggiani .

We also applied our rules to $`k`$-hop exclusion models Binder  ( with LK added), which are a generalization of the TASEP with stationary product measures and asymmetric current-density relations. Due to this fact shocks appear, which are non-symmetric with respect to $`\rho =1/2`$. MC simulations are in full accord with our predictions Willmann\_unpublished .

## IV KLS model with Langmuir kinetics

A much studied one-dimensional driven diffusive system with interactions between the particles is the following variant of the KLS model Mukamel2 ; Hager ; Popkov :

In the interior, particles at site $`i`$ move to site $`i+1`$, provided it is empty, with a rate that depends on the state of sites $`i1`$ and $`i+2`$.

$`01000010`$ with rate $`1+\delta `$
$`11001010`$ with rate $`1+ϵ`$
$`01010011`$ with rate $`1ϵ`$
$`11011011`$ with rate $`1\delta `$
At site 1 particles can enter the lattice provided the target site is empty. The rate depends on the state of site 2. Similarly, particles can leave the system at site $`L`$ with a rate depending on the state of site $`L1`$. The boundaries mimic the action of reservoirs with densities $`\rho _{}`$ and $`\rho _+`$. For $`\rho _{}=\rho _+`$ the stationary state is that of an one-dimensional Ising model with boundary fields. The current-density relation can be calculated exactly using transfer matrix techniques Hager . It turns out that for strong enough repulsion between the particles ($`ϵ0.9`$) a current-density relation with two maxima arises (see Fig. 3). The parameter $`\delta `$ determines the skewness of $`j(\rho )`$ with respect to the vertical line $`\rho =1/2`$. For $`\delta =0`$, the system has particle-hole symmetry resulting in $`j(\rho )`$ being symmetric with respect to 1/2. For simplicity we consider this case in the rest of the paper.

The phase diagram of this family of models with strong particle repulsion is known to exhibit 7 different phases, among them two maximal-current and one minimal-current phase. The phase diagram is determined by the interplay of diffusion, branching and coalescence of shocks Popkov .

When equipping these models with Langmuir kinetics one expects that a very rich phase diagram with many more than the original 7 phases will appear. We will not attempt to give this full phase diagram here, but instead present two new features, which cannot be observed in systems without a concave region in the current-density relation: localized downward shocks and double shocks.

### IV.1 Localized downward shocks

In the regime where the current-density relation of the KLS model exhibits two maxima at densities $`\rho _1^{}`$ and $`\rho _2^{}`$, where $`\rho _1^{}<\rho _2^{}`$ and a minimum at $`\rho =1/2`$ (at $`\delta =0`$) there is a region where downward shocks are stable according to Ref. Hager, ; Popkov,  (and rule (C)). These are characterized by $`\rho _l(0.5,\rho _2^{})`$ and $`\rho _r(\rho _1^{},0.5)`$. This suggests that localized downwards shocks may appear when introducing the kinetic rates. In deed, in the KLS model with Langmuir kinetics for certain values of the boundary densities $`\rho _{}`$ and $`\rho _+`$ ,which strongly depend on the kinetic rates $`\mathrm{\Omega }_a`$ and $`\mathrm{\Omega }_b`$, one gets a stable downward shock according to rules (A,B,C). We give an example for this case on Fig. 4.

One can see that employing the general theory described above yields a stationary profile with a localized downward shock, which coincides with the MC results up to finite size effects, while a simple mean field approach would fail as it would not be able to capture the difference between the KLS model with $`ϵ>0`$ and the TASEP (KLS with $`ϵ=0`$).

### IV.2 Localized double shocks

Let $`\rho _{1,2}^{}`$ be defined as the inflection points of the current-density relation ($`\rho _1^{}<\rho _2^{}`$). As is known from the studies of the KLS model Popkov ; Hager , if we start an infinite system from a step-like initial density profile with $`\rho _{}(\stackrel{~}{\rho }_1,\rho _1^{})`$ on the left and $`\rho _+(\rho _2^{},\stackrel{~}{\rho }_2)`$ on the right, we get a time-dependent solution having two shocks: One of these has negative mean velocity, while the other has positive and in the middle there is an expanding region with $`\rho =1/2`$ (for $`\delta =0`$) which corresponds to the minimal current phase in a system with open boundariesHager ; Popkov .

This leads us to the conjecture that introducing the kinetic rates for certain values of $`\rho _{},\rho _+,\mathrm{\Omega }_a,\mathrm{\Omega }_d`$ one may achieve a stable double shock structure. In Fig. 6 we present an example for such a case. Application of rules (A,B,C), which is presented in detail in Appendix A, yields the same double shock structure as the MC up to finite size effects. Note, that a simple mean field approach could not predict a double shock.

## V Conclusions

In this work we present a hydrodynamic equation which, together with some rules treating the discontinuities, correctly describes the stationary states of one-dimensional driven diffusive systems with Langmuir kinetics and open boundaries. It captures both systems without correlations in a steady state ( as e.g. the TASEP and the $`k`$-hop exclusion models) and systems with correlations as the KLS model. For the latter the two new phenomena of a stationary localized downward shock and a localized double shock (corresponding to phase separation to three distinct regions) were presented which a mean field approach would not reproduce. The exact currect of driven diffusive systems without LK enters the hydrodynamic description since the bulk has sufficient time to relax between subsequent annihilation/creation events. An interesting, paradoxical feature of these phenomena is that fluctuating shocks get localized due to extra noise (LK), which is highly unexpected.

###### Acknowledgements.

We wish to thank the authors of Ref. Parmeggiani,  for communicating their results prior to publication and M. Salerno for useful discussions. AR acknowledges financial support by Deutsche Forschungsgemeinschaft. ABK acknowledges the support of the Camille and Henry Dreyfus New Faculty Awards Program (under Grant NF-00-056) and the hospitality of Forschungszentrum Julich.

*

## Appendix A Double shock density profile from the rules (A)- (C)

Here we demonstrate how one determines the stationary density profile using the rules (A), (B) and (C) from the section II. As an example we take the parameters which yield a double (localized) shock structure in the KLS model ($`\rho _{}=0.23`$, $`\rho _+=0.745`$, $`\mathrm{\Omega }_a=0.03`$ and $`\mathrm{\Omega }_d=0.01`$). The KLS-model parameters are: $`\delta =0,ϵ=0.9`$ (see section IV).

First suppose that there is a boundary layer at $`x=0`$. According to rule (B) it is stable only if $`\rho _{}^{}>1\rho _{}=0.77`$. If this is the case then in the bulk there is no allowed jump since these trajectories of the flow-field (see Fig. 8) stay always above $`\rho =0.75`$ (rules (A) and (C)) which yields $`\rho _+^{}>0.75`$. But then the boundary layer at $`x=1`$ does not satisfy rule (B). This contradiction shows that there is no boundary layer at $`x=0`$. One can use the same argument to show that there is no boundary layer at $`x=1`$ either.

Now one can see that the stationary density profile close to the left boundary follows the line of the flow-field for which $`\rho (x=0)=\rho _{}=0.23`$. Since there is no boundary layer at the right end it is clear that somewhere in the bulk it has to make a jump.

Note that this trajectory crosses the line $`\rho =\stackrel{~}{\rho }_1`$ at $`x=x_1`$. Suppose that the jump takes place before at $`x<x_1`$. In this case, according to rule (A) it would jump over $`\stackrel{~}{\rho }_2=1\stackrel{~}{\rho }_1`$ which would result in a boundary layer at $`x=1`$ which is not allowed. If the jump takes place at $`x>x_1`$ then $`\rho _1^{}<\rho _r<0.5`$ and since from this region there is no allowed jump it would end up at $`\rho _1^{}<\rho _+^{}<0.5`$ resulting again in an unstable boundary layer on the right side. This shows that the jump is located at $`x=x_1`$ and from here the density profile follows the trajectory which starts at $`x=x_1`$ with the value $`\rho =0.5+0`$.

One can easily see that we need another jump to connect this trajectory with the one which ends at $`x=1`$ with $`\rho =\rho _+`$. Applying rule (A) (continuity of the current) we can get the point $`x_2`$ where the second jump is located.


text:
# Fluctuational transitions through a fractal basin boundary

## Abstract

Fluctuational transitions between two co-existing chaotic attractors, separated by a fractal basin boundary, are studied in a discrete dynamical system. It is shown that the mechanism for such transitions is determined by a hierarchy of homoclinic points. The most probable escape path from the chaotic attractor to the fractal boundary is found using both statistical analyses of fluctuational trajectories and the Hamiltonian theory of fluctuations.

The mechanism of fluctuational escape from a chaotic attractor (CA) through a fractal basin boundary (FBB) represents one of the most challenging unsolved problems in fluctuation theory Kautz:87 ; Beale:89 ; Grassberger:89 ; Hamm:91 . The unpredictable and highly complex stochastic behavior of such systems arises in part from the presence of limit sets of complex geometrical structure, and in part from the fractality of the basin boundary Guken ; Ott . For this reason, the central question – whether or not there exists a generic mechanism for fluctuational transitions through the FBB – has remained unanswered. More specifically, it has remained unclear: (i) if boundary conditions can be found both on the CA and on the FBB; (ii) if there exits a unique escape path from the CA to the FBB; (iii) whether this path can be determined using the Hamiltonian theory of fluctuations; (iv) if there is any deterministic structure involved in the transition through the FBB itself; and (v) what effect is exerted by the noise intensity. If transitions across FBBs are characterised by general features, a knowledge of them could considerably simplify analyses of both stability and control for chaotic dynamical systems, which are problems of broad interdisciplinary interest Fradkov:98 ; Boccaletti:00 .

A promising approach to the solution of this problem is based on the analysis of fluctuations in the limit of very small noise intensity. In this limit, a stochastic dynamical system fluctuates to remote states along certain most probable deterministic paths Onsager:53 ; Dykman:92a ; Luchinsky:97 , corresponding to rays in the WKB-like asymptotic solution of the Fokker-Planck equation Freidlin:84 . The possibility of extending such an approach to chaotic systems, both continuous and discrete, was established earlier Kautz:87 ; Beale:89 ; Grassberger:89 ; Hamm:91 . It was shown also that the presence of homoclinic tangencies, causing the fractalization of the basins, causes a decrease in the activation energy Soskin:01 .

In this Letter we show that a generic mechanism of fluctuational transition between co-existing CAs separated by an FBB does exist, that it is determined by a hierarchy of homoclinic original saddles forming the homoclinic structure, and that there is a unique most probable escape path (MPEP) from the CA that approaches an accessible orbit on the fractal boundary.

To demonstrate the existence of this escape mechanism, we take as an example the two-dimensional map introduced by Holmes Holmes:79 . The properties of this map, including the structures both of its CA and of its locally disconnected FBB, are generic for a wide class of maps and flow systems Cartwright:51 ; Yorke:85a . It is this fact, taken with the results of our investigations of escape in other systems, that allows us to conclude that the escape mechanism described below is indeed a typical one. The Holmes map is

$`x_{n+1}`$ $`=`$ $`y_n`$ (1)
$`y_{n+1}`$ $`=`$ $`bx_n+dy_ny_n^3+\xi _n,`$
where $`\xi _n`$ is zero-mean, white, Gaussian noise of variance $`D`$. Due to symmetry, the noise-free system (1) with $`b=0.2`$ and $`2.0d2.745`$ has pairs of co-existing attractors, the basins of which are separated by a boundary that may be either smooth or fractal, depending on the choice of parameter values. We choose for our studies $`b=\mathrm{\hspace{0.17em}0.2}`$ and $`d=\mathrm{\hspace{0.17em}2.7}`$, which corresponds to there being two co-existing CAs separated by a locally disconnected FBB (see Fig. 1). The fractal dimension of the boundary has been determined numerically (dim = 1.84472) by using the “uncertainty exponent” technique introduced in Grebogi:85 . The chaotic attractors in (1) appear as the result of a period-doubling cascade and, for the parameters chosen, each of them consists of two disconnected parts.

We have modelled the system (1) numerically, exciting it with weak noise, and have collected trajectories that include escape paths from one CA to the other, and also the corresponding realisations of noise that induced these transitions. By ensemble averaging a few hundred such escape trajectories and noise realisations, we have obtained the optimal escape path and the corresponding optimal force, which are shown in Fig. 2. The results of this statistical analysis allow us both to determine the boundary conditions near the CA and the FBB, and to demonstrate the uniqueness of the MPEP. It can be seen in particular that, in leaving the CA, the system (1) falls into a small neighbourhood of the saddle point of period 1 (S1) located between its two disconnected parts and having the multipliers $`\rho _1=0.118975`$ and $`\rho _2=1.681025`$. Its stable manifolds separate the parts of the CA, while the unstable ones belong to the CA. The system makes a few iterations in the neighbourhood of S1 (initial plateau in Fig. 2(a)) and then moves to the FBB in three steps, crossing it at a saddle point of period 3 (S3) with multipliers $`\rho _1=0.001016`$ and $`\rho _2=7.875512`$. Calculations have shown that, for the chosen parameter values, S3 lies within the FBB. Moreover, its stable manifold (solid black line in Fig. 1) is dense in the FBB and detaches the open neighborhood, including an attractor, from the FBB itself. This allows us to classify it as an accessible boundary point Grebogi:87 .

An analysis of the structure of escape paths inside the FBB has shown that the homoclinic saddle points play a key role in its formation. In the system (1), we observe an infinite sequence of saddle-node bifurcations of period $`3,4,5,6,7\mathrm{}`$, which occur at parameter values $`d_3<d_4<d_5<d_6<d_7\mathrm{}`$ and are caused by tangencies of the stable and unstable manifolds of the saddle point O at the origin. The homoclinic orbits appearing as a result of these bifurcations were classified earlier as original saddles, and it was also shown that their stable and unstable manifolds cross each other in a hierarchical sequence Grebogi:87 . To characterize this hierarchical relation between original saddles, it is reasonable to introduce a parameter $`\mu `$ equal to the ratio $`\mathrm{ln}(\rho _1(P))/\mathrm{ln}(\rho _2(P))`$, where $`\rho _1`$ and $`\rho _2`$ are the multipliers of a saddle point $`P`$. Calculations have shown that, for the original saddles of period $`3,4,5,6,7,8\mathrm{}`$ in (1), the following hierarchical sequence of index $`\mu `$ values occurs: $`\mu _3=3.339,\mu _4=3.08,\mu _5=2.999,\mu _6=2.339,\mu _7=1.958,\mu _8=1.539`$. Moreover, the values of index $`\mu `$ corresponding to the other homoclinic saddle cycles are close to zero. Correspondingly the probability of finding the system in their neighbourhood tends to zero.

These results allows us to infer features of a fluctuational transition through a locally disconnected FBB that are probably generic, as follows: (i) it always occurs through a unique accessible boundary point; and (ii) the original saddles forming the homoclinic structure of the system play a key role in the formation of the paths inside the FBB, the difference in their local stability defining the hierarchical relationship between them. Thus, we may claim that complicated structure of escape trajectories, caused by the thin homoclinic structure and their randomness, has in many respects a deterministic nature.

Having now understood the mechanism of escape, we can seek the MPEP. According to the Hamiltonian theory of fluctuations Kautz:87 ; Beale:89 ; Grassberger:89 ; Hamm:91  the MPEP is the path which minimizes the energy

$$S=\frac{1}{2}\underset{n=1}{\overset{N}{}}\xi _n^T\xi _n,$$
(2)
of the possible realizations of noise $`\{\xi _n\}`$ inducing a transition of the system (1) from the CA (with the initial condition on S1) to the FBB (with the final condition on the accessible orbit S3). The Lagrangian of the corresponding variational problem can be found following Grassberger:89  (cf. Dykman:90 ) in the form

$$L=\frac{1}{2}\underset{n=1}{\overset{N}{}}\xi _𝐧^𝐓\xi _𝐧+\underset{n=1}{\overset{N}{}}\lambda _𝐧^𝐓(𝐱_{𝐧+\mathrm{𝟏}}𝐟(𝐱_𝐧)\xi _𝐧),$$

where (1) is taken into account using the Lagrange multiplier $`\lambda _n`$. Varying $`L`$ with respect to $`\xi _n`$, $`\lambda _n`$, and $`x_n`$, the following area-preserving map is obtained:

$`x_{n+1}`$ $`=`$ $`y_n`$
$`y_{n+1}`$ $`=`$ $`bx_n+dy_ny_n^3+\lambda _n^y`$ (3)
$`\lambda _{n+1}^x`$ $`=`$ $`(d3x_{n+1}^2)\lambda _n^x/b\lambda _n^y/b`$
$`\lambda _{n+1}^y`$ $`=`$ $`\lambda _n^x`$
Equations (Fluctuational transitions through a fractal basin boundary) are supplemented by the following boundary conditions

$$\underset{n\mathrm{}}{lim}\lambda _n^y=0,(x_n^0,y_n^0)S1,(x_n^1,y_n^1)S3.$$
(4)
The MPEP is the minimum-energy heteroclinic trajectory connecting S1 to S3 in the phase space of (Fluctuational transitions through a fractal basin boundary). The solution of this boundary value problem is in general complicated, because of the presence of multiple local energy minima Luchinsky:02a  induced by the complex geometrical structure of the unstable manifolds of S3 (see e.g. Graham:84a  for a discussion). The solution of the boundary value problem involves a parameterization of the structure of the multiple local minima requiring, in turn, a proper parameterization of the unstable manifold in the vicinity of the initial conditions Beri:03 . The MPEP found by this method is shown in Fig. 2. It can be seen that, within the range shown by the vertical dotted lines in Fig. 2(a), the theoretical MPEP closely coincides with the path obtained by statistical analysis of escape trajectories in the Monte Carlo simulations. Note that no further action is required to bring the system to the other attractor once it has reached the accessible orbit of the FBB, i.e. once it has reached the points numbered 8 in Fig. 2(a) and (b); correspondingly, the optimal force measured in the numerical simulations (inset) falls back to zero.

The existence of an almost deterministic mechanism of transition across the FBB raises important questions about the effect of noise on this mechanism, and on the structure of escape trajectories inside the FBB. We have therefore used randomly chosen initial conditions in a very small neighborhood of the accessible point S3 through which escape occurs (see Fig. 2(b)). By definition, any arbitrarily small neighborhood of S3, lies within the FBB, and must contain points belonging to the basins of both attractors. Therefore the system can cross the FBB starting from a very small neighborhood of S3, even in the absence of noise. By collecting all such successful escape paths, we have calculated the probabilities for the system to pass via small neighborhoods of different original saddle cycles during its escape, both in the presence and absence of noise. As can be seen from Fig. 3, the corresponding probabilities demonstrate the same hierarchical interrelationship in both cases, which is determined by the value of index $`\mu `$ defined above. This structure is evidently robust with respect to noise-induced perturbations. The addition of noise causes a slight broadening of the distribution in Fig. 3 and a small increase of the probability for the system to escape via original saddles of larger period.

In conclusion, we have described the mechanism by which noise-induced escape occurs through a locally disconnected FBB. We have found the (unique) most probable escape path from a chaotic attractor to the fractal boundary, using both statistical analyses of fluctuational trajectories, and the Hamiltonian theory of fluctuations. We have shown that the original saddles forming the homoclinic structure play a key role in effecting the transition through the FBB itself. In particular, their local stability defines the hierarchical relationship between the probabilities for the system to pass via small neighborhoods of different original saddle cycles during its escape, both in the presence and in the absence of noise. We emphasize that the escape mechanism we have revealed must be applicable to the broad class of two dimensional maps and flows Holmes:79 ; Cartwright:51 ; Yorke:85a  that exhibit the same type of FBB. For instance, one possible application of our results is to the development of an energy-optimal control scheme for the CO2 laser, a discrete model of which demonstrates the type of FBB considered above Grigorieva .

The authors would like to thank Ulrike Feudel, Igor Khovanov and Suso Kraut for fruitful discussions. The work was supported by the Engineering and Physical Sciences Research Council (UK), the Russian Foundation for Fundamental Science, and INTAS.


text:
# On Einstein Causality and Time Asymmetry in Quantum Physics

## 1 Introduction

G. C. Hegerfeldt discovered some interesting features of the structure of quantum physics appertaining to microphysical causality . In particular, he showed that an initially localized particle immediately develops “infinite tails”, unless it remains localized for all times. This is a very broad and general result in that it does not depend on the details specific to a particular quantum system such as the form of the Lagrangian and type of interaction. Only the existence of a positive selfadjoint Hamiltonian, and therewith the symmetry of time translations, is assumed.

In its final form, the discovery of Hegerfeldt is encapsulated in the following theorem:

###### Theorem 1.1.

 Let $`H`$ be a selfadjoint operator, positive or bounded from below, in a Hilbert space $``$. For given $`\psi _0`$, let $`\psi _t,t`$, be defined as $`\psi _t=e^{iHt}\psi _0`$. Let A be a positive operator111It seems to us that Hegerfeldt’s proof holds only when $`A`$ is a bounded positive operator, though this is not explicitly stated in . in $``$ and $`p_A(t)`$ be defined as $`p_A(t)=\psi _t,A\psi _t`$. Then, either $`p_A(t)0`$ for almost all $`t`$ and the set of such $`t`$’s is open and dense, or $`p_A(t)=0`$ for all $`t`$.

It is a classical result due to Stone that every selfadjoint operator $`H`$ leads to a strongly continuous one parameter group of unitary operators $`U(t)=e^{iHt}`$. If the spectrum of the operator $`H`$ is bounded from below, then the mapping $`tU(t)`$ admits an analytic extension into the open lower half of the complex plane, $`^{}`$, i.e., the mapping $`zU(z)=e^{iHz}`$ is strongly analytic for every $`z`$ with $`\mathrm{}z<0`$. Further, since $`z_1+z_2U(z_1+z_2)=U(z_1)U(z_2)`$ and $`U(z)1`$, the analytic mapping $`zU(z)`$ furnishes a representation of the semigroup $`^{}`$ (under addition) by contractions in $``$. The proof of Theorem 1.1 is anchored in this analytic extension of $`e^{iHt}`$ into the semigroup $`e^{iHz},z^{}`$.

In the next section, we shall show that a result analogous to Hegerfeldt’s theorem holds for certain quantum systems the time evolution of which is given by a semigroup of bounded normal operators. Specifically, we shall prove that $`p_A(t)=\psi _t,A\psi _t0`$ for almost all $`t0`$, unless $`p_A(t)0`$ for all $`t0`$. The condition $`t0`$, rather than $`t`$ as in Hegerfeldt’s result, is a consequence of the fact that the time evolution is now furnished by a semigroup –not a unitary group as required for Theorem 1.1.

Semigroup time evolutions are generally understood as representing a time asymmetry or an irreversibility at the quantum physical level. The main theoretical question given rise to by our mathematical result is whether such asymmetry in time translations implies or is consistent with Einstein causality. The immediate emergence of “infinite tails” for an initially localized particle, a prediction of Hegerfeldt’s theorem, remains a feature of the particular semigroup time evolution developed in the next section. Thus, insofar as particle localization is concerned, if such localization for a relativistic particle is meaningful in the first place, it appears that asymmetry in time evolution and causality are quite distinct notions –in particular, the former does not ensure the latter.

However, to further examine the issues of Einstein causality for irreversible processes, we must also consider the space translations consistent with the semigroup time evolution (i.e., the set of spacetime translations resulting from all possible boosts of the time evolution semigroup). The characterization of relativistic quasistable states by irreducible representations of a particular subsemigroup of the Poincaré group, studied in Section 3, provides a concrete example with such a semigroup of spacetime translations. Within the context of this example, we discuss aspects of Einstein causality inferred by the main mathematical result of this paper.

## 2 A Variant of Hegerfeldt’s Theorem

The following variant of Theorem 1.1 holds:

###### Theorem 2.1.

Let $``$ be a Hilbert space, and $`H`$, a normal operator in $``$. Suppose the following conditions hold on $`\sigma (H)`$, the spectrum of $`H`$:

$`\underset{\lambda \sigma (H)}{sup}(\lambda _y)`$ $`=`$ $`k_0<\mathrm{}`$ (2.1)
$`\underset{\lambda \sigma (H),\lambda _x0}{sup}\left({\displaystyle \frac{\lambda _y}{\lambda _x}}\right)`$ $`=`$ $`k_1<0`$ (2.2)
$`\underset{\lambda \sigma (H),\lambda _x0}{inf}\left({\displaystyle \frac{\lambda _y}{\lambda _x}}\right)`$ $`=`$ $`k_2>0`$ (2.3)
where $`\lambda =\lambda _x+i\lambda _y`$.
Let $`e^{iHt},t[0,\mathrm{})`$, be the semigroup generated by $`iH`$, and let $`p_A(t)`$ be defined as

$$p_A(t)=\psi _t,A\psi _t$$
(2.4)
where $`A`$ is a bounded selfadjoint operator in $``$, and $`\psi _t=e^{iHt}\psi _0`$ for given $`\psi _0`$. Then,
either

$$p_A(t)0\text{for almost all}t[0,\mathrm{})$$
(2.5)
and such $`t`$’s are open and dense in $`(0,\mathrm{})`$
or

$$p_A(t)=0\text{for all}t[0,\mathrm{})$$
(2.6)
Proof. It is well known  that if $`\sigma (H)`$ satisfies (2.1), the operator $`(iH)`$ generates a strongly continuous one parameter semigroup of bounded normal operators $`e^{iHt},t[0,\mathrm{})`$. Conditions (2.2) and (2.3) ensure that the strongly continuous mapping $`te^{iHt}`$ admits an extension into the domain $`D=\{z=t+iy:t>0,k_1<y<k_2\}`$ such that it is strongly continuous on $`D`$ and strongly analytic on $`D(0,\mathrm{})`$.

The dual semigroup $`e^{iH^{}t}`$ admits an extension into $`D^{}=\{z=t+iy:t>0,k_2<y<k_1\}`$. As before, the extension is continuous on $`D^{}`$ and analytic on $`D^{}(0,\mathrm{})`$ in the strong operator topology. Therefore, the function $`p_A(t)=\psi _t,A\psi _t=\psi _0,e^{iH^{}t}Ae^{iHt}\psi _0`$ has an extension into $`DD^{}=\{z=t+iy:t>0,k<y<k\}`$, where $`k=\mathrm{min}\{|k_1|,k_2\}`$.

This extension, which we denote by $`p_A`$, is continuous on $`DD^{}`$ and analytic on $`(DD^{})(0,\mathrm{})`$. Further, since $`\overline{p_A(\overline{z})}=p_A(z)`$ for every $`z(DD^{})(0,\mathrm{})`$ and $`p_A(t)`$ is real for $`t(0,\mathrm{})`$, by the Schwarz principle of reflection, $`p_A`$ is analytic in $`(DD^{})`$. Since $`(0,\mathrm{})`$ is in the domain of analyticity of $`p_A`$, either (2.5) or (2.6) must hold. $`\mathrm{}`$

Remark: Notice that unlike in Theorem 1.1, the positivity of $`A`$ is not necessary for our proof.

## 3 Application to Relativistic Decaying States

Wigner’s pioneering work established that a correspondence exists between the unitary irreducible representations of the Poincaré group $`𝒫`$ and relativistic (stable) particles by way of their mass and spin. In the Hilbert space $`L_{mj}^2(^3)`$ of momentum wavefunctions for a particle of mass $`m`$ and spin $`j`$, the relevant unitary representation can be realized by

$$(U(\mathrm{\Lambda },a)f)(\stackrel{}{p},j_3)=e^{ip_\mu a^\mu }\underset{j_3^{}}{}D_{j_3^{}j_3}^j(W(\mathrm{\Lambda },p))f(\stackrel{}{\mathrm{\Lambda }^1p},j_3^{})$$
(3.1)
where $`(\mathrm{\Lambda },a)𝒫`$, $`\{\stackrel{}{p}\}=^3`$ and $`p^2=m^2`$. The $`W(\mathrm{\Lambda },p)`$ are the well known Wigner rotations. The $`D^j`$ matrices provide a $`(2j+1)`$-dimensional representation of the quantum mechanical rotation subgroup.

The unitary irreducible representation for the particle of mass $`m`$ and spin $`j`$ can be realized also in the Hilbert space $`L_j^2(^3)`$ of “four-velocity wavefunctions”:

$$(U(\mathrm{\Lambda },a)f)(\stackrel{}{u},j_3)=e^{imu_\mu a^\mu }\underset{j_3^{}}{}D_{j_3^{}j_3}^j(W(\mathrm{\Lambda },u))f(\stackrel{}{\mathrm{\Lambda }^1u},j_3^{})$$
(3.2)
where $`u_\mu =\frac{p_\mu }{m}`$. The crucial property that makes (3.2) possible is that Wigner rotations $`W(\mathrm{\Lambda },p)`$ are functions only of the quotients $`\frac{p_\mu }{m}=u_\mu `$, but not of the $`p_\mu `$ themselves . That is, $`W(\mathrm{\Lambda },p)=W(\mathrm{\Lambda },u)`$. The inner product under which the velocity wavefunctions form the Hilbert space $`L_j^2(^3)`$ is

$$(f,g)=\underset{j_3}{}\frac{d^3\stackrel{}{u}}{2u_0}\overline{f(\stackrel{}{u},j_3)}g(\stackrel{}{u},j_3)$$
(3.3)
where $`\{\stackrel{}{u}\}=^3`$ and $`j_3=j,j+1,\mathrm{},j`$.

It is generally believed that relativistic resonances and unstable particles are associated with simple poles of an analytic $`S`$-matrix. The location of the pole, $`𝗌_R=(Mi\mathrm{\Gamma }/2)^2`$, contains information about the “mass” and “width” of the particular quasistable state222We use the term “quasistable state” to refer to resonances and unstable particles collectively. while the specific partial wave in which the quasistable state occurs as a pole defines its spin $`j`$. Much like the unitary irreducible representations (3.2) of $`𝒫`$ which characterize isolated stable particles, we may define a representation for a quasistable state of (complex) square mass $`𝗌_R`$ and spin $`j`$ in the Hilbert space $`L_j^2(^3)`$ by

$$(U(\mathrm{\Lambda },a)f)(\stackrel{}{u},j_3)=e^{i\sqrt{𝗌_R}u_\mu a^\mu }\underset{j_3^{}}{}D_{j_3^{}j_3}^j(W(\mathrm{\Lambda },u))f(\stackrel{}{\mathrm{\Lambda }^1u},j_3^{})$$
(3.4)
where, as in (3.2), $`\{\stackrel{}{u}\}=^3`$ and $`j_3=j,j+1,\mathrm{},j`$.

Notice that the operators $`U(\mathrm{\Lambda },a)`$ defined by (3.4) are bounded in $`L_j^2(^3)`$ if and only if $`(\mathrm{\Lambda },a)𝒫_+`$, where

$$𝒫_+=\{(\mathrm{\Lambda },a):(\mathrm{\Lambda },a)𝒫,a_00,a^20\}$$
(3.5)
The subset $`𝒫_+𝒫`$ is a subsemigroup of $`𝒫`$ as it remains invariant under the product rule $`(\mathrm{\Lambda }_1,a_1)(\mathrm{\Lambda }_2,a_2)=(\mathrm{\Lambda }_1\mathrm{\Lambda }_2,a_1+\mathrm{\Lambda }_1a_2)`$ for $`𝒫`$. Further, the inverse mapping $`(\mathrm{\Lambda },a)(\mathrm{\Lambda },a)^1=(\mathrm{\Lambda }^1,\mathrm{\Lambda }^1a)`$ transforms every element $`(\mathrm{\Lambda },a)𝒫_+`$ (except $`(I,0)`$) out of $`𝒫_+`$. The $`𝒫_+`$ is the semidirect product of the group of orthochronous Lorentz transformations and the semigroup $`T_+`$ of spacetime translations into the forward light cone:

$$T_+=\{a:a_00,a^20\}$$
(3.6)
The operators defined by (3.4) provide a continuous representation of $`𝒫_+`$ by contractions in $`L_j^2(^3)`$. This representation is characterized by a real, discrete spin value $`j`$ and a complex square mass value $`𝗌_R=(Mi\mathrm{\Gamma }/2)^2`$, where $`M`$ and $`\mathrm{\Gamma }`$ may be interpreted as the mass and width of the quasistable state, respectively. Notice further that, as in (3.2), the Lorentz subgroup of $`𝒫_+`$ is unitarily represented by (3.4). It must be mentioned that complex mass representations with real velocities of the Poincaré transformations have been considered in the past .

The Hamiltonian $`H=P_0`$ for (3.4) acts in $`L_j^2(^3)`$ as

$$(P_0f)(\stackrel{}{u},j_3)=\sqrt{𝗌_R}u_0f(\stackrel{}{u},j_3),u_0=(1+\stackrel{}{u}^2)^{1/2}$$
(3.7)
This shows that $`P_0`$ is a normal operator whose spectrum $`\sigma (P_0)=\{(Mi\mathrm{\Gamma }/2)u_0:u_0[1,\mathrm{})\}`$. Thus, $`\sigma (P_0)`$ satisfies the conditions demanded in Theorem 2.1, and we have a semigroup of time evolution

$$(U(t)f)(\stackrel{}{u},j_3)=e^{i\sqrt{𝗌_R}u_0t}f(\stackrel{}{u},j_3)=e^{iM(1+\stackrel{}{u}^2)^{1/2}t}e^{\frac{\mathrm{\Gamma }}{2}(1+\stackrel{}{u}^2)^{1/2}t}f(\stackrel{}{u},j_3),t0$$
(3.8)
Now, if $`f_t=e^{iP_0t}f`$ for any $`fL_j^2(^3)`$ and $`A`$ is any bounded selfadjoint operator in $`L_j^2(^3)`$, then by Theorem 2.1, $`p_A(t)=f_t,Af_t`$ is either almost never zero for $`t0`$ or identically zero. In particular, if we choose for $`A`$ a projection operator $`N(V)`$ representing the localization of the quasistable particle in a finite volume $`V`$, a calculation can be carried out in complete analogy to that in  to show that there appear “infinite tails” for any $`t0`$. That is, the asymmetric time evolution of the quasistable state does not preclude the propagation of its survival probability into all space for arbitrarily small positive times. However, unlike Hegerfeldt’s theorem, ours does not yield non-zero probabilities for $`t<0`$ since the time evolution is now given by a semigroup with $`t0`$. The mathematical result notwithstanding, it is important here to mention that the spatial localization of a relativistic quasistable state, as given by the above projection operator $`N(V)`$, may not be a physically meaningful notion. It has been argued  that, even for a relativistic stable particle, the localization by way of the projection $`N(V)`$ considered in  does not have much physical content.

What may be relevant and more meaningful as a directly measurable quantity (unlike survival probability) is the decay probability of the quasistable state into its possible decay products. Thus, if the operator $`A`$ taken to be the projection representing the measurement of certain decay products, then (2.4) and (2.5) infer that the probability to detect decay products becomes non-zero immediately after the quasistable state is created at $`t=0`$. This property, however, does not say anything about the Einstein causality for the spacetime propagation of decay probability, a problem which requires that we consider the whole semigroup (3.6) of spacetime translations. To that end, recall first that the operators $`e^{iP_\mu a^\mu }`$ of (3.4) are bounded in $`L_j^2(^3)`$ only when $`aT_+`$. For such $`aT_+`$,

$$(U(a)f)(\stackrel{}{u},j_3)=e^{i\sqrt{𝗌_R}u_\mu a^\mu }f(\stackrel{}{u},j_3)=e^{i\sqrt{𝗌_R}(u_0t\stackrel{}{u}.\stackrel{}{x})}f(\stackrel{}{u},j_3)$$
(3.9)
where $`a=(t,\stackrel{}{x})`$.

Consider now a certain fixed $`\stackrel{}{x}`$. By (3.9), the $`U(a)`$ is a bounded operator only when $`t|\stackrel{}{x}|`$. Further, the mapping $`aU(a)`$ given by (3.9) admits an extension into $`D=\{z=t+iy:t>|\stackrel{}{x}|,y<\frac{\mathrm{\Gamma }}{2M}(t|\stackrel{}{x}|)\}`$ such that it is strongly continuous in $`D`$ and strongly analytic in $`D(|\stackrel{}{x}|,\mathrm{})`$. Starting from this observation, the proof technique of Theorem 2.1 can be invoked to conclude that, for any bounded operator $`A`$ in $`L_j^2(^3)`$ and given $`fL_j^2(^3)`$, the quantity $`p_A(t,\stackrel{}{x})=U(a)f,AU(a)f`$ is either almost never zero for $`t|\stackrel{}{x}|`$ or identically zero.

Either of these two possibilities is clearly causal, and so what remains to be examined is $`p_A(t,\stackrel{}{x})`$ in spacelike regions. Ideally, we expect a causal theory to predict vanishing expectation values $`p_A(t,\stackrel{}{x})=0`$ in all regions of spacetime outside the forward semi-lightcone $`T_+`$. While the theory presented here does not yield this ideal result, certain properties of the transformation operators $`U(a)`$ for $`aT_+`$ motivate the plausible argument that the theory is mathematically meaningful only for $`aT_+`$. In particular, the operators $`U(a)`$ are unbounded for $`aT_+`$, and for such $`a`$, the following results are true:

* There exist elements $`f`$ of $`L_{mj}^2(^3)`$ such that $`U(a)f`$ is infinite.
* For any $`fL_{mj}^2(^3)`$ and $`\alpha >0`$, there exists some $`aT_+`$ such that $`U(a)f>\alpha `$. That is, the function $`aU(a)f`$ grows without bound for any $`fL_{mj}^2(^3)`$.

If we want to maintain the general idea that the vectors $`f`$ and $`U(a)f`$ represent the same physical system as viewed by two different observers, then the transformations $`U(a)`$ for $`aT_+`$ need to be excluded from the theory. In particular, if $`U(a)`$ for $`aT_+`$ are admitted, then it implies that for every normalized state $`f`$ one can find an observer, spacelike separated from the first, to whom the quantity $`U(a)f`$ can be arbitrarily large. That is, there is no natural way to attribute a probability interpretation for $`U(a)f`$ when $`aT_+`$. Likewise, for a state $`f`$ that is initially normalized, i.e., $`f=1`$, the quantity $`p_A(t,\stackrel{}{x})=U(a)f,AU(a)f`$ cannot be interpreted as the expectation value of the observable $`A`$ in an arbitrary state $`fL_j^2(^3)`$ at $`a=(t,\stackrel{}{x})`$, if $`a𝒫_+`$. It may be possible333This was pointed out by a referee. that finite values can be restored if the expectation values are defined as $`p_A(t,\stackrel{}{x})=\frac{U(a)f,AU(a)f}{U(a)f,U(a)f}`$. This, however, implies that the state vector is re-normalized for each $`aT_+`$, and how such an operation can be given an unambiguous physical meaning is not clear.

## 4 Concluding Remarks

Theorem 1.1 proved by Hegerfeldt seems to imply, on the face of it, that if nature admits only positive (or bounded from below) energies, then its behavior is non-causal at the microphysical level. It may be possible to resolve this conflict, as the author suggests, once field theoretic concepts such as vacuum fluctuations and weak causality are introduced. Mathematically, the centrally significant component of Hegerfeldt’s result is the semi-boundedness of the Hamiltonian and the resulting unitary group which determines the time evolution of stationary quantum systems represented in a Hilbert space. It is this reversible group evolution which brings about the non-vanishing probabilities for (almost) all $`t`$, both positive and negative.

On the other hand, it has been argued that certain quantum mechanical processes such as decay exhibit asymmetric, irreversible time evolutions . Such asymmetric time evolutions can be described by semigroups. The main technical result we reported in this paper is a variant of Hegerfeldt’s theorem that applies for certain asymmetric time evolutions given by semigroups. The characterization of unstable particles by the complex mass representations of the semigroup $`𝒫_+`$ provides a concrete example where the new theorem can be applied. The semigroup representation (3.4) characterized by $`𝗌_R`$ and $`j`$, leads to non-zero probabilities $`p_A(t)`$ for (almost) all $`t0`$. Negative times are excluded from the prediction by virtue of the semigroup time evolution defined only for $`t0`$.

To examine the problem of Einstein causality for the spacetime evolution of the expectation values $`p_A(t,\stackrel{}{x})`$ of an observable $`A`$, it is necessary to consider, in addition to the time evolution semigroup, the set of space translations consistent with this semigroup (i.e., those obtained by boosting the asymmetric time translations). It was seen in Section 3 that $`p_A(t,\stackrel{}{x})`$ may be given a clean, unambiguous interpretation only in non-spacelike regions if the spacetime translations of the quantum system are defined by (3.4). Thus, while we have not obtained Einstein causality in the stronger sense, i.e., $`p_A(t,\stackrel{}{x})=0`$ for all $`aT_+`$, in the foregoing weaker sense, the characterization of the quasistable states by irreducible representations of the Poincaré semigroup appears to infer Einstein causality for the expectation values of decay products.

## Acknowledgment

We acknowledge the financial support from the Welch Foundation.


text:
# Comment on: ”Radiation-Induced ”Zero-Resistance State” and the Photon Assisted Transport”

## Abstract

Shi and Xie  predicted the N-shape current-voltage characteristic $`(CVC)`$ for a 2DEG in a magnetic field where the zero-resistance state has been observed in recent experiments. However it is known that in the absence of a magnetic field the zero resistance state (zero differential resistance state) is achieved in a system with the S-shape $`CVC`$. The difference in the behaviour of systems with N- and S-shape $`CVC`$ was studied more than three decades ago (see for example the review  and references therein) and is briefly explained in this Comment. At present it is not clear whether the N-shape $`CVC`$ may lead to the zero resistance state in a 2DEG in a magnetic field.

An interesting effect has been observed in recent papers : the resistance of a 2-DEG subjected to microwave irradiation drops to zero in some interval of an applied magnetic field $`B.`$ Possible mechanisms for this phenomenon are discussed in a number of papers . It was shown in Refs.  that in the presence of irradiation and a sufficiently strong magnetic field the conductivity $`\sigma _{xx}`$ may become negative in a weak electric field $`E_x`$. In addition, using a simple model, the authors of Ref.  have calculated the current-voltage characteristic $`I(V)`$ of an irradiated system with a density-of-states periodic in energy $`\epsilon `$. They have demonstrated that not only regions with negative conductance $`G`$, but also regions with positive $`G`$, may have negative differential conductance $`G_d=dI/dV`$ on the current-voltage characteristic ($`CVC`$) curve.

In the present Comment I would like to note that effects in systems with negative differential conductance (or resistance) depend in a crucial way on the type of $`CVC`$, i.e. whether it has the N- or S-shape form (see, for example, the review  and references therein). It was established long ago that the states corresponding to regions with negative $`G_d`$ on $`CVC`$ are unstable with respect to nonhomogeneous fluctuations. However the type of instability and final nonhomogeneous state which is formed as a result of the instability depend on the type of $`CVC.`$

In the case of the N-shape $`CVC`$ (three voltages $`V_i`$ correspond to one current $`I`$) a domain with a strong electric field (as happens, for example, in the Gunn effect) arises as a result of an instability in a homogeneous initial state (it is assumed that the total voltage is fixed). This domain moves in the direction of the applied electric field $`E=V/L`$ ($`L`$ is the length of the sample). In the presence of the high field domain the $`CVC`$ changes drastically: an almost flat region (plateau) appears on the $`I(V)`$ curve. The differential conductance $`G_d`$ (but not the resistance $`R_d=dV/dI`$) is zero on this plateau. The current $`I_{pl}`$ corresponding to this flat part of $`CVC`$ depends on the particular form of the $`CVC`$. In the case of a ”symmetric” $`I(V)`$ curve with absolute negative conductance (i.e. min{$`dI/dV`$} corresponds to $`V=0`$ ) this current is equal to zero (that is the total conductance is zero). This means that when the bias voltage varies in some limits, the current remains equal to zero. The $`CVC`$ obtained in Ref.  on the basis of a toy model can be regarded as a chain of N-shape $`I(V)`$ curves.

In the case of the S-shape $`CVC`$ (three currents $`I_i`$ correspond to one voltage $`V`$) small perturbations increase in the direction transverse to the current $`I`$ if this current corresponds to the part of $`CVC`$ with negative $`R_d`$. As a result of this instability, a domain of a strong current density arises in the sample (it is assumed that total current is fixed, otherwise the system goes over to a state on the stable part of the $`CVC`$). In the presence of this current domain (or a filament with a higher current density) the $`CVC`$ is modified so that an almost vertical part appears on the $`I(V)`$ curve ($`R_d`$ is close to zero). The electric field $`E_{ver}`$ corresponding to this part of the $`CVC`$ is again determined by the shape of $`CVC`$. All these statements concerning the S-type of $`CVC`$ are based on the study of the so called superheating mechanism of the S-shape $`I(V)`$ curve  . In the case of a ”symmetric” S-shape $`CVC`$ with negative $`G`$ (i.e. min{$`dV/dI`$} corresponds to $`I=0`$) the field $`E_{ver}`$ is zero (when the bias current varies in some limits, the voltage across the sample is zero ). Similar conclusions were made in a recent paper on a simple phenomenological model .

Of course the behaviour of a system with the N- or S-shape $`CVC`$ in the absence or presence of a magnetic field $`B`$ is different. In the latter case the Hall field $`E_y`$ arises for example if the transverse current $`I_y`$ is zero. The form of the $`CVC`$ of a homogeneous 2DEG (in the absence of domains of the electric field or current) depends on whether the Hall current $`I_y`$ or the Hall field $`E_y`$ is zero. For instance, if the $`I_x(V_x)`$ dependence corresponds to the N-shape $`CVC`$ in the absence of the Hall field $`E_y`$ $`(I_y0)`$, this correspondence changes completely in the absence of the Hall current $`I_y`$ ($`E_y0`$). It acquires a complicate form and can not be assigned to the N- or S-shape type of $`CVC`$ (several values of $`I_x`$ correspond to one $`E_x`$ and several values of $`E_x`$ correspond to one $`I_x`$). Therefore one has to analyse the type of instability and the form of the final state of the system separately for these two cases ($`I_y=0`$ or $`E_y=0`$).

To summarize, one can conclude that the absolute negative conductance obtained in Refs.  is not sufficient to explain the observed zero resistance (or differential resistance) states. In a system with negative conductance at low dc fields $`E,`$ the $`CVC`$ may have either the N- or S-shape form, and the behaviour of the system will be different in these both cases. I also would like draw attention to the fact that negative conductance was predicted long ago in Refs. (2DEG and 3DEG in a quantizing magnetic field in the presence of microwave irradiation) and in Refs. (superlattices under ac irradiation).

I would like to thank Sh.M.Kogan as well as B.Huckestein, N.Rick and J.Shi for useful discussions and comments.


text:
# The Discovery of Submillimeter Galaxies

## 1. Introduction

In 1997 the SCUBA camera was commissioned at the 15-m JCMT in Hawaii, providing 2.5-arcmin-wide images at wavelengths of 450 and 850 $`\mu `$m with resolutions of 7 and 15′′respectively. The MAMBO 1.25-mm camera at the 30-m IRAM telescope provides a similar capability (see Dannerbauer et al. 2002), while the 350-$`\mu `$m SHARC-II and 1.1-mm BOLOCAM cameras being commissioned at the 10-m Caltech Submillimeter Observatory (CSO) should provide significantly enhanced performance within a year. In the future, BOLOCAM is expected to observe from the 50-m mm-wave GTM/LMT telescope under construction in Mexico, a larger-format 8$`\times `$8-arcmin2 camera SCUBA-II is being designed in the UK, and developments of MAMBO are planned for APEX – a new 12-m telescope at the 5000-m ALMA site in Chile. APEX will join the 10-m-class Japanese ASTE submm telescope at the ALMA site for which a large bolometer camera is being designed.

SCUBA was the first submm camera able to survey fields large enough to detect the redshifted thermal dust emission from previously unknown galaxies (Smail, Ivison & Blain 1997). The peak of the emission from galaxies, typically at 60–200 $`\mu `$m in their restframe, corresponding to a range of dust temperatures between about 60 and 20 K, is redshifted into SCUBA’s observing bands. The steep submm Rayleigh–Jeans slope of the dust emission ensures that distant galaxies with similar bolometric luminosities and spectral energy distributions (SEDs) would produce similar flux densities at all redshifts from $`z0.5`$ to $`z5`$, assuming the same spectral energy distribution (SED). The detectability of galaxies does depend on the details of the SED, in general being greater for cooler temperatures at a fixed luminosity and redshift (see Blain et al. 2002).

More than 300 high-$`z`$ submm galaxies (SMGs) have now been detected by SCUBA and MAMBO, while BOLOCAM (Glenn et al. 1998) could detect new examples at a rate of order 1 per hour. Various types of surveys have been made: in narrow fields to exploit the gravitational lensing effect of clusters of galaxies (Smail et al. 2002; Chapman et al. 2002a; Cowie, Barger & Kneib 2002); in the Hubble Deep Field (Hughes et al. 1998); and in wider shallower surveys (Eales et al. 1999; Borys et al. 2002; Scott et al. 2002) covering a total of about 0.25 deg2.

These SMGs are responsible for a significant fraction of the star-formation/ AGN activity at $`z>1`$ (Blain et al. 1999a). It is vital to understand their individual properties if we are to understand the formation of galaxies as a whole. Their inferred space density is similar to that of giant Elliptical galaxies in the local Universe, and it has been suggested that they are high-$`z`$ formation events of these rare galaxies (Eales et al. 1999), presumably in the most massive dark-matter halos with the lowest specific angular momenta. It is more likely that they reflect a more common, short-lived phase involving the formation of galactic bulges in perhaps episodic mergers (Blain et al. 1999b). The test of these ideas is to measure the redshifts, clustering properties and mass distribution of a representative sample of SMGs.

## 2. Finding redshifts and studying astrophysics

The coarse ($``$15′′) resolution of submm images ensures that there are many possible faint optical counterparts to the detected galaxy. Hence, while some SMGs have bright optical counterparts (Ivison et al. 1998, 2000), most remain unidentified at optical wavelengths. Final confirmation of a correct identification in both position and redshift is provided by detecting (sub)mm-wave CO molecular line emission with a suspected redshift in the narrow spectral window of existing (sub)mm-wave spectrographs (Frayer et al. 1998, 1999).

Mm-wave continuum interferometer images of the fields have reduced positional uncertainties to $`1`$ arcsec for several SMGs (Downes et al. 1999; Frayer et al. 2000; Dannerbauer et al. 2002), but are very expensive in observing time. Deep near-IR imaging to $`K>22`$ tends to reveal plausible red counterparts for many SMGs (Smail et al. 1999; Frayer et al. 2003) by their $`JK`$ and $`IK`$ colors. Deep optical/near-IR spectroscopy of these galaxies can then be attempted to find redshifts. The most efficient technique for determining redshifts, however, appears to be to exploit very deep radio images. For reasonable SEDs, SMGs with 850-$`\mu `$m flux densities of about 5–10 mJy (with luminosities $`10^{13}`$ L) should be detectable in $`10`$$`\mu `$Jy RMS 1.4-GHz VLA images out to $`z3`$, if they lie on the far-IR–radio correlation observed for local galaxies (Condon 1992): (see Barger, Cowie & Richards (2000) and Chapman et al. (2001, 2002b). The wide field (0.25 deg2) and accurate sub-arcsec astrometry of VLA images, coupled with the low surface density of the faintest radio sources as compared with optically-selected galaxies yield accurate positions for a large fraction ($`70`$%) of the SMGs brighter than 5 mJy at 850 $`\mu `$m. This provides an opportunity for efficient multi-object optical spectroscopy. In March 2002, about 25 spectra were obtained using the Keck-LRIS spectrograph (Chapman et al. 2003). These will be subject to CO molecular line spectroscopy to confirm the identifications and to study both their gas masses (via velocity dispersions) and excitation conditions (via line–line and line–continuum ratios), using mm-wave interferometers in the 2002-2003 Northern winter. Only a handful of reliable redshifts were available for SMGs in 2001: now it is likely that a luminosity function of these radio-selected SMGs should be available in 2003.

## 3. Photometric redshifts

The key target of investigating the SMGs is now to find their physical properties, especially their masses. However, just obtaining a reasonably complete redshift distribution is important for fixing their form of evolution and ensuring their fractional contribution to the energy emission of all galaxies is correctly accounted for (Eales et al. 1999; Blain et al. 1999a).

As redshift surveys were generally unsuccessful until (Chapman et al. 2003), photometric techniques have been proposed to provide redshift information (Eales et al. 1999; Hughes et al. 1998, 2002). The key information available is radio (Carilli & Yun 1999), submm and mid-/far-infrared photometry, typically from VLA, SCUBA and ISO respectively. The IRAS survey is not sufficiently deep to detect SMGs; ISO data is deeper, but covers only a small area and is only useful for the closest (Soucail et al. 1999) or brightest examples (Ivison et al. 1998). The SIRTF space telescope will be ideal for finding far-IR SEDs.

Fitting a thermal spectrum to a galaxy at uncertain redshift leads to an unavoidable degeneracy between the inferred dust temperature $`T`$ and redshift $`z`$. The peak wavelength of the SED is determined in the observers frame, but this is shifted in exactly the same way by either a fractional increase in $`T`$ or a corresponding fractional decrease in $`(1+z)`$. This makes any far-IR/submm-based photometric redshift only as accurate as the knowledge of the temperature (Blain, Barnard & Chapman 2003), and not to $`\mathrm{\Delta }z0.5`$ as claimed by Hughes et al. (2002). Despite non-thermal radio emission being due to an entirely different process, the submm–radio properties of galaxies on the far-IR–radio correlation conspire to produce a similar $`T`$$`(1+z)`$ degeneracy if $`T<60`$ K (Blain 1999). If a reliable link exists between dust temperature and luminosity, then it is possible to break this degeneracy; however, the accuracy of the result is then determined by the scatter in the LT relation, which is likely to be at least 30%, implying at least this great an uncertainty in redshift (Blain et al. 2003).

The addition of $`K`$-band near-IR data (Dannerbauer et al. 2002) could also help, but first the intrinsic scatter in the ratio between the $`K`$-band and far-IR luminosity of the galaxies must be known. Based on observations of fairly complete samples of SMGs (Ivison et al. 1998; Frayer et al. 2003), the K-band magnitudes are certainly scattered by as much as $`\mathrm{\Delta }K2`$ mag. Spectroscopic redshifts remain essential for accurate study of SMGs.

## 4. Spectroscopic mm-wave redshifts

Correct identifications of SMGs via either samples of faint radio galaxies that feed targets to multi-object spectrographs or deep near-IR imaging and spectroscopy, must be confirmed and verified using mm-wave interferometers or single-antenna telescopes to detect CO molecular rotational line emission at integer multiples of 115 GHz in the galaxy’s rest frame.

There are other possibilities for obtaining spectroscopic redshifts, especially the direct detection of molecular or atomic fine-structure spectral lines at mid-/far-IR and submm wavelengths. The key is to obtain wide-band spectral coverage at these wavelengths, in order to allow searches for redshifts. At present, mm-wave spectrographs cover a total frequency range of only $`\mathrm{\Delta }\nu /\nu 0.05`$, and so a redshift must be known to well within 5% before confirmation can be made.

Powerful correlators for the 100-m GBT will allow searches for CO(1-0) lines (rest frequency 115 GHz) from high-redshift galaxies in the 22 & 44-GHz radio bands. At shorter submm wavelengths new very stable 230/345-GHz spectral line receivers at the CSO will have $`10`$ GHz bandwidths (Rice et al. in prep.). Dispersive techniques may allow very-wideband, low-resolution spectrographs at millimeter or far-infrared wavelengths. These include the ZSPEC and WaFIRS waveguide/grating concepts for space-borne and ground-based applications (Bradford et al. in prep).

### Acknowledgments.

I thank Vicki Barnard, Matt Bradford, Scott Chapman, Dave Frayer, Andrew Harris, Tom Phillips, Naveen Reddy, Frank Rice, Chip Sumner and Jonas Zmuidzinas. Work on a wide-band correlator for the CSO is being supported by the Caltech Tolman endowment and the Research Corporation.

## References

Barger, A.J., Cowie, L.L., & Richards E.A., 2000, AJ, 119, 2092
Blain, A.W. 1999, MNRAS, 309, 955
Blain, A.W., Barnard, V.E., & Chapman, S.C. 2003, MNRAS, 338, 733
Blain, A.W., Smail, I., Ivison, R.J., & Kneib, J.-P., 1999a, MNRAS, 302, 632
Blain, A.W. et al. 1999b, MNRAS, 209, 715
Blain, A.W., Smail, I., Ivison, R.J., Kneib, J.-P., & Frayer, D.T. 2002, Physics Reports, in press (astro-ph/0202228)
Borys, C., Chapman, S.C., Halpern, M., & Scott, D. 2002, MNRAS, 330, 63
Carilli, C.L., & Yun, M.S. 1999, ApJ, 513, L13
Chapman, S.C., et al. 2001, ApJ, 548, L147
Chapman, S.C., Scott, D., Borys, C., & Fahlman, G.G. 2002a, MNRAS, 330, 92
Chapman, S.C., et al. 2002b, ApJ, in press
Chapman, S.C., et al. 2003, Nature, in press
Condon, J.J., et al. 1992, ARA&A, 30, 575
Cowie, L.L., Barger, A.J., & Kneib, J.-P., 2002, AJ, 123, 2197
Dannerbauer, H., et al. 2002, ApJ, in press (astro-ph/0201104)
Downes, D. et al. 1999, A&A, 347, 809
Eales, S., et al. 1999, ApJ, 515, 518
Frayer, D.T., et al. 1998, ApJ, 506, L7
Frayer, D.T., et al. 1999, ApJ, 514, L13
Frayer, D.T., Smail, I., Ivison, R.J., & Scoville, N.Z. 2000, AJ, 120, 1668
Frayer, D.T., et al. 2003, AJ, submitted
Glenn, J., et al. 1998, Proc. SPIE, 3357, 326
Hughes, D.H., et al. 1998, Nature, 394, 241
Hughes, D.H., et al. 2002, MNRAS, 335, 871
Ivison, R.J., et al. 1998, MNRAS, 298, 583
Ivison, R.J., et al. 2000, MNRAS, 315, 209
Scott, S.E., et al. 2002, MNRAS, 331, 817
Smail, I., Ivison, R. J., & Blain, A. W. 1997, ApJ, 490, L5
Smail, I., et al. 1999, MNRAS, 308, 1061
Smail, I., et al. 2002, MNRAS, 331, 495
Soucail, G., et al. 1999, A&A, 343, L70


text:
# The influence of hydrogen adsorption on magnetic properties of Ni/Cu(001) surface

## I Introduction

Ultrathin ferromagnetic films grown on nonmagnetic substrate show peculiar magnetic behavior bluegel . One of the unique phenomena which is observed in the ultrathin Ni films on Cu(001) substrate is the spin direction reorientation transition (SRT) from an in-plane to out-of-plane magnetization direction when the Ni-overlayer thickness is increased babe ; kirs99 ; kirs00 . It is very important for the spintronic magnetic device applications and its microscopic understanding attracted recently both experimental and theoretical uibe99 ; shic99  interest. The gas adsorption on magnetic film provides the way to monitor the surface magnetic properties resulting in strong influence on the SRT critical thickness $`d_c`$. For Ni/Cu(001) films, the hydrogen adsorption is observed to reduce $`d_c`$ by $``$ 4 monolayers (ML) from its value of $``$ 11 ML for the vacuum/Ni/Cu-films kirs99 ; kirs00 .

The key quantity which drives the SRT in ultrathin magnetic films is the magneto-crystalline anisotropy energy (MAE). It determines the preferred film magnetization orientation by minimizing the free energy of the system. For the free-standing Ni/Cu(001) films, the $`d_c`$ is found to be determined uibe99 ; shic99  by competition of the uniaxial MAE due to the tetragonal distortion of Ni-films (so called “volume” MAE) and the surface MAE. The adsorption of a gas produces non-trivial changes in geometrical structure as well as magnetic properties of Ni films wimm85 ; wein85  influencing both volume and surface MAE, and results in change of the SRT critical thickness.

In recent years it has become possible to make use of ab initio density functional theory to predict the MAE in ultra-thin films jansen . From theoretical and computational point of view, the MAE calculations are extremely difficult due to the very high energy resolution in the range of few $`\mu `$eV which is required. For the H/Ni/Cu(001) films, the problem becomes even more complex due to the need of accurate account of structural relaxation. It requires us to use the most accurate total energy full-potential linearized augmented plane wave method (FP-LAPW) wimm85  in order to account simultaneously on equilibrium geometrical and magnetic structures.

## II Method and Results of Calculations

The experiments show that H2 adsorbs dissociatively on Ni(001) surface in fourfold hollow sites sten85 . Therefore, as a structural model for ab-initio calculations we use free-standing Ni-films (H/Nid/H, $`d`$ = 1-11) and Ni-films on Cu-substrate (H/Nid/Cu7/Nid/H, $`d`$ = 1–6) with ordered $`p(1\times 1)`$H adsorbate overlayer as shown in Fig. 1.

The in-plane experimental lattice constant of Cu $`a_{Cu}=3.615`$ Å LB  is used which remains unchanged in the calculations. All interlayer spacings $`d_{ij}`$ (see Fig. 1) are relaxed to the equilibrium values. In order to determine the equilibrium slab geometry we employ the FP-LAPW method in FLAIR implementation (unpublished improved and rewritten version of the original FLAPW codes wimm85 ). The scalar-relativistic atomic force technique is employed for the total energy minimization. Here, the Perdew and Wang pw91  approximation for exchange-correlation potential is used with the plane wave energy cutoff $`E_{cut}`$ of 13 Ry, and the 28 special k-points in the 1/8th irreducible part of 2-dimensional Brillouin zone (2D BZ) are used for the BZ integrations. The convergence better than 1 $`\times `$ 10-6 e/(a.u.)3 is achieved for charge/spin densities, and better than 0.2 mRy/a.u. for the atomic force acting on individual atom.

The relaxed interlayer distances for different free-standing films (H/Nin/Cu7/Nin/H, $`n`$ = 3–6) are shown in Fig. 2 and compared with the relaxation of tetragonal $`p(1\times 1)`$H Ni(100) surface. Clearly, we can distinguish the outward relaxation of top Ni layer as well as the strong influence of Ni/Cu interface on inter-layer distances in the close vicinity of the interface. The relaxation of Ni/Cu subsurface is different only for 5 ML Ni. The deeper interlayer distances (thicker film) approach the bulk value.

The equilibrium values of inter-layer distances for two model systems: eleven-layer tetragonal Ni-film covered with H and seven-layer Cu-film with overlayer of six Ni layers covered with H are shown in Table 1. Both interfaces (H/Ni, Ni/Cu) influence strongly the interlayer distances in the overlayer slab. Note, that the calculated bond length $`d_{HNi}`$ of 3.48 a.u. is slightly shorter than the bond length for the bulk nickel hydride (3.52 a.u.). The outward relaxation of top Ni-layer is found, the H-ML influences strongly only the interlayer distance between the Ni-surface and Ni-subsurface layer. Below the interface, the interlayer distances oscillate around the bulk value for strained tetragonal Ni bulk ($`d_{}=3.20`$ a.u.). The electron screening in the metal is responsible for fast damping of these oscillation.

The layer-resolved spin magnetic moments $`M_s`$ for these systems are shown in Table 2. These values correspond to the magnetic moments in the Ni “muffin-tin” spheres ($`R_{MT}=2.2`$ a.u.). It is seen that there is strong reduction of the spin magnetization for the top Ni-layer due to the interaction with the H-adlayer. The strong hybridization of the H s state with the filled majority Ni d band changes the band structure and the surface density of states (SDOS) leading to the decrease of spin-majority and increase of the spin-minority SDOS. Away from the interface, the spin magnetic moments are slowly converging to their bulk values. The layer-resolved magnetic moments for different H/Nin/Cu7-films ($`n`$ = 3–6) on Cu-substrate are shown in Fig. 3. We note that for very thin H/Nin/Cu7-films with n=1,2 the Ni-local moments disappear and the system becomes non-magnetic. With the increase of the Ni-film thickness, the local Ni $`M_s`$ of 0.24 $`\mu _B`$ is formed at the Ni-interface and then increases away from the H/Ni interface. When approaching the Ni/Cu interface, the Ni-atom magnetic moments start to decrease again and become $``$ 0.45 $`\mu _B`$ for Ni/Cu interface layer (depending of the Ni-film thickness).

The anisotropic energy density of a tetragonal ferromagnetic film is written as Landau :

$`E/V=K_1^vm_z^2K_2^vm_z^4K_3^vm_x^2m_y^2`$
$`{\displaystyle \frac{2}{d}}(K_1^sm_z^2+K_2^sm_z^4+K_3^sm_x^2m_y^2),`$
where, $`K^v`$ terms are 2nd and 4th-order volume-type anisotropy constants, $`K^s`$ terms are surface/interface-type 2nd and 4th-order anisotropy constants, $`m_{x,y,z}`$ are magnetization cosines with respect to the crystal axes, and $`d`$ is a thickness of magnetic film. We assume that the 4th-order terms in Eq.(1) are significantly smaller than 2nd-order uniaxial anisotropy constant $`K_1^v`$, and neglect them. The MAE then can be characterized with the difference in the total energy when magnetization is oriented along  (in-plane, $``$) and  (out-of-plane, $``$) axes (MAE = $`E[100]E[001]`$).

We use the relativistic version flapwso  of the FP-LAPW method to solve self-consistently the Kohn-Sham-Dirac equations with spin-orbit coupling included to obtain the ground state charge and spin densities for the magnetization directed along -axis. The fully optimized interlayer distances as discussed above are used for free-standing Ni-films and Ni-films on Cu-substrate with H-adlayer. The MAE is obtained by applying the force theorem to the spin - axis rotation Igor : from the self-consistent ground state charge and spin density obtained for the  spin axis, a calculation of the band structure for  spin axis orientation is performed, and difference of the single particle eigenvalue sums is then taken to be the MAE. For the MAE calculations the k-points mesh equivalent to 6400 k-points in the full 2D BZ is used guaranteeing the MAE convergence better than 10 $`\mu `$eV.

The MAE as a function of the Ni-film thickness for the H/Nid/H films with $`d`$ = 3,5,7,9 ML is shown in Fig. 4. For the $`d=3`$ we found small and negative MAE (-0.02 meV) keeping the Ni-film magnetization in  plane. It becomes positive for the $`d=5`$ (0.432 meV) resulting in the out-of-plane magnetization switching. The MAE is positive with further increase of $`d`$ and shows pronounced oscillations. The linear interpolation MAE = $`dK_V+\mathrm{\hspace{0.33em}2}K_I`$ yields the estimates for the “volume” $`K_V`$= 0.106 meV/atom and “interface” $`K_I`$=-0.127 meV/atom MAE contributions. The calculated “volume” MAE agrees well with $`K_V`$= 0.0835 meV/atom calculated for free-standing unrelaxed Ni-films without H-adlayer shic99  and extrapolated to T=0 K experimental value of 0.072 meV/atom of Ref. farle . The H/Ni “interface” $`K_I`$ is calculated to be substantially smaller than vacuum/Ni “surface” $`K_S`$ = -0.447 meV/atom shic99  for free-standing Ni-films and extrapolated to T=0 K experimental value of -0.7 meV/atom kirs99 . We note that for the very thin Ni-films considered here, validity of the linear fit for the MAE and its separation into “volume” and “surface” contributions is not well justified due to the strong dependence of the Ni-film magnetic properties on the thickness of the film, especially in the presence of the H-adlayer.

In order to estimate the $`d_c`$ critical thickness for SRT, we take into account the shape anisotropy energy (SAE) due to the magnetic dipole interaction, which provides additional in-plane anisotropy. This anisotropy is estimated using the relation SAE = $`2\pi M^2`$ to the spin magnetization density M (in CGS units). It yields for the SAE the values of -0.016 meV, -0.032 meV, -0.055 meV and -0.073 meV corresponding respectively to the H/Nid/H films with $`d`$ = 3,5,7,9 layers (which are consistent with experimentally derived SAE/Ni atom of -0.0075 meV kirs99 ). This additional negative SAE slightly shifts the MAE in Fig. 4 downwards yielding the $`d_c`$ = 4.

Indeed, the H/Nid/H model is a way too simple to describe quantitatively the SRT in in Ni/Cu layers with H-adlayer. As it was already mentioned in Ref. kirs99 ; uibe99  the Ni/Cu interface can play an important role. In order to evaluate the influence of Ni/Cu interface on SRT we perform the calculations for H/Nid/Cu5/Nid/H ($`d`$ = 3,4) films, where the 5 Cu layers play a role of the substrate. The use of 5 ML of Cu layers instead of 7 ML allows to reduce computational effort without producing any significant impact on the magnetic properties and the MAE of H/Ni/Cu-films, since they are originated from Ni magnetic film and not the Cu non-magnetic substrate, and the relaxed Cu-interlayer distance below the Ni-interface is found to be very close to its bulk value of 3.427 a.u. For the $`d=3`$ we found negative MAE of -0.192 meV and SAE of -0.036 meV keeping the Ni-film magnetization in  plane. Already for $`d=4`$ the MAE becomes positive 0.543 meV while the SAE is small and negative -0.052 meV. Again as in the case of free-standing H/Nid/H films we get the SRT critical thickness $`d_c=4`$.

To summarize, we found that hydrogen adsorption for the Ni films on Cu substrate yields the reduction of the SRT critical thickness. The calculated $`d_c=4`$ for both H/Nid/H and H/Nid/Cu5/Nid/H films agrees well with the experimental $`d_c`$ of 7 ML kirs99 . We show that the MAE has strong and oscillatory dependence on the Ni-film thickness which deviates substantially from the linear fit. We attribute the decrease of the $`d_c`$ due to the H adsorption to originate from strong reduction of the magnitude for the Ni “surface” MAE contribution. In turn, this is caused by strong decrease of the exchange splitting at the H/Ni interface due to the strong hybridization of the H $`s`$ state with the Ni-bands.

### Acknowledgment

The financial support was provided by the Academy of Sciences of the Czech Republic (Grant No. A1010214), CMS Vienna (GZ 45.504), and by the RTN project ”Computational Magnetoelectronics” of the European Commission (HPRN-CT-2000-00143).