text:
# Total Cross-sections

hep-ph/0205197

IISc/CTS/04-02

Total Cross-sections111Talk given by G.Pancheri at Photon2001,Ascona, Switzerland
R.M. Godbole1, A. Grau2, and G. Pancheri3

1. Centre for Theoretical Studies, Indian Institute of Science, Bangalore, 560012, India. E-mail: rohini@cts.iisc.ernet.in

2. Centro Andaluz de Física de Partículas Elementales and Departamento de Física Teórica y del Cosmos, Universidad de Granada, Spain. E-mail igrau@ugr.es

3. INFN Frascati National Laboratories, Via E. Fermi 40, I00044 Frascati, Italy. E-mail:pancheri@lnf.infn.it

Abstract

We examine the energy dependence of total cross-sections for photon processes and discuss the QCD contribution to the rising behaviour.

>

## Abstract

We examine the energy dependence of total cross-sections for photon processes and discuss the QCD contribution to the rising behaviour.

A look at total cross-sections for the processes $`pp,p\overline{p},\gamma p,\gamma \gamma hadrons`$ immediately raises a number of questions, like: what gives the energy dependence of total cross-sections? Are photon data properly normalized? Are the predictions from factorization, quark counting and VMD, consistent with the complete set of data available in the same energy range?

In this talk we describe work in progress towards a QCD Description of the energy dependence of total cross-sections. The issue has both a theoretical and a practical interest, as it is necessary to have a reliable model to predict total hadronic cross-sections from $`\gamma \gamma `$ collisions, which form a bulk of the hadronic backgrounds at the Linear Colliders, in order that these are properly evaluated. Indeed, convoluting the photon spectrum with various predictions for $`\gamma \gamma hadrons`$, one finds that those for $`e^+e^{}e^+e^{}hadrons`$ differ by 30$`40`$%. In order to reduce this uncertainty, it is necessary to drastically reduce the range of variability present in $`\gamma \gamma `$ collisions, where models can differ by more than a factor two in their predictions for the total cross-section. These differences are due to those in the absolute normalization and the slope with which the total cross-section rises in these models, all being consistent with the current data.

In general the task of describing the energy behaviour of total cross-sections can be broken down into three parts: i) the rise, ii) the initial decrease, iii) the normalization. The rise alone can be obtained

* in the Regge-Pomeron model, with $`\sigma _{total}=Xs^ϵ+Ys^\eta `$, through $`s^ϵ`$, although it does not seem that the same power $`ϵ`$ fits protons and photons: one finds $`ϵ_{pp}=0.08,ϵ_{\gamma \gamma }=0.10.2`$. To overcome this problem, it has been suggested to add more power terms, thus increasing the number of free parameters.
* from factorization, but there remain the problem of getting the proton-proton cross-section from first principles
* using the QCD calculable contribution from the parton-parton cross-section, whose total yield increases with energy
* a combination of the above two

In the Minijet Model, the rise is driven by the LO QCD contribution to the integrated jet cross-section

$$\sigma _{jet}=_{p_tmin}\frac{d^2\sigma _{jet}}{d^2\stackrel{}{p}_t}d^2\stackrel{}{p}_t=\underset{partons}{}_{p_{tmin}}d^2\stackrel{}{p}_tf(x_1)𝑑x_1f(x_2)𝑑x_2\frac{d^2\sigma ^{partons}}{d^2\stackrel{}{p}_t}$$

which depends on the densities and very dramatically on $`p_{tmin}`$, the minimum transverse momentum cut-off. To ensure unitarity, the mini-jet cross-sections are embedded into the eikonal formulation, which gives the Eikonal Minijet Model in LO QCD (EMM)

$$\sigma _{pp(\overline{p})}^{\mathrm{inel}}=2d^2\stackrel{}{b}[1e^{n(b,s)}],\sigma _{pp(\overline{p})}^{\mathrm{tot}}=2d^2\stackrel{}{b}[1e^{n(b,s)/2}cos(\chi _R)]$$

In the EMM, one puts $`\chi _R=0`$. To proceed further, one can separate the non perturbative from the perturbative behaviour, with $`n(b,s)=n_{NP}(b,s)+n_P(b,s)`$, and then factorize b vs. s behaviour. The simplest model has $`n(b,s)=A(b)[\sigma _{soft}+\sigma _{jet}]`$.

Taking the matter distribution A(b) to be the convolution of the Fourier transform of the form factors of the colliding particles, the s-dependence is then entirely contained in $`\sigma _{soft}`$, parametrized so as to reproduce the low-energy data, and $`\sigma _{jet}`$, which is given by the LO QCD jet cross-sections.

The consistency between $`\gamma p`$ and $`\gamma \gamma `$ can be studied by applying the EMM model with same set of parameters to the relevant data. The total cross-section predictions for photon processes in the EMM model include the probability $`P_{had}`$ for the photon to behave like a hadron, a probability expressed through a parameter obtained using VMD, $`P_{had}1/240`$. With the EMM for the $`\gamma p`$ total cross-sections, using for A(b) the convolution of proton (dipole) and pion-like (monopole with scale $`k_0`$) formfactor, one obtains a band of values symmetrically encompassing all the data. We then apply the same formalism and the same parameters to the $`\gamma \gamma `$ case, and find the band shown in Fig.(1), which spans all the data, but with the lower bound slighly below the data, especially at low energies. We also see that the Aspen model prediction, obtained using factorization, is clearly lower than the data.

The comparison between $`\gamma p`$ and $`\gamma \gamma `$ indicates the existence of a problem in the normalization of $`\gamma \gamma `$ data, first noticed in. Indeed one can see that using VMD and Quark Counting to put proton and photon data on same scale, $`\gamma p`$ falls in place, $`\gamma \gamma `$ data remain higher than the rest, basically the same result suggested by the EMM model. From these considerations, it would appear that data for $`\gamma \gamma `$ total x-section are overestimated by about 10%. We also notice that the normalization problem can confuse the issue of the rise.

Further refinements of the minijet model are possible, using soft gluon summation to include initial state acollinearity among partons. The model proposed to do this introduced an energy dependence in the impact parameter distribution, namely

$$n(b,s)=A_{soft}(b)\sigma _{soft}+A_{PQCD}(b,s)\sigma _{jet}^{LO}$$

with $`A_{PQCD}(b,s)`$ given by the Fourier transform of the transverse momentum distribution of the initial parton pair, due to initial state soft gluon radiation. Using the QCD resummation techniques, this leads to

$$A_{PQCD}(b,s)\frac{e^{h(b,s)}}{d^2\stackrel{}{b}e^{h(b,s)}},withh(b,s)=_{k_{min}}^{k_{max}}d^3n_{gluons}(k)[1e^{i\stackrel{}{k}_t\stackrel{}{b}}]$$

$`k_{max}`$, which is energy dependent, can be taken to be the kinematic limit, averaged over the parton densities, while $`k_{min}=0`$. The difficulty in using $`k_{min}=0`$ stems from our ignorance on $`\alpha _s(k_t)`$ as $`k_t0`$ . To proceed further one needs to make models for this behaviour. Our model uses a singular but integrable parametrization for $`\alpha _s`$ in the infrared limit. This introduces a strong energy dependence in the impact parameter distribution, physically understandable as follows. As the energy increases, one probes smaller and smaller $`k_t`$ values. The more singular $`\alpha _s`$ is, the more is the emission of soft gluons making the initial partons more acollinear resulting in loss of parton luminosity and a decrease in the jet cross-section. This effect is what one might call the taming of the rise. We show in Fig.(2) a preliminary result with GRV densities and $`p_{tmin}=2GeV`$.


format: text

text:
# Current-induced highly dissipative domains in high 𝑇_𝑐 thin films

## I introduction

It has been observed in the early 70’s that in one-dimensional (1D) classical superconductors quantum effects can dominate the dissipation and give rise to spectacular phenomena such as voltage steps in the I-V curves . The steps are due to the sequential development of phase-slip centers (PSC) , which are dissipative regions, where the superconducting order parameter oscillates. Each time it drops to zero the phase of the order parameter ’slips’ by 2$`\pi `$ . The PSC formation could be observed only in a very narrow interval of temperature below $`T_c`$ . At lower temperature, sharp features in the I-V characteristics were also observed, but caused by a different mechanism: thermal instabilities resulting in the development of localized normal regions where the temperature exceeds $`T_c`$ (hot spots) .

Do high $`T_c`$ compounds behave similarly? Signatures of phase slippage as been reported near $`T_c`$ in YBa2Cu3O7-δ grain boundary junctions , but also at low temperature in YBa2Cu3O7-δ films  and ceramics , suggesting that the region where the dissipation is dominated by the nucleation of PSC is not restricted to the vicinity of $`T_c`$. In addition, these experiments would imply, if the PSC interpretation is correct, that the PSC picture, developed for a 1D situation, is still valid in a 2D superconductor.

Nevertheless only few publications interpret the steps in the I-V curves at low temperature in terms of PSC. Discontinuities in the I-V characteristics are generally attributed to heating effects (or ’bolometric’ effects) because thermal instabilities can produce very fast voltage variations in thin films. Indeed the thermal response time, which is equal to the heat capacity of the film divided by the thermal conductance of the film-to-substrate boundary, is typically a few nanoseconds.

Some reports however do not interpret the voltage jumps induced by the application of a high current as hot spots, but the proposed explanations are all different from the PSC nucleation scenario.

At low temperature, the steps observed in Nd2-xCexCuOy were explained in terms of a strongly energy dependent density of states . In this picture, unlike the PSC model, these steps do not reflect the nucleation of a spatially limited domain surrounded by low dissipative regions, but are caused by an oscillation of the density of state as a function of energy in the superconducting mixed state. This model is valid only in the clean or superclean limit and therefore cannot be a general model for all superconductors. In contrast, steps attributed to spatially limited domains were reported in narrow YBa2Cu3O7-δ bridges and interpreted as an effect of edge pinning .

Close to $`T_c`$, the application of a high current density has been seen to induce large voltage jumps bringing the sample close to the normal state resistance in YBa2Cu3O7-δ  and in Bi2Sr2CaCu2O8+δ  thin films. These discontinuities in the I-V characteristics were attributed to flux flow instabilities, an effect predicted 25 years ago by Larkin and Ovchinnikov , involving the density of states in the vortex core.

It is likely that a single mechanism is at the origin of most of the abrupt voltage variations observed in high $`T_c`$ superconducting films at high current density. In an effort to provide an overall picture of this extreme regime, we have performed time resolved transport measurements on two different compounds, Bi2Sr2CaCu2O8+δ (BSCCO) and Nd1+xBa2-xCu3O7-δ (NBCO) thin films, over a wide temperature range.

Based on our measurements, we show that a very systematic and general phenomenon appears in HTS compounds at a certain critical current $`j^{}(T)`$, causing a rapid development of highly dissipative domains. We argue that, unlike what is concluded for conventional superconductors, even at low temperature these domains do not result from a thermal instability. Instead they behave just like 1D PSC’s observed in metallic superconductors near $`T_c`$.

We also show that the current induced breakdown of superconductivity seems to result from a large scale propagation starting with an elementary PSC. A full understanding of the way the current induces a switch into the normal state is important in the development of certain applications, in particular in the use of HTS thin films as fault current limiters , bolometric detector  or ultrafast superconductive switches .

## II experimental details

The samples used for the present study are BSCCO and NBCO epitaxial films, 200 nm thick, deposited by rf magnetron sputtering on (100) oriented MgO substrates. The NBCO film are Nd-rich (the nominal target composition is Nd1.15Ba1.85Cu3O7-δ). This prevents the formation of screw dislocations, resulting in good cristallinity , but $`T_c`$ decreases slightly . MgO was chosen for its high thermal conductivity (about 25 times larger than SrTiO3 and LaAlO3) in order to optimize the heat removal during the measurements. As revealed by x-ray analysis, all films were $`c`$-axis oriented. 60 $`\mu `$m wide bridges were patterned by standard photolitography into a four-probes configuration, with a distance of 500 $`\mu `$m between the voltage contacts. A gold layer was deposited on all the contact pads to reduce their resistance. The critical temperatures after the patterning, taken at mid-transition, are 80 K for the NBCO and 76 K for the BSCCO thin film.

To limit the Joule heating, the current is applied by single pulses. The purpose-built pulsed current source could deliver square or triangular pulses with an amplitude up to 1 Amp and a duration ranging from 10 $`\mu `$s to 1 s. During square pulses, the current varied by less than 1$`\%`$. Current stability is crucial since any spike above the critical value can result in an abrupt resistance change, leading to wrong interpretations. The voltage is recorded as a function of time by a Tektronics digital oscilloscope. The current is measured simultaneously on the second channel of the oscilloscope by picking the voltage across a 1 $`\mathrm{\Omega }`$ resistor in series with the sample.

Resistive measurements were performed from 4 K up to $`T_c`$. The sample is glued by a GE heat conducting varnish onto a copper block and surrounded by 100 mbar of helium gas. The sample temperature is controlled by a heater and monitored by a carbon-glass probe placed inside the copper block. Once a pulse is applied, the heat front crosses the substrate and reaches the copper block in about 40 $`\mu `$s. In the dissipative state (normal or mixed state), the heat produced in the film causes a smooth resistance increase, consistent with a simple simulation of heat diffusion through the substrate to the copper block. The thermometer probe is however too far to record the real time temperature increase of the film. Therefore we use the evolution of the film resistance to obtain information about its temperature during a pulse. A magnetic field generated by a superconducting magnet in a persistent mode can be applied parallel to the $`c`$-axis.

## III results

### A critical current

The currents applied in this experiment are well above the critical current which is traditionally defined at the onset of resistance. At high current density (above 106 A/cm2), the measured voltage across the bridge exhibits abrupt variations as a function of current or temperature. Two kinds of voltage jumps are observed: small steps where the resistance increase is less than 1 % of $`\rho _n`$ ($`\rho _n`$ being the normal resistivity just above $`T_c`$), and large jumps where the resistance reaches a value close to $`\rho _n`$. The voltage increase during a small jump is instantaneous at our time resolution, i.e. the voltage variation is faster than the microsecond rise time of our differential amplifier. With constant current pulses, we show in this section that both kinds of voltage steps are due to a single mechanism, which occurs when the current reaches a critical value that we call $`j^{}(T)`$.

A typical set of results for identical square pulses is shown in fig. 1 for the NBCO thin film at 0.1 T. The applied current density is 1.9 106 A/cm2 and each curve corresponds to a different sample temperature, measured just before the pulse application, varying from 60.0 K to 62.1 K. Each curve exhibits two clear steps before the total breakdown of superconductivity occurs, with the resistance rising up to $`\rho _n`$. Note that each particular jump occurs at a given resistance independent of the initial temperature.

The smooth increase between the steps is due to Joule heating but cannot be recorded in real time by the temperature probe, located in the sample holder. The probe will indicate a temperature corresponding only to the initial resistance appearing at the beginning of the pulse. We use these initial values to obtain the $`R`$ versus $`T`$ function at a particular current value. Knowing this calibration function we can relate, at any time during the pulse, the resistance to the sample temperature, and in particular we can determine the temperature at which each jump occurs.

We observe that for a given current, each particular jump (small or large) always occurs at the same temperature, that we call $`T^{}`$, indicated on the right hand side of fig. 1. $`T^{}`$ depends on the current ($`T^{}(j)`$ decreases when the current increases), and by varying $`j`$, we observe the steps over a wide temperature range, from 4 K up to a temperature close to $`T_c`$ (0.7 $`T_c`$), in both NBCO and BSCCO films. Above a certain temperature, which turns out to be near $`T_c`$, the uniform resistance presumably caused by vortex motion is so high that the sample temperature ’runs away’ and the abrupt transition is replaced by a smooth upturn of the voltage as a function of time.

The data of fig. 1 indicate that in one sample for a given current density, there are several temperatures $`T^{}`$, each one corresponding to a different location along the bridge. The local character of the dissipation was established by adding a third voltage contact in the middle of the bridge: within this configuration, when measuring simultaneously the voltage on both sides of the bridge, the small jumps were detected on only one segment, never on both. As a consequence, while the large jump corresponds to a transition of the whole bridge, each small step is produced by the sudden formation of a localized highly dissipative domain (HDD). Fig. 1 and 2 show that these domains can persist over relatively long times without affecting the superconducting state in the rest of the bridge. The extreme reproducibility of the HDD formation (even when the films were warmed up to room temperature and cooled down again) suggests that the $`T^{}`$ spreading for a given current is due to slight spatial inhomogeneities along the bridge.

One of the crucial points, demonstrated in fig. 2, is that the formation of an HDD is also the origin of the global sample transition characterized by the large voltage jump. In this experiment performed on the NBCO thin film, each curve was measured with a different current (explaining the variation of the jump temperature $`T^{}`$). Above a certain current, the small step indicated by the arrow is replaced by a global transition.

We now move to the shape of $`T^{}`$ in the $`jT`$ plane. To obtain experimentally the reciprocal function $`j^{}(T)`$, we have to detect the formation of a particular HDD over the whole temperature range. The easiest way is to look at the large voltage jump. Fig. 3 illustrates the procedure to obtain one point of $`j^{}(T)`$: we fix the current, here 3.1 106 A/cm2, and then we increase the initial sample temperature, until the global transition occurs immediately at the beginning of the pulse. In this example the sample switches at 47.8 K, i.e. $`j^{}(47.8K)=3.1`$ 106 A/cm2. By applying the same procedure for small steps, we have verified that the critical current associated with the small steps has the same temperature dependence than the large jumps on some limited intervals (we could not follow a particular small jump from 4 K up to $`T_c`$). This confirms the conclusion drawn from fig. 2, namely that the small and the large steps have the same origin.

Fig. 4 shows the results for NBCO (filled circles) and BSCCO (filled triangles) films. It is remarkable that for both compounds the temperature dependence of $`j^{}`$ is the same. $`j^{}(0)`$, the extrapolated value at $`T=0`$, which represents the highest current density that can be carried by the film while remaining in the superconducting state, is also found to be similar in the two compounds. $`j^{}(0)`$ in the NBCO film was in fact slightly larger than in BSCCO (5.8 106 A/cm2 instead of 4.4 106 A/cm2), possibly due to the higher $`T_c`$ of the former.

It should be noted that upon increasing the magnetic field $`j^{}(T)`$ is globally reduced but the shape of the $`j^{}`$ versus $`T`$ curve remains the same (see fig. 4 for the 1 T case). The zero-field data are not shown here because they do not exhibit significant differences compared with the 0.1 T case.

Hence, the temperature dependence of $`j^{}`$ appears to be a general characteristics of HTS compounds. It can be well fitted by:

$$\frac{j^{}(T)}{j^{}(0)}=\left(1\left(\frac{T}{T_c}\right)^p\right)^{3/2},$$
(1)
with $`p`$ ranging from 2 to 2.5. It is worth noting that the temperature dependence described by Eq. 1 also fits other reported measurements of high current-induced phenomena. Fig. 5 displays the critical current for PSC formation found by Jelila et al.  and the critical current giving rise to a voltage instability near $`T_c`$ measured by Xiao et al. , both in YBa2Cu3O7-δ films. Like for our measurments, the temperature dependence of these data agrees with expression 1.

The picture that emerges so far is that a HDD nucleates somewhere along the bridge when the temperature reaches a critical value $`T^{}(j)`$. In some cases, it rapidly propagates to the whole sample, causing the complete switching into the normal state. The propagation speed of a highly dissipative region at high current densities has been recently measured in large YBa2Cu3O7-δ films at 77 K : it ranges from 20 to 120 m/s. In our experiment, because of the shortness of the bridge, the voltage variation during the propagation seems instantaneous. In conclusion we find that the small steps and the total transition have a common origin. To understand this origin we will focus, in the next section, on small voltage jumps.

### B I-V characteristics

To investigate the nature of the HDD, we present here the I-V characteristics of these domains at low temperature. The aim is to see the HDD signature by sweeping the current, instead of the temperature as in the previous section, across the critical value $`j^{}`$ associated with one particular domain. In this study, we chose to look at the BSCCO film, where numerous small steps could be observed.

In the following experiment the voltage response is recorded during a current ramp. Unlike the experiment using square pulses, we need that the temperature does not increase irreversibly during the application of the current, i.e that the measured voltage is history-independent (except for the possible hysteresis in the HDD formation and disappearance). This absence of heating was verified by applying upward and downward current ramps (see fig. 6 (a) and (b)). The ramp duration was reduced so that both curves matches together (see fig. 6 (c)). The only difference is the value of $`i^{}`$ for the second step. Alternatively, we observed that the I-V curves remain the same when the ramping rate is varied. For the example of fig. 6, the same characteristics are found when the ramping rate is raised from 20 A/s to 500 A/s.

Fig. 7 shows that the raw I-V characteristics (curve A) are different below and above $`i_1^{}`$, the critical current corresponding to the first step. Below $`i_1^{}`$, the I-V curve can be well fitted by a power law: $`Ui^q`$ and represents the homogeneous dissipation in the whole bridge due to vortex motion (curve B represents the extrapolation of this voltage to high currents). Above $`i_1^{}`$, the voltage comes from the uniform dissipation of vortices (curve B) and the contribution of the HDD. Therefore, to obtain the voltage of a single HDD, we have to subtract curve B from the raw data.

The resulting I-V characteristics (curve C) are non-ohmic: the I-V have piecewise linearity but the extrapolated zero-voltage current is non-zero. The differential resistance of the first domain can however be estimated: $`dU/di=0.6`$ $`\mathrm{\Omega }`$. Other domains are successively created when the current is raised (see fig. 8). A very noticeable feature is that after each jump, the slope of the linear portion increases.

It should be noted that the characteristics displayed in fig. 7, performed at 30 K, are representative of the BSCCO behaviour at low temperature, because we do not see a significant change below 30 K. For instance, as one can see in fig. 4, $`j^{}(T)`$ does not vary significantly at low temperature.

## IV discussion

Since in metallic superconductors, resistive domains result, at low temperature, from thermal instabilities, we will first examine if the voltage steps can be attributed to hot spots creation. By ”hot spot” we mean a localized domain with a temperature above $`T_c`$. A hot spot can result from a thermal instability that takes place when the heat generated in the system (which can be the film) cannot be evacuated into the heat sink. Such an instability produces necessarily domains above $`T_c`$ .

In Fig. 9 we plot the resistance of a single HDD as a function of the current, at $`T=10`$ K. Here, the resistance is defined as the voltage divided by the total current, which is the appropriate quantity for a normal domain. The HDD resistance is obtained after subtracting the background, as in fig. 7, and therefore it is the resistance we would obtain by placing local voltage probes just across the HDD.

If the HDD is a hot spot, then its length can be deduced directly from the ratio of its resistance over the total bridge resistance just above $`T_c`$, i.e 31 $`\mathrm{\Omega }`$. This yields a length of 1.3 $`\mu `$m. As a consequence, right after the nucleation the HDD would be a thin line (1.3 x 60 $`\mu `$m) cutting the bridge perpendicularly to the current direction.

The main point comes directly from the inspection of fig. 9. Since the HDD, when created, does not have an ohmic behavior, it cannot be a normal zone above $`T_c`$, unless the domain temperature or size changes as a function of the applied current. However, above $`i_1^{}`$, the I-V curve are the same with an ascending or a descending current ramp. This indicates that the application of the current does not result into an irreversible heating or expansion of the HDD. Therefore, if these changes occurs, they do in a stationnary manner.

Gurevich and Mints  have shown that, in a stationary situation, the hot spot size must decrease with an increasing current, in order to evacuate the excess of heat generated in the domain. Hence, a change in the domain size cannot explain the resistance increase of fig. 9.

Furthermore, the resistance increase is too large to be explained by a temperature increase alone. Indeed, if the value just after the nucleation ($`R_{nucl}=0.07`$ $`\mathrm{\Omega }`$) corresponds to the normal state at 90 K, then when the resistance reaches 0.26 $`\mathrm{\Omega }`$ (on the right side of fig. 9), the HDD would be over room temperature (334 K). Furthermore, in this sample, I-V measurement were performed up to 0.3 Amp, and the resistance continued to increase, corresponding to a hypothetical hot spot of 700 K. At these temperature, the HDD would irreversibly damage the film, in particular causing a loss of oxygen, and the I-V would not be reproducible.

Thus, the hot spot model cannot explain the observed voltage steps. This conclusion is also supported by the following quantitative argument.

Let us estimate the temperature increase $`\mathrm{\Delta }T`$ in a hypothetical hot spot of 1.3 $`\mu `$m, just after its creation, when we neglect the heatlink with the helium gas. Since all the power produced in excess must be evacuated toward the substrate, we have: $`R_{nucl}i_1^2=S\lambda \mathrm{\Delta }T`$, where $`S`$ is the area of the interface with the substrate and $`\lambda `$ the thermal conductance between the film and the substrate. Taking $`\lambda =10^7`$ WK-1m-2, a value generally obtained for the boundary conductance between YBa2Cu3O7-δ and MgO , the temperature increase in the HDD right after its formation would be around 2 K. This rough estimate yields that the power generated in the HDD is not sufficient to maintain the temperature above $`T_c`$.

Plotted as the voltage versus current, as in fig. 7, we see that the only deviation from an ohmic behaviour is a non-zero intersect at zero voltage. Therefore the HDD can be understood as a normal region where the current that generates dissipation is just a fraction of the total current. In other words, in addition to the normal current, a supercurrent would be flowing across the HDD.

In fact, with respect to the three features, piecewise linearity, positive zero-voltage intersect and hysteresis, the I-V curves obtained here resemble the measurements performed in tin whiskers  and microbridges . In both systems, dissipation was well explained by the standard model of PSC nucleation of Scokpol, Beasley and Tinkham (SBT) . According to this model, only part of the total current contributes to the dissipation of a PSC and this current produces a normal state resistance. Therefore, the I-V characteristics of a single PSC is given by:

$$V=2\mathrm{\Lambda }\rho (jj_0),$$
(2)
where $`\rho `$ is the normal resistivity, $`\mathrm{\Lambda }`$ the diffusion length of the quasiparticle generated in the PSC core, and $`j_0`$ the time average supercurrent. In the SBT description, the suppression of the order parameter occurs in a small region of size $`\xi `$, but the measured voltage drop is governed by the diffusion of quasiparticles over a much larger distance $`\mathrm{\Lambda }`$.

Several experiments carried out in Sn films show that the ”PSC kind” of I-V, i.e. a series of linear section separated by sharp steps is also visible in wide samples . It suggests that a 2D system exhibits the same behaviour as a narrow superconductor, but the 1D core is replaced by a thin line crossing the bridge width (i.e. a phase-slip line) or a channel of rapidly moving vortices . This extension to a 2D case is also supported by theory .

Therefore, we can evaluate $`\mathrm{\Lambda }`$, the extension of the HDD, using Eq. 2: the differential resistance of a single HDD is $`2\mathrm{\Lambda }\rho /\sigma `$, where $`\sigma `$ is the bridge cross section. This determination requires an estimate of the normal carriers resistivity $`\rho `$ at a temperature $`T<T_c`$. Taking into account the reduction of the scattering rate at low temperature, one can assume a resistivity of the form $`\rho (T)=\rho _d+\rho _s(T)`$ where $`\rho _d`$ is due to defects and becomes dominant at low temperature  and $`\rho _s(T)`$ is the linear component measured above $`T_c`$. Qualitatively such a behavior for the normal state resistivity is observed under high magnetic fields . Since we see little change below 30 K in the I-V characteristics, we estimate $`\rho _d`$ as $`\rho _nT_d/T_{onset}`$ where $`T_d`$ is the saturation temperature we estimate to be $`30`$ K and $`T_{onset}`$ is 80 K , the temperature just above the transition. This estimate yields 30 $`\mu \mathrm{\Omega }`$cm for the resistivity entering in expression (2). Therefore, the first HDD shown in fig. 7, with differential resistance of 0.6 $`\mathrm{\Omega }`$, has a characteristic length $`2\mathrm{\Lambda }=27\pm 10`$ $`\mu `$m, much larger than in the hot spot case. The error allows for the large uncertainty in $`\rho `$. Note that the first HDD could include several PSC and therefore this estimate represents an upper limit of the PSC size.

Coming back to the temperature dependence of $`j^{}`$, it is interesting to note that such a temperature dependence is what is most likely expected for the depairing current $`j_d(T)`$ in high $`T_c`$ superconductors. Indeed, using Ginzburg-Landau relations for $`j_d`$, we have:

$$j_d(T)\frac{1}{\lambda ^2(T)}\frac{1}{\xi (T)}\left(\frac{1}{\lambda ^2}\right)^{3/2},$$
(3)
where the $`1/\lambda ^2(T)`$ factor is the temperature dependence coming from the superfluid density and $`1/\xi (T)`$ is the contribution of the critical velocity. The last proprotionality is obtained by assuming that the Ginzburg-Landau parameter $`\kappa =\lambda /\xi `$ is only weakly temperature dependent. The temperature dependence of $`j_d`$ will be the one displayed on the right hand side of eq. 1, because the temperature dependence of the superfluid density measured from zero to $`T_c`$ is well described by $`1(T/T_c)^p`$ with $`p`$ ranging from 1.7 to 2.9 .

However, the amplitude of $`j^{}`$ is well below the depairing current (about two orders of magnitude). In a perfect superconducting filament, in contrast, a PSC would develop when the supercarriers reaches the critical velocity. The fact that a HDD appears at a current representing a small fraction of $`j_d`$ implies that the energy cost to create a phase-slip line or a channel of fast moving vortices is well below the one needed to destroy completely superconductivity in the film. A theoretical prediction of $`j^{}`$ would require an understanding of the exact mechanism causing the development of the PSC. This is presently beyond the scope of this article. Nevertheless, we can deduce from the similarity of $`j^{}`$ in the different compounds that this mechanism is not dependent on the nature of the defects, and therefore not related to vortex pinning.

## V conclusion

In every high $`T_c`$ thin film of BSCCO and NBCO we measured, the high current part of the I-V characteristics exhibits some sharp voltage variations. These are robust features and are observed over a wide temperature range, from 4 K up to the vicinity of $`T_c`$. This behavior of superconducting films traversed by a high density of current seems very systematic, since it has been reported by other groups in YBa2Cu3O7-δ  or Nd2-xCexCuOy , with various interpretations. We believe that a very general mechanism occurs in the superconducting state at high current density independently of the defects, such as pinning centers, etc. By means of time resolved resistance measurements, we have shown that this mechanism accounts for both small voltage steps and a more dramatic event like the current-induced breakdown of superconductivity. The small voltage steps are associated with the nucleation of local dissipative domains. We have brought new arguments showing that these domains do not result from thermal instabilities, like hot spots (the situation is quite different in bulk superconductors, where Joule heating may dominate). Another mechanism, directly related to superconductivity must therefore be at the origin of the formation of these domains. Even though the samples cannot be considered as one-dimensional, the transport properties are in good agreement with a dissipative state dominated by phase-slip centers, in which the superconductor is locally driven far from its equilibrium point. If this interpretation is correct, transport measurements at high current density allow to access intrinsic properties of the superconductor. The fundamental character of this dissipation process is confirmed by the temperature dependence of the critical current which is similar in different high $`T_c`$ compounds and turns out to be close to the one expected for the depairing current.

###### Acknowledgements.

We thank P. Martinoli for stimulating discussion and careful reading of the manuscript. This work was supported by the Swiss National Science Foundation.


text:
# 1 Introduction

## 1 Introduction

The random percolation problem in 2-d  can be extended to a number of interesting variations, e.g. the ‘bootstrapping model’ , where certain sites are culled leading to a change in percolation behaviour or the ‘diagenesis model’ , with growing as well as culling. In the present model we show that an introduction of simultaneous growing and culling processes, followed by a randomisation after each time step gives very interesting behaviour.

The model may be applicable to a real situation such as condensation/evaporation or adsorption/desorption of a layer of ‘molecules’ at different surface temperatures. We develop two versions of the model in subsequent sections . In model I, we assume ‘growth’, that is condensation at sites where there are a large number of occupied nearest neighbours and ‘evaporation’ at sites with most neighbouring sites vacant. This algorithm mimics an attractive interaction between the ‘molecules’. We find that a certain coverage ($`p_{inv}`$) of the surface is invariant, depending on the exact algorithm employed. An initial concentration of occupied sites $`p>p_{inv}`$ leads to complete coverage, whereas a lower starting point $`p<p_{inv}`$ leads to zero coverage, i.e. complete evaporation. In section 2 we discuss the algorithm for growing and culling and in section 3 a hypothetical case with a symmetric rule for activated growing and culling is discussed. In sections 4 and 5, more realistic models II and III for the condensation/evaporation process as function of vapor pressure and surface temperature are presented. In the last section conclusions and future plans are discussed.

## 2 Model I

We start with a square lattice in two-dimensions and sites are randomly occupied by particles with a probability $`p`$. Then we grow new particles at vacancies, or cull particles at occupied sites, according to one of the rules below.

Rule 1
(a) occupied sites having all four nearest neighbour (nn) sites vacant and occupied sites with only one nn site occupied are culled. (b) A new particle grows, at a vacant site with four nn positions occupied, and also at a vacant site where only three nn sites are occupied.

Rule 2
(a) Only particles with four vacant nn are culled. (b) Only sites with four occupied nn ‘grow’ a new particle.

Rule 3
(a)Only particles with four vacant nn are culled. (b) A new particle grows, at a vacant site with four nn positions occupied, and also at a vacant site where only three nn sites are occupied.

Rule 4
(a) Occupied sites having all four n.n. sites vacant and occupied sites with only one nn site occupied are culled. (b) Only sites with four nn occupied ‘grow’ a new particle.

In all the above rules particles or vacant sites with 2 nn occupied are left undisturbed. Figure 1 illustrates the rules pictorially.

After completing the grow-cull operations, the concentration of particles changes from an initial value $`p_i`$ to a final value $`p_f`$. We now randomise the positions of the remaining sites over the whole lattice, and repeat the grow-cull operations with $`p_f`$ becoming the new $`p_i`$. The results of the above procedure for the four different rules outlined above can be determined analytically.

For each of the four rules we find an initial coverage $`p_{inv}`$ which remains invariant after repeating the steps of growing-culling followed by randomisation. $`p_{inv}`$ corresponds to the coverage for which the probability of growth equals the probability of culling for the particular rule employed. For some other starting concentration say $`p_0>p_{inv}`$, we would get a new $`p=p_1`$ which is $`>p_0`$, since the growth probability exceeds culling probability. The randomisation that follows makes $`p_1`$ the new initial coverage, which after growing-culling gives a still larger $`p_2`$, and so on. So the coverage approaches 1.0 or 100 %, as shown in figure 2.

For $`p_0<p_{inv}`$, on the other hand $`p_1<p_0`$ and the system evolves towards $`p=0`$, i.e. zero coverage, $`p_{inv}`$ is an unstable fixed point for the system, whereas $`p=0`$ and $`p=1`$ are stable fixed points.

The processes described above can be studied by computer simulation as well as analytically. Using a parallel algorithm for growth and culling gives a result exactly in agreement with the calculated result, while a sequential algorithm gives different results as expected.

The stable fixed point can be determined as follows. For a coverage $`p`$, the probability of growing and culling $`P_{gr}`$ and $`P_{cl}`$ for the four different rules can be written as follows. For Rule 1

$$P_{cl}=p(1p)^4+4p^2(1p)^3$$
(1)
$$P_{gr}=p^4(1p)+4p^3(1p)^2$$
(2)
In eq.1 the 1st term on the right is the probability for an occupied site having 4 vacant nn, in a random distribution. since $`p`$ is the probability of a site being occupied, and $`(1p)`$ is the probability of being vacant. The other terms can be written down similarly with proper weight factors.

For the other rules, we have analogous relations. For rule 2

$$P_{cl}=p(1p)^4$$
(3)
$$P_{gr}=p^4(1p)$$
(4)
For Rule 3

$$P_{cl}=p(1p)^4$$
(5)
$$P_{gr}=p^4(1p)+4p^3(1p)^2$$
(6)
and Rule 4

$$P_{cl}=p(1p)^4+4p^2(1p)^3$$
(7)
$$P_{gr}=p^4(1p)$$
(8)
$`p_{inv}`$ is easily obtained for any of the rules by setting $`P_{cl}=P_{gr}`$ and solving for $`p`$.

For Rules 1 and 2, which are symmetric, we find $`p_{inv}=0.5`$, while for Rule 3 and Rule 4 $`p_{inv}`$ has the complementary values 0.32 and 0.68 respectively. The same values are obtained for a steady coverage by carrying out a computer simulation of the processes described by the rules.

## 3 Activated growth and culling

From physical considerations it is expected that - the probability of condensation on a surface decreases with increasing surface temperature, and probability of evaporation increases . If we assume that rule 3 holds at high temperatures and rule 4 at low temperature, we may introduce a hypothetical temperature dependent growth and culling rule as follows.

$$P_{cl}=p(1p)^4+4exp(E_0/kT)p^2(1p)^3$$
(9)
$$P_{gr}=p^4(1p)+4(1exp(E_0/kT))p^3(1p)^2$$
(10)
The temperature dependence enters through the exponential factor with $`E_0`$ as an energy characteristic of the system.

The new rule can be seen to reduce to Rule 3 for $`T\mathrm{}`$, and to Rule 4 for $`T0`$. The fixed point $`p_{inv}`$ goes accordingly from 0.32 to 0.68 as temperature is lowered from a very high value compared to a characteristic temperature $`T_0=E_0/k`$. $`p_{inv}`$ can be determined by equating the growth and culling probabilities as before at different temperatures.

The exponential temperature dependence introduced here for evaporation is realistic, but the rule for condensation is purely hypothetical, introduced to make the equation symmetric. In the next section we discuss a more realistic approach. The final coverage in this case is either $`0`$ or $`1`$ depending on the temperature and the initial coverage.

## 4  Condensation and Evaporation

In this section we study the condensation/evaporation problem using the approach developed in the preceding section. We find that computer simulation of the problem is instructive, and may lead to a reinterpretation of some existing ideas.

The process of condensation/evaporation or adsorption/desorption is described traditionally by two different sets of models – one for a mobile adsorption layer and one for an immobile layer. Probabilities for sticking and evaporation on a two-dimensional monolayer are specified according to the physics behind the model.These are functions of the temperature, superincumbent pressure and the existing coverage. At equilibrium, the sticking and evaporation probabilities are set equal and the resulting equation is solved to get the equilibrium coverage at that temperature and pressure.

The simplest mobile layer model is a two-dimensional ideal gas, and the improved versions include interaction between particles, similar to a two-dimensional Van der Waal gas. The so-called ‘immobile layer’ models introduce a sticking probability, depending on how long a molecule in the vapor above the surface is in contact with a surface site. The simplest ‘immobile model’ is the Langmuir equation derived as follows.

For vapor condensation the particle flux, i.e. the number of particles deposited per unit time per unit surface area is equal to:

$$c=\frac{P\lambda }{h}$$
(11)
where $`P`$ is the pressure,

$$\lambda =\sqrt{h^2/(2\pi mkT)}$$

is de Broglie length, $`h`$ is Planck constant.

Evaporation probability from saturated surface (at $`p=1`$) may be approximated by  :

$$d=\frac{kT}{h\lambda ^2}\mathrm{exp}(E_e/kT)$$
(12)
Condensation probability is set equal to evaporation probability, giving equilibrium:

$$\frac{P\lambda }{h}(1p)=\frac{kT}{h\lambda ^2}\mathrm{exp}(E_e/kT)p$$
(13)
or

$$\frac{P\lambda ^3}{kT}(1p)=\mathrm{exp}(E_e/kT)p$$
(14)
and we have the simple Langmuir equation

$$Pb=\frac{p}{1p}$$
(15)
where

$$b=\frac{\lambda ^3}{kT}\mathrm{exp}(E_e/kT),$$

$`E_e`$ is activation energy for evaporation.

In this approximation it is assumed that evaporation energy is the same for all configurations.

A modification of this model is the Fowler-Guggenheim model , one form of this is given below. Here the evaporation probability for a particle depends on the number of occupied neighbors. Each of the neighbors exerts an attractive force on the particle, which must be overcome for evaporation. We introduce the following parameters for convenience, - $`f`$ is a pressure parameter given by

$$f=P\mathrm{\Lambda }^3/gkT$$

$$\mathrm{\Lambda }=\sqrt{T^{}}\lambda $$

and

$$g=exp(E_e/kT)=exp(T_0/T)=exp(1/T^{})$$

where $`T_0=E_e/k`$ is a characteristic evaporation temperature, $`T^{}=T/T_0`$ is a reduced temperature.

The condition for condensation rate to equal evaporation rate is now

$$f(1p)/(T^{})^{5/2}=p^5g^4+4p^4g^3(1p)+6p^3g^2(1p)^2+4p^2g(1p)^3+p(1p)^4$$
(16)
The successive terms on the RHS of equation (16) represent probabilities for a particle to have 4, 3, 2, 1 and 0 nearest neighbours in a random distribution On rearrangement this reduces to

$$f(1p)/(T^{})^{5/2}=gp(1+pgp)^4$$
(17)
## 5 Condensation/evaporation from a different viewpoint: Models II and III

In the above models the question of mobility of the adsorption layer is not introduced explicitly, though the FG model is classified as an immobile model. Moreover the probabilities for evaporation considered in eq.(16). are valid only as long as the distribution is random. For a strictly immobile layer the distribution ceases to be random once the site-dependent evaporation starts. We formulate the problem so that the lateral mobility if present is introduced explicitly, and we can look at both mobile and immobile situations within the same framework. As in the previous section, we start with a two-dimensional lattice with a certain fraction occupied randomly by the particles. The adsorbate particles are also present as vapor above the surface and the pressure and temperature have a key role to play.

We visualise the condensation/evaporation (or adsorption/desorption) as a two-step process, there is one characteristic time for the sticking and evaporation and another for lateral diffusion of the molecules over the surface. For the ‘mobile interface’ situation, the timescale for lateral motion is very small, and each condensation-evaporation step is followed by a complete randomisation. The other extreme is the ‘immobile interface’, here we drop the randomisation process altogether. It is also possible to consider intermediate situations where the two characteristic times are comparable.

According to this picture, setting the sticking probability equal to the evaporation probability at equilibrium, is valid for the mobile situation only, not the immobile case. This is because, if the surface molecules cannot move laterally, after one round of growth and evaporation, the expressions for condensation/evaporation probability are no longer valid, because the distribution is no longer random. So for the mobile case, we have a result similar to the Fowler-Guggenheim formalism, but with a different interpretation. Our solution for the final coverage depends on the initial coverage, besides temperature and pressure. The coverage isotherm of the system may follow a different path during increasing and decreasing pressure showing hysteresis. This mobile case, can be worked out by analytical calculation as well as computer simulation.

The ‘immobile surface’ case, cannot however, be calculated analytically, as we do not know the condensation/evaporation probabilities after one round of growth and culling, since the system has lost its random distribution. But we can still simulate the system on a computer. We have two cases - firstly the Fowler-Guggenheim equation in the form of eq.(16), we call this model II. Secondly we take a situation where only particles with one occupied neighbor or none at all (isolated particles) are allowed to evaporate (as in rule 1a, section 2), we call this model III. The corresponding equation is

$$f(1p)/(T^{})^{5/2}=4p^2(1p)^3+p(1p)^4$$
(18)
Temperature and pressure dependence of sticking probabilities are assigned as in FG model eq. (16).

### 5.1 Mobile interface layer – analytical study

Here, there is a complete randomisation after each growing-culling sequence, so it is meaningful to equate the probability of growth to the probability of culling and solve for the invariant coverage. We can look at the process as an iteration of the following two steps

$$p_f=p_i+P_{gr}P_{cl}$$
(19)
and

$$p_i=p_f$$
(20)
here $`p_i`$ is the initial coverage and $`p_f`$ the final. Also $`P_{gr}=P\lambda (1p)/h`$ and $`P_{cl}=wT^^2exp(1/T^{})`$. Here $`w`$ is a parameter, with a suitably chosen arbitrary value. This iteration done numerically gives the same result as solving analytically

$$P_{gr}=P_{cl}$$
(21)
and also agrees with an explicit computer simulation of the process. At certain values of temperature and pressure, there are three solutions for the coverage. The middle one is an unstable fixed point(UFP), and the other two stable fixed points. So a starting coverage above the upper fixed point leads to the coverage stabilising at the upper fixed point, whereas if we start with a coverage below UFP we end up at the lower fixed point. The computer simulation shows exactly the same behaviour. Figures \[3a\] and \[4a\] show the results for eqs (17) and (18). The difference with the standard thermodynamic treatment is as follows. If we start from a very low coverage say at a temperature $`0.5`$ and pressure $`10^5`$ and gradually increase the pressure, according to earlier FG model the coverage increases as shown in fig3a until a phase transition to the upper fixed point takes place according to the equal area Maxwell’s rule . In the present model however, stable coverage depends on initial $`p`$ as well as $`T^{}`$ and $`f`$. If we start with a very low $`p`$ at low pressure, the system follows the same path as FG initially but undergoes the phase transition later as shown in fig.3a (pt. B) at a pressure where the lower fixed point meets the UFP. It then follows path BC. While decreasing the pressure from a higher value the system follows the different path shown in the figure, with a phase transition where the upper fixed point coincides with UFP (pt.D). So here we have a marked hysteresis. Hysteresis in adsorption/desorption is usually attributed to presence of pores , but we see here that it may also have a different origin. In model III (see fig.4a) the reverse path starting from high pressure and complete coverage shows no decrease in $`p`$ but continues to low pressures with $`p=1`$. This is because in this model there is no evaporation unless some sites are vacant. Model II looks more realistic under normal conditions. It must be noted that this is quite different from the usual form of Fowler-Guggenheim, where one would not expect hysteresis, and the present version considers a mobile layer.

### 5.2 Computer simulation of the mobile case

In the mobile case, the adsorbed molecules can move laterally on the surface. In model II, the physical situation simulated is described by eq. (16). A two-dimensional square lattice of unit spacing and size $`300\times 300`$, is occupied randomly with an initial coverage $`p_{initial}`$. Every occupied site is assigned the value $`1`$, empty sites are assigned the value $`0`$. The occupied sites are then culled paralally with a probability determined by the number of their occupied nearest neighbours. The site having n occupied neighbours, has the culling probability $`p^{n+1}g^n`$, where $`p`$ is the occupation probability and

$$g=exp(1/T^{})$$

The vacant sites are filled with a probability $`f/(T^{})^{5/2}`$ where $`f`$ is the pressure parameter defined earlier. After one round of growth and culling is complete, the concentration of the occupied sites $`p_{final}`$ is calculated.

In the next time-step, the $`p_{final}`$ of the previous time-step becomes the new $`p_{initial}`$. The square lattice is then randomly occupied afresh with this $`p_{initial}`$. A complete time-step begins with the random occupation of all sites with a $`p_{initial}`$ and ends with the assignment of the $`p_{final}`$ to the $`p_{initial}`$ of the next time-step. This iterative process stops when $`p_{final}=0`$ or $`1`$ or $`p_{final}`$ saturates with increasing time to a definite value. In the simulation, we checked upto 50,000 time-steps. The ‘mobility’ of the molecules is simulated by the randomization of the concentration $`p_{initial}`$ in the beginning of every time-step.

In model III described by eq.(18), the same iterative process described in model II is carried out, except for the condition of culling. In this case an occupied site is culled with the probability $`g`$ if the sum of its nearest neighbours is less than two.

All the sites in both the models are updated parallelly and periodic boundary conditions are applied both along the x- and y- directions. Both the models are studied over an effective temperature range of $`0.5`$ to $`3`$ and with $`f`$ varying from $`10^4`$ to $`10^2`$, for the entire range of $`p_{initial}`$ from $`0`$ to $`1`$. The simulation results agree to within $`10^3`$ of the numerical results. These are presented in the figures ($`3a,3b,4a,4b`$ and $`5`$).

A significant difference between models II and III is evident from figs.3a and 4a, showing isotherms for different temperatures. In model II there is a critical temperature above which there is no phase transition, but in model III there is always a phase transition. Fig.5 shows how the coverage evolves with time for $`T^{}=0.5`$ and $`f=0.02`$ for $`p_{initial}`$ varying from 0 to 1.

### 5.3 Immobile interface layer

Let us now make the adsorbed molecules immobile. In this case obviously, analytical calculation is not feasible. If we start with an initial random configuration with a certain coverage, as soon growth and culling at preferential locations starts we can no longer calculate probabilities for further evolution exactly. So here we resort to computer simulation. The simulation of the ‘immobile’ case of the models II and III follow the same iterative process, except for the randomization of the sites with the given probability at the beginning of every time-step. Here a new round of growth and culling commences on the geometrical configuration reached at the end of the last time-step. Since in the preceding section we saw that computer simulation results agree with calculations, we are confident that the simulation gives reliable results, with the system size and algorithm used.

We find that the results are quite different from the $`mobile`$ case. Whereas in the mobile case, we found definite fixed points where the coverage converged regardless of the exact starting $`p`$ (see fig.5), here for certain ranges of the temperature and pressure parameters, we get a different stable coverage for even closely spaced initial $`p`$ values. Figs.(6a & 6b) show the time evolution of the coverage for the immobile system for typical situations under quite different conditions of pressure. In fig.6a. for $`p_{initial}`$ above 0.3 the coverage always goes to 1 at rates depending on the initial coverage, but for lower starting points it stabilises to different values for each $`p_{initial}`$. In a similar study of the $`mobile`$ case one always ends up at either the upper or the lower fixed point (fig.5). But fig.6a shows that for the immobile case, initially culling dominates, for large $`p_{initial}`$ there is a rapid drop in coverage within a few hundred time-steps. After that sticking catches up, and for $`p>0.3`$ (approximately) there is a steady increase, and $`p1`$ linearly with a slope decreasing as $`p_{initial}`$ decreases. For $`p_{initial}`$ lower than $`0.3`$ $`p`$ apparently saturates, showing no variation upto 50,000 time-steps. But when the same situation is studied at a much lower $`f`$ ($`10^{10}`$), (fig.6b), $`p_{final}`$ always saturates to different values depending on $`p_{initial}`$, over the time-scales checked. For example, in fig.6a., for $`p_{initial}=0.9`$, $`p_{final}=1`$ whereas the same $`p_{initial}`$ in fig.6b. saturates to $`p_{final}`$ of 0.8979. This suggests that even for the same $`p_{initial}`$ there must be a definite combination of pressure and temperature where there is a change in the behaviour – $`p_{final}`$ reaching the stable fixed point $`1`$ or a constant value $`<1`$. Further investigations probing the exact phase diagram is in order and will be done in the near future.

All the above figures are the results on model III. Similar studies on model II show no significant qualitative difference with model III except for the exact numerical values.

## 6 Conclusions

We have investigated the behaviour of a two-dimensional surface layer explicitly allowing or forbidding lateral motion of the particles. We find that the cases for a mobile and an immobile layer show quite different characteristics, moreover none agrees fully with the Langmuir or FG model. The fixed points for our model II with ‘mobile’ particles, are the same as the FG solutions, which are supposed to be valid for an immobile layer. Details of the evolution of the surface coverage from low to high pressure are however different. For the immobile layer the results are again different. These results are similar to what was reported in a previous work attempting to simulate diagenetic processes of restructuring in sandstones . However the two-dimensional character is more appropriate in case of the condensation/evaporation problem. It will be interesting if experimental support for the behaviour described here can be demonstrated. Hysteresis has been observed for adsorption on porous surfaces , but here we see that competition between site-specific sticking and evaporation may also cause hysteresis.

We can simulate intermediate behaviour between mobile and immobile layers, by allowing the particles to execute a random walk for some time on the surface before the next grow-cull operation. Making the time for the walk extremely large will correspond to the complete randomisation done here. We have plans to study this in future.

Recent work on adsorption/desorption problems \[8-13\] show that considerable theoretical and experimental studies are being done on such problems. They are of interest in modern devices which use surface properties extensively and also environment relate problems such as explaining the ozone hole in polar regions . Complicated mathematical and computer simulation methods such as density functional theory and molecular dynamics simulations are being used. In this work we show that very simple monte carlo simulations also reveal some interesting features and may turn out to be quite useful in shedding some light on such problems. In conclusion the growth-culling-randomisation extension of the standard 2-d percolation problem promises to yield more new and interesting results.

## 7 Acknowledgement

DST, Govt. of India and Ministry of Ukraine for Education and Science, are gratefully acknowledged for grant of an Indo-Ukrainian collaboration project. Authors thank S S Manna and J K Bhattacharjee for suggestions and discussion on this subject. Authors are grateful to S N Bose National Centre for Basic Sciences for extending their computer facilities.

figure Captions :

Figure 1. The growth and culling processes are illustrated, the rules are implemented as follows : rule 1 - (a), (b) (d) and (e), rule 2 - (a) and (e), rule 3 - (a), (d) and (e), rule - 4 (a), (b) and (e).

Figure 2. This is a schematic diagram showing how the system evolves towards the fixed points, here the unstable fixed point is p = 0.5, valid for rules 1 and 2.

Figure 3a. Coverage vs. pressure for the mobile case (model II) with $`T^{}=0.5`$ If we start with a very low coverage ($`p`$) and gradually increase pressure, $`p`$ increases along the curve as shown upto B. After this it undergoes a phase transition to $`p=1`$ along BC. On decreasing pressure, it follows the path CDA showing a well marked hysterisis.

Figure 3b. Plot of final coverage vs. pressure for different $`T^{}`$.

Figure 4a. Plot of coverage vs. pressure for $`T^{}=0.5`$ for model III. Here the final coverage depends on initial coverage. For example, if we start at a point with coverage $`<B`$ the final coverage reaches point B. For any starting coverage between points B and C, the final coverage is B. For starting coverage above C, the final coverage is point D ($`p=1`$). Starting with a low coverage, if pressure is gradually increased, the coverage proceeds along the curve OBE where it undergoes phase transition to F.

Figure 4b. Plot of final coverage vs. pressure for different $`T^{}`$ (model III).

Figure 5. Plot of coverage vs. time for the mobile case. For initial coverages $`>0.68`$ (dotted lines) the final coverage is $`p=1`$. For coverages $`<0.68`$ (solid lines) the final $`p`$ is the lower fixed point.

Figure 6a. Plot of coverage vs. time for the immobile case in model III for$`T^{}=4`$ and $`f=10^4`$. For initial coverage $`<0.35`$, the final coverage saturates to different values. But for larger initial coverage, the system approaches $`p=1`$. The numbers on the curves indicate the initial coverage. The y- axis is shifted from $`t=0`$ to show clearly the initial drop of the coverage at early times.

Figure 6b. Plot of coverage vs. time for the immobile case in model III for $`T^{}=0.3`$ and $`f=10^{10}`$. Here the final coverage saturates to different values for all the initial coverages studied over 50,000 time-steps.


text:
# Multifractal Analysis of the Dynamical Heterogeneity in a Two-Dimensional Supercooled Liquid [1]

## Abstract

The dynamical heterogeneity in supercooled liquids measured by a molecular dynamics simulation has been quantified on the basis of the multifractal formalism. The singularity spectrum becomes broader as the glass transition is approached. This behavior is similar to that observed in the case of the Anderson-localization transition. Keywords: supercooled liquid, glass, dynamical heterogeneity, multifractal

The dynamical heterogeneity in supercooled liquids has been studied intensively in recent years. It can be quantified by introducing the concept of bond-breakage and the distribution of the broken bonds near the glass-transition is found to be similar to that of the spin state of the Ising model near the critical point. The bond-breakage corresponds to the collective jump motion of particles among cages.

On the other hand, we have analyzed the dynamical heterogeneity of the critical spin state of the Ising model on the basis of the multifractal formalism. In this Short Note we try to combine above two approaches and quantify the dynamical heterogeneity observed in a molecular dynamics simulation of a supercooled liquid in terms of multifractality. The so-called cooperatively rearranging regions (CRR) can be identified with weakly bonded regions with various sizes which form multifractal structure.

We use the same model as in Ref. 1, the binary mixture of soft core particles. Two species of particles, 1 and 2, with the size $`\sigma _1`$ and $`\sigma _2`$, and the mass $`m_1`$ and $`m_2`$ are considered. The interaction potential is given by

$$v_{\alpha \beta }(r)=ϵ(\frac{\sigma _{\alpha \beta }}{r})^{12},$$
(1)
where $`r`$ is the distance between two particles and $`\sigma _{\alpha \beta }=(\sigma _\alpha +\sigma _\beta )/2`$ with $`\alpha ,\beta =1,2`$.

The molecular dynamics simulation is done using the second order symplectic integrator and the microcanonical ensemble. The temperatures in the following are determined by the kinetic energy.

The bond is introduced as in Ref. 1 when

$$r_{ij}(t_0)1.1\sigma _{\alpha \beta },$$
(2)
where $`r_{ij}(t_0)`$ is the distance between $`i`$-th particle of species $`\alpha `$ and $`j`$-th particle of species $`\beta `$ at time $`t_0`$. The position of the bond is represented by the midpoint of the two particles at this time. The bonds are unambiguously defined owing to the sharpness of the first peak of their pair correlation function. The bond is regarded as surviving at time $`t_0+t`$ when

$$r_{ij}(t_0+t)1.6\sigma _{\alpha \beta }.$$
(3)
The number of surviving bonds $`N_{\mathrm{suv}}(t_0+t)`$ is fitted by

$$N_{\mathrm{suv}}(t_0+t)=N_{\mathrm{suv}}(t_0)\mathrm{exp}(t/\tau _\mathrm{b}),$$
(4)
introducing the life-time of bonds $`\tau _\mathrm{b}`$. In the following the broken bonds are determined for $`t=0.05\tau _\mathrm{b}`$.

In our simulation the followings are fixed. Two species of particles have the same particle number, $`N_1=N_2=80000`$. The size and mass ratios are chosen as $`\sigma _2/\sigma _1=1.4`$ and $`m_2/m_1=2`$ respectively. The time step of the integration $`\mathrm{\Delta }t`$ is chosen as $`\mathrm{\Delta }t=0.005\tau _1`$ with $`\tau _1=(m_1\sigma _1^2/ϵ)^{1/2}`$. The potential is truncated at $`r=4.5\sigma _1`$. The particle density is chosen as $`n=0.8/\sigma _1^2`$. The simulations were performed on our personal computer with Pentium-4 CPU.

In Fig. 1 the snapshot of the broken bonds is shown where the position of the bond is defined at time $`t_0`$ as mentioned in the above. In Figs. 1(a) and 1(b) the temperatures are chosen as $`T=2.54`$ and $`T=0.859`$ respectively in the unit of $`ϵ/k_\mathrm{B}`$. These two states are realized starting from an initial state with the temperature 9.69 and cooled in a stepwise manner as in Ref. 1. In the stepwise cooling process the temperatures in the first, second, third and fourth stages are 9.69, 2.54, 1.43 and 0.859 respectively. For the former state the temperature was changed from 9.69 to 2.54 after spending 5000$`\mathrm{\Delta }t`$ at the first stage and the broken bonds are observed at the second stage with $`t_0=5000\mathrm{\Delta }t`$. For the latter state the temperature was changed successively from 9.69, 2.54, 1.43 to 0.859 after spending 5000$`\mathrm{\Delta }t`$ at the first, second and third stages each and the broken bonds are observed at the fourth stage with $`t_0=60000\mathrm{\Delta }t`$. The times $`t_0`$ and $`t`$ are measured from the instance of the temperature change to the stage of the observation. The life-times of bonds were measured to be $`\tau _\mathrm{b}=4.09\times 10^3\mathrm{\Delta }t`$ for $`T=2.54`$ and $`\tau _\mathrm{b}=2.43\times 10^5\mathrm{\Delta }t`$ for $`T=0.859`$.

Qualitatively the snapshot shows a multifractal behavior and the degree of intermittency increases as the temperature is decreased. In the following we quantify this observation on the basis of the multifractal formalism.

The multifractality is quantified by the singularity spectrum, $`f(\alpha )`$, and the spectrum is obtained by the following formulae,

$$\alpha (q)=\underset{k=1}{\overset{N_L}{}}\mu _k(q,l)\mathrm{ln}\mu _k(1,l)/\mathrm{ln}(l/L),$$
(5)
and

$$f(\alpha (q))=\underset{k=1}{\overset{N_L}{}}\mu _k(q,l)\mathrm{ln}\mu _k(q,l)/\mathrm{ln}(l/L),$$
(6)
with

$$\mu _k(q,l)=\{\mu _k(l)\}^q/\underset{k^{}=1}{\overset{N_L}{}}\{\mu _k^{}(l)\}^q.$$
(7)
Here the measure, $`\mu _k(l)`$, is defined as the probability to find a broken bond in the $`k`$-th square box of the linear dimension $`l`$ in a snapshot of broken bonds. The total number of the boxes $`N_L=(L/l)\times (L/l)`$ where $`L`$ is the linear dimension of the simulation region with periodic boundary condition.

In Fig. 2 the singurarity spectra for $`T=2.54`$ and $`T=0.859`$ are shown. The dots represent the data set $`(\alpha (q),f(\alpha (q)))`$ determined by the above formulae and should be tangential to $`f(\alpha )=2`$ and $`f(\alpha )=\alpha `$. The maximum $`\alpha _{\mathrm{max}}`$ and minimum $`\alpha _{\mathrm{min}}`$ are determined by the intercept of the spectrum extrapolated to the $`\alpha `$-axis and correspond to the singularity exponent for the box with the lowest density of broken bonds and the one with the highest density respectively.

From the data in Fig. 2 we can quantify the degree of intermittency. The singularity spectrum becomes broader as the temperature is decreased approaching the glass transition. At the same time $`\alpha _{\mathrm{max}}`$ becomes larger and $`\alpha _{\mathrm{min}}`$ becomes smaller. Namely $`\alpha _{\mathrm{max}}=2.74`$ and $`\alpha _{\mathrm{min}}=1.66`$ for $`T=2.54`$ and $`\alpha _{\mathrm{max}}=2.97`$ and $`\alpha _{\mathrm{min}}=1.57`$ for $`T=0.859`$. This finding is similar to that reported in the case of the Anderson localization where the singularity spectrum becomes broader as the strength of the impurity potential is increased and the localization transition is approached. In both cases the degree of intermittency corresponds to that of localization.

In conclusion we have quantified the dynamical heterogeneity in a model supercooled liquid measured by a molecular dynamics simulation on the basis of the multifractal formalism. In supercooled liquids various sizes of CRR are formed and constitute multifractal structures. The degree of intermittency increases approaching the glass transition. Since our simulation parameters are still far from those for the glass transition, we need further study to elucidate the nature of the glass transition.

The authors are grateful to Jun Matsui for enlightening discussions on the molecular dynamics simulation. This work was supported in part by a Grand-in-Aid for Scientific Research from the Ministry of Education, Culture, Sports, Science and Technology of Japan.


text:
# Minimal Brownian Ratchet: An Exactly Solvable Model

(Last modified: June 19, 2003. Printed: )

## Abstract

We develop an analytically solvable three-state discrete-time minimal Brownian ratchet (MBR), where the transition probabilities between states are asymmetric. By solving the master equations we obtain the steady-state probabilities. Generally the steady-state solution does not display detailed balance, giving rise to an induced directional motion in the MBR. For a reduced two-dimensional parameter space we find the null-curve on which the net current vanishes and detailed balance holds. A system on this curve is said to be balanced. On the null-curve, an additional source of external random noise is introduced to show that a directional motion can be induced under the zero overall driving force.

The Brownian ratchet and pawl system was first correctly explained by Smoluchowski smoluchowski1912  and later revisited by Feynman feynman1963  – this has inspired much activity in the area of Brownian ratchets, despite flaws in Feynman’s analysis of the thermal efficiency of the ratchet engine parrondo1996  and detailed balance abbott2000 .

Interest has revived because molecular motors astumian1994  have been described in terms of Brownian ratchet doering1995 ; magnasco1993  models. Another area of interest has been in Parrondo’s paradox harmer1999a  where losing strategies cooperate to win. This can be illustrated in terms of games that lose when played individually, but win when alternated – this has been shown to be a discrete-time Brownian ratchet harmer1999b , otherwise known as a ‘Parrondian game’. Parrondo’s games have significantly sparked recent interest in the areas of lattice gas automata meyer02 , spin models moraal00 , random walks and diffusions key02 ; pyke02 ; cleuren02 , biogenesis davies01 , molecular transport kinderlehrer02 ; heath02 , noise induced patterns buceta02 , stochastic control allison01a ; kocarev02 , stochastic resonance allison01b  and quantum game theory flitney02 ; lee02 . Recently, Reimann reimann02  has performed an extensive review of the ratchet field.

Jarzynski et al. jarzynski1999  developed an exactly solvable Brownian ratchet that can be operated as a heating system or refrigerator, depending on the parameters between two heat reservoirs of different temperatures. However this is treated as a six state system and solution is via matrix inversion of coupled linear equations. The derivation is somewhat complex, so the physical picture and key ingredients of the observed properties are obscured.

Westerhoff et al. westerhoff1986  have analyzed enzyme transport using a four-state model. In this paper, for the first time, we develop a three-state discrete-time Brownian ratchet model that can be solved analytically. We call it the Minimal Brownian Ratchet (MBR) footnote\_01 . By setting up and solving the steady-state solution of the corresponding master equations we obtain the probability current and the null-surface, in the parameter space, of the noisy and noise-free MBR. The obtained solution does not show any critical behavior and can be suitably explained in terms of non-singular behaviors.

The minimal ingredients of a Brownian ratchet are an asymmetric potential and random noise. In Fig. 1 we show the state diagram of the MBR. The MBR has three states, $`\{S_0,S_1,S_2\}`$, where the transition probabilities between states are asymmetric. The transition probability that a random walker in state $`S_k`$ steps in the positive direction is $`p_k`$. The probability of a shift in the negative direction is $`\stackrel{~}{p_k}`$. This is true for $`k\{0,1,2\}`$. We define the positive direction as counterclockwise. The condition of normalization, $`p_k+\stackrel{~}{p_k}=1`$, is automatically enforced by our choice of symbols. These ingredients comprise a three-state random walk model with generalized asymmetric potential and we call it a noise-free MBR.

It is straight forward to set up the following difference equations for the probability distributions of the noise-free MBR model:

$$P_k(t+1)=P_{k+1}(t)\stackrel{~}{p}_{k+1}+P_{k+2}(t)p_{k+2}$$
(1)
for all cyclic (modulo-3) state indices $`k`$. $`P_k(t)`$ is the probability for the random walker at time $`t`$ to be on the state of $`S_k`$. This can be written in matrix form as $`𝐏_{t+1}=𝐏_tB`$, where $`𝐏_t`$ is the time varying probability (row) vector at time $`t`$ and $`B`$ is the transition probability matrix. We can write:

$$\left[B_{i,j}\right]=\left[\begin{array}{ccc}0& p_0& \stackrel{~}{p}_0\\ \stackrel{~}{p_1}& 0& p_1\\ p_2& \stackrel{~}{p_2}& 0\end{array}\right].$$

The steady-state probability, after a sufficiently long time, $`lim_t\mathrm{}𝐏_𝐭=𝐏`$ is simply given as

$$𝐏=𝐏B$$
(2)
which is a characteristic value problem. A partial probability current, $`I`$, can be defined as

$$I=P_kp_kP_{k+1}\stackrel{~}{p}_{k+1}.$$
(3)
If $`I=0`$ there is no net current and detailed balance onsager\_1931  is satisfied, otherwise there exists a net current and the system will assume a non-equilibrium steady-state.

Solving the Eq. 2 together with the normalization condition, $`P_0+P_1+P_2=1`$ , is again straightforward. Using the standard methods for characteristic value or eigenvalue problems, we obtain

$$P_k=(\stackrel{~}{p}_{k+1}+p_{k+1}p_{k+2})/D$$
(4)
for all $`k`$. The denominator $`D`$ is given as

$$D=2+p_0p_1p_2+\stackrel{~}{p_0}\stackrel{~}{p_1}\stackrel{~}{p_2}.$$
(5)
These expressions are consistent with the results of Pearce pearce2000 . It is easy to check that they are the solution to Eq. 2 by direct substitution.

We can substitute the results from Eqs. 4 into Eq. 3 to solve for the net current, $`I`$,

$$I=(p_0p_1p_2\stackrel{~}{p_0}\stackrel{~}{p_1}\stackrel{~}{p_2})/D.$$
(6)
The condition for detailed balance $`I=0`$ is then

$$p_0p_1p_2=\stackrel{~}{p_0}\stackrel{~}{p_1}\stackrel{~}{p_2},$$
(7)
which is the equation of a two-dimensional surface in the three dimensional parameter space, $`\{p_0,p_1,p_2\}`$. Note that the Eq. 6 is independent of state index $`k`$ as is required by its definition given in Eq. 3.

The second part of the MBR is to introduce additional random noise to the system. To the noise-free MBR, we add more noise, controlled by the parameter $`\gamma `$, to the MBR as follows. With a probability of $`\stackrel{~}{\gamma }`$, a random walker follows the dynamic rule of the noise-free MBR otherwise, with the probability of $`\gamma \left(=1\stackrel{~}{\gamma }\right)`$, the walker randomly takes a right or left step with the equal probability of a half. For $`\gamma =0`$ the model is exactly same with the noise-free MBR. In the other limit, for $`\gamma =1`$, the randomizing process dominates and the system reduces to a simple unbiased random walk where the net current remains zero. It is important to note that $`\gamma `$ influences the level of noise in the ratchet but is not identical with the noise itself. We refer to $`\gamma `$ as a ‘noise parameter’. With this modification the transition probability matrix $`B`$ changes as

$$\left[B_{i,j}\right]=\left[\begin{array}{ccc}0& \stackrel{~}{\gamma }p_0+\gamma /2& \stackrel{~}{\gamma }\stackrel{~}{p_0}+\gamma /2\\ \stackrel{~}{\gamma }\stackrel{~}{p_1}+\gamma /2& 0& \stackrel{~}{\gamma }p_1+\gamma /2\\ \stackrel{~}{\gamma }p_2+\gamma /2& \stackrel{~}{\gamma }\stackrel{~}{p_2}+\gamma /2& 0\end{array}\right].$$

From the transition matrix, we know that adding the random noise with parameter $`\gamma `$ effectively changes the existing parameters as

$$p_k\stackrel{~}{\gamma }p_k+\gamma /2$$
(8)
and same holds for the $`\stackrel{~}{p_k}`$’s. The steady state solution and net current for noisy MBR can be obtained by exchanging all the $`p_k`$’s in Eqs. 4 and 6 according to Eq. 8. The expression for current is given as

$$I_\gamma =[\stackrel{~}{\gamma }^3A_{}+\stackrel{~}{\gamma }(\gamma /2)(1\gamma /2)B]/D_\gamma $$
(9)
where

$`D_\gamma `$ $`=`$ $`2+\stackrel{~}{\gamma }^2A_++\stackrel{~}{\gamma }(\gamma /2)(1+\gamma /2)`$
$`A_\pm `$ $`=`$ $`p_0p_1p_2\pm \stackrel{~}{p_0}\stackrel{~}{p_1}\stackrel{~}{p_2}`$ (10)
$`B`$ $`=`$ $`p_0+p_1+p_2\stackrel{~}{p_0}\stackrel{~}{p_1}\stackrel{~}{p_2}.`$
It is possible to further restrict the choices of $`\{p_0,p_1,p_2\}`$ without losing the important properties of the ratchet. Parrondo’s original definition imposed the further constraints $`p_0=q`$ and $`p_1=p_2=p`$. This reduced the parameter space to a two dimensional space with parameters $`\{p,q\}`$. In two dimensional $`\{p,q\}`$ parameter space, the condition of detailed balance, i.e., Eq. 7, gives the equation for a curve that we call the null-curve:

$$q=\frac{\stackrel{~}{p}^2}{p^2+\stackrel{~}{p}^2}=\frac{(1p)^2}{p^2+(1p)^2}.$$
(11)
The null-curve is a special case of the more general null-surface or null-hypersurface, in higher dimensions.

Fig. 2 shows the ‘positive’ and ‘negative’ net current regions of the noise-free MBR. Note that as expected from the symmetry of the system the curve is invariant under the transformations $`q(1q)`$ and $`p(1p)`$. This also apparent from a consideration of Eq. 7.

On the null-surface, we introduce additional random noise to the system by controlling the value of $`\gamma `$. For $`\gamma =0`$ the model is exactly same as the noise-free MBR and the net current remains zero since we are on the null-surface. In the other limit, for $`\gamma =1`$, the randomizing process dominates and the system reduces to a simple unbiased random walk where the net current is also zero. However, counter-intuitively, for $`0<\gamma <1`$ non-zero current is induced by introducing random noise controlled by $`\gamma `$.

In Fig. 3 we show the current versus noise parameter, $`\gamma `$, for different values of parameters $`p`$ and $`q=\stackrel{~}{p}^2/(p^2+\stackrel{~}{p}^2)`$. As $`\gamma `$ is increased from zero, the current increases to a maximum and then falls off, which has the form of stochastic resonance allison01b . The position of this extremum can be obtained from the condition, $`I_\gamma /\gamma =0`$. $`\gamma `$ varies a little from $`0.408`$ for $`p=0.1`$ to $`0.423`$ for $`p=0.5`$.

Fig. 4 shows the net current, $`I_\gamma `$, versus noise $`\gamma `$ when the system is not on the null-surface any more. In the off-balance region, $`q\stackrel{~}{p}^2/(p^2+\stackrel{~}{p}^2)`$, the net current is not zero for $`\gamma =0`$ but still should be zero for $`\gamma =1`$ and the intermediate behavior is qualitatively the same as the balanced behavior. The actual values of $`p`$ and $`q`$ for the various curves in Fig. 4 are in linear increments of $`0.01`$ for $`p`$ and $`q`$. The top curve has parameters $`p=0.77`$ and $`q=0.12`$. The bottom curve has parameters $`p=0.73`$ and $`q=0.08`$.

We generalize the MBR by introducing a bias into the added noise. The walker takes a right step with probability of $`0.5ϵ`$ and a left step with probability of $`0.5+ϵ`$. For $`ϵ0`$ this noise introduces non-zero net current. The new parameter, $`ϵ`$, is essentially a measure of the bias in the added noise.

We can generalize this model to a system of size $`N`$ by repeating the unit cell of modulo-3 $`N`$ times with a periodic boundary condition. In this case, the periodic potential ensures $`p_k(t)=p_{k+3n}(t)n=0,\pm 1,\pm 2,\mathrm{}.`$ Because of the normalization condition, $`_{k=1}^NP_k(t)=1`$, the current will be reduced by factor of $`N`$. Otherwise, the corresponding master equations and solutions are exactly same as the minimal model.

For different moduli, in principle, we can also set up the master equations and solve them exactly by matrix inversion for the set of linear equations. It can be shown that these results have qualitatively the same statistical behavior as the 3-state MBR. Note that for even number moduli there are oscillatory non-stationary solutions.

The transformation in Eq. 8 tells us effectively that the MBR gives the same results as a biased random walk, where the transition probability is not symmetric but biased. Although this analogy can be used to investigate the characteristics of MBR, it is absolutely impossible to determine whether the system is itself a biased random walk or an MBR by analyzing the result of measurements, without prior knowledge that the model is a combination of a balanced unbiased random walk and added random external noise. This makes the MBR valuable for understanding the minimal features of the discrete-time Brownian ratchet. The MBR has applicability in discrete-time processes where the transition probabilities do not fluctuate in time, such as in game and computation theory where transitions occur at precisely defined times.

Funding from GTECH Australasia, Sir Ross and Sir Keith Smith Fund and the Australian Research Council (ARC) is gratefully acknowledged.


text:
# BNL–69052 March 1, 2002 Kaon Physics at BNL To be published in the Proceedings of the The Fifth KEK Topical Conference - Frontiers in Flavor Physics; KEK, Tsukuba, Japan; November 20-22, 2001; S. Hashimoto and T.K. Komatsubara, Eds.

## 1 INTRODUCTION

The AGS has had a broad and rich program in kaon physics. With the recent successful commissioning of the Relativistic Heavy Ion Collider (RHIC), the primary role of the AGS has shifted to become an injector of heavy ions for RHIC. Nevertheless, the AGS remains the highest intensity proton synchrotron in the world and is designed to be available for $``$20 hours/day when not filling RHIC, and as such retains an important role in the US high energy physics (HEP) program. DOE and BNL have approved and agreed to fund one new HEP experiment to run at the AGS between RHIC fills: the E949 experiment which seeks to make a precise measurement of the branching ratio B($`K^+\pi ^+\nu \overline{\nu }`$ ).

### 1.1 $`𝑲\mathbf{}𝝅𝝂\overline{𝝂}`$

The kaon physics program at the AGS is centered on the two golden modes: $`K^+\pi ^+\nu \overline{\nu }`$ and $`K_L^{}\pi ^{}\nu \overline{\nu }`$ . These modes are interesting as there is essentially no theoretical ambiguity in extracting fundamental CKM parameters from measurements of the branching ratios. The theoretical uncertainty in B($`K^+\pi ^+\nu \overline{\nu }`$ ) is $``$7% and is even smaller in B($`K_L^{}\pi ^{}\nu \overline{\nu }`$ ), only $``$2%; in both cases the hadronic matrix element is extracted from the $`K^+\pi ^{}e^+\nu _e`$ ($`K_{e3}`$) decay rate.

The unitarity of the CKM matrix can be expressed as

$$V_{us}^{}V_{ud}+V_{cs}^{}V_{cd}+V_{ts}^{}V_{td}=\lambda _u+\lambda _c+\lambda _t=0$$

with the three vectors $`\lambda _iV_{is}^{}V_{id}`$ converging to form a very elongated triangle in the complex plane. The length of first vector $`\lambda _u=V_{us}^{}V_{ud}`$ is precisely determined from $`K_{e3}`$ decay. The height of the triangle, $`Im\lambda _t`$, can be measured by $`K_L^{}\pi ^{}\nu \overline{\nu }`$ and the length of the third side $`\lambda _t=V_{ts}^{}V_{td}`$ is measured by $`K^+\pi ^+\nu \overline{\nu }`$ . Measurements of the two $`K\pi \nu \overline{\nu }`$ modes, along with the well known $`K_{e3}`$, will completely determine the unitarity triangle.

Comparison of CKM parameters as measured from the golden $`K\pi \nu \overline{\nu }`$ , $`B_d^{}\psi K_S^{}`$ and $`\mathrm{\Delta }M_{B_d}/\mathrm{\Delta }M_{B_s}`$ modes, provide the best opportunity to over-constrain the unitary triangle and to search for new physics. In particular, comparisons of

* $`|V_{td}\text{ }|`$ from $`K^+\pi ^+\nu \overline{\nu }`$ and from the ratio of the mixing frequencies of $`B_d`$ and $`B_s`$ mesons $`\mathrm{\Delta }M_{B_d}/\mathrm{\Delta }M_{B_s}`$  ,
* $`\beta `$ from B($`K_L^{}\pi ^{}\nu \overline{\nu }`$ )/B($`K^+\pi ^+\nu \overline{\nu }`$ ) and from the time dependent asymmetry in the decay $`B_d^{}\psi K_S^{}`$

offer outstanding opportunities to explore the Standard Model (SM) picture of CP-violation.

The SM prediction for the $`K\pi \nu \overline{\nu }`$ branching ratios are B($`K^+\pi ^+\nu \overline{\nu }`$ ) = $`(0.75\pm 0.29)\times 10^{10}`$ and B($`K_L^{}\pi ^{}\nu \overline{\nu }`$ ) = $`(0.26\pm 0.12)\times 10^{10}`$ . In addition, an upper limit on B($`K^+\pi ^+\nu \overline{\nu }`$ ) can be derived in a theoretically unambiguous way from the current limit on $`B_s`$ mixing; this limit is B($`K^+\pi ^+\nu \overline{\nu }`$ ) $`<1.32\times 10^{10}`$ .

## 2 SUMMARY OF RECENT RESULTS

A large number of new results are available from several recently completed rare kaon decay experiments at BNL (running during 1995–98). There were three experiments running during this period — E787: Search for $`K^+\pi ^+\nu \overline{\nu }`$ , E865: Search for $`K^+\pi ^+\mu ^+e^{}`$ and E871: Search for $`K_L^{}\mu e`$ .

### 2.1 Lepton Flavor Violation Searches

Two of these experiments, E865 and E871, made use of the tremendous kaon flux available at the AGS to push the sensitivity of searches for physics beyond the SM (BSM), in particular, to search for lepton flavor violating decays. These experiments have, or will soon, reach single event sensitivities of $``$$`10^{12}`$ and have excluded or severely limited many possible extensions to the SM. A summary of limits from these searches for BSM physics is provided in Table 1. The E871 experiment is finished with all data analysis, but new results on $`K^+\pi ^+\mu ^+e^{}`$ and other modes are expected from E865 in the near future.

### 2.2 Chiral Perturbation Theory

A number of other modes have been observed for the first time, or with significantly larger statistics than was previously available. These modes are of interest as a testing ground for Chiral Perturbation Theory (ChPT). A summary of these measurements is provided in Table 2.

A byproduct of the E871 search for $`K_L^{}\mu e`$ , and a demonstration of the capabilities of the experiment is the observation of the $`K_L^{}e^+e^{}`$ decay, with the smallest BR ever measured for an elementary particle decay. Along with $`K_L^{}e^+e^{}`$ , a very high statistics measurement of the until recently considered ‘rare’ decay $`K_L^{}\mu ^+\mu ^{}`$ was made. This mode is very interesting due to its historical significance, and to the fact that the short distance component of this decay is sensitive to $`\rho `$ or $`Re\lambda _t`$. The dominant contribution to this decay is from two real photons and is well known from QED and the measured value of B($`K_L^{}\gamma \gamma `$ ). It may be possible, with additional theoretical work  and/or new data (for example, from $`K_L^{}\mu ^+\mu ^{}e^+e^{}`$ from KTeV at FNAL) to determine the long distance dispersive contribution to this mode, although there is some controversy on this point . If all long distance contributions can be calculated, this measurement will be able to provide an independent constraint on $`\lambda _t`$. A summary of the $`K_L^{}\mathrm{}^+\mathrm{}^{}`$ modes can be found in Table 2.

### 2.3 CKM Matrix and CP-violation

The third experiment, E787, was designed to search for the decay $`K^+\pi ^+\nu \overline{\nu }`$ . This experiment made use of the tremendous flux of low energy kaons to measure this very rare mode (branching ratio $`10^{10}`$) with only one detectable particle in the final state. A summary of $`\pi ^+\nu \overline{\nu }`$ results is provided in Table 3.

The first $`K^+\pi ^+\nu \overline{\nu }`$ signal was observed in the 1995 data sample of the E787 experiment . No new events were seen in the data sample from 1996–97 , and with a background of $`0.08\pm 0.02`$ events and a signal of one event a branching ratio of $`1.5_{1.2}^{+3.4}\times 10^{10}`$ was measured. That event was in fact in a very clean region of the signal box with a standard model signal to background ratio of 35. An analysis of the final E787 data sample from the last run in 1998 has recently been reported . With a measured background of $`0.066_{0.025}^{+0.044}`$, one new event was observed. The final plot of range vs. energy from the combined E787 1995–98 data sample for events passing all other cuts is shown in Figure 1.

The branching ratio, as determined from these two events is $`1.57_{0.82}^{+1.75}\times 10^{10}`$. This branching ratio is a factor of two higher than expected in the SM and is higher than allowed by the current limit on $`B_s`$ mixing. Of course, the uncertainty on the BR measurement is large due to limited statistics and new data from the E949 experiment are eagerly awaited.

The new event found in the 1998 data sample is in a relatively clean region of the accepted signal region: the SM signal to background ratio for this event is 3.6. An event display is shown in Figure 2.

From the 90% CL limits $`0.56\times 10^{10}<B(K^+\pi ^+\nu \overline{\nu }\text{ })<3.9\times 10^{10}`$ a limit on the branching ratio of the neutral mode can be derived

$`B(K_L^{}\pi ^{}\nu \overline{\nu }\text{ })`$ $`<`$ $`4.4\times B(K^+\pi ^+\nu \overline{\nu }\text{ })`$
$`<`$ $`1.7\times 10^9(90\%\mathrm{CL}).`$
Limits on $`|V_{td}|`$ and $`\lambda _t`$ can be obtained (these are 1-$`\sigma `$ limits except for $`Im\lambda _t`$ which is 90% CL),

$`0.007<`$ $`|V_{td}|`$ $`<0.030,`$
$`2.9\times 10^4<`$ $`|\lambda _t|`$ $`<1.2\times 10^3,`$
$`0.88\times 10^3<`$ $`Re\lambda _t`$ $`<1.2\times 10^3,`$
$`Im\lambda _t`$ $`<1.1\times 10^3.`$
Even with the large statistical error, this new measurement of B($`K^+\pi ^+\nu \overline{\nu }`$ ) provides a non-trivial contribution to global fits of the CKM parameters . The constraints on $`\lambda _t`$ from this result are shown in Figure 3. The constraints from the other golden B modes, $`\mathrm{\Delta }M_{B_d}/\mathrm{\Delta }M_{B_s}`$ and $`B_d^{}\psi K_S^{}`$ are shown on the same plot. One can immediately see that the allowed region is tightly constrained to a narrow crescent by $`K^+\pi ^+\nu \overline{\nu }`$ and $`\mathrm{\Delta }M_{B_d}/\mathrm{\Delta }M_{B_s}`$ .

New data from the successor to E787, E949 will make a significant contribution to our knowledge of the CKM parameters.

In addition, E787 has searched for the decay $`K^+\pi ^+\nu \overline{\nu }`$ in the pion kinematic region below the $`K^+\pi ^+\pi ^{}`$ ($`K_{\pi 2}`$) peak . This region contains more of the $`K^+\pi ^+\nu \overline{\nu }`$ phase space, but is complicated by a significant background from $`K^+\pi ^+\pi ^{}`$ decays with the $`\pi ^+`$ scattering in the scintillating fiber target and down-shifting its kinematics into the search region. The data from the 1996 run of E787 has been analyzed ($``$20% of the entire E787 data sample). One event was observed in the search region, consistent with the background estimate of $`0.73\pm 0.18`$. This implies an upper limit on B($`K^+\pi ^+\nu \overline{\nu }`$ ) $`<4.2\times 10^9`$ (90% C.L.), and is consistent with the 2 events observed above the $`K_{\pi 2}`$ peak and the SM spectrum. Some additional reduction of the background levels in the remaining E787 data may be possible, but the major focus will shift to the new E949 experiment, which has significantly enhanced photon veto capabilities that will further suppress this background. In addition, the next experiment after E949, CKM at FNAL, will be essentially free of this background since there is no stopping target.

## 3 CURRENT PROGRAM

The current high energy physics program at BNL consists of one experiment: E949 — A measurement of the branching ratio B($`K^+\pi ^+\nu \overline{\nu }`$ ). This is the first of the AGS-2000 experiments to be approved and is to be funded by DOE for 6000 hours of running beginning in 2002. Unfortunately, the FY02 run will fall substantially short of the planned running time, so the expectation is that E949 will need to run beyond FY03.

### 3.1 E949

E949 is an upgraded version of the E787 experiment, planning to capitalize on the full AGS beam to collect $`K^+\pi ^+\nu \overline{\nu }`$ data at 14 times the rate of the E787 run in 1995. The new detector has substantially upgraded photon veto capabilities, enhanced tracking, triggering, monitoring, and DAQ capability, and will run at a higher AGS duty factor and a lower kaon momentum (with an increased fraction of useful, stopped kaons). It has been designed to reach a sensitivity of at least 5 times beyond E787 and observe of 5–10 SM events. The background level for E949 measurement of B($`K^+\pi ^+\nu \overline{\nu }`$ ) above the $`K_{\pi 2}`$ peak is reliably projected from E787 data to be $``$10% of the Standard Model signal.

E949 should see up to 10 SM events (or 20 events at the branching ratio measured by E787) within the next couple of years. This is an exciting opportunity to make a significant contribution to quark mixing and CP-violation that should be fully exploited. A history of the search for $`K^+\pi ^+\nu \overline{\nu }`$ is shown in Figure 4.

The next step towards a precision measurement of B($`K^+\pi ^+\nu \overline{\nu }`$ ) will be the CKM experiment at FNAL. CKM has been given scientific (Stage–1) approval by FNAL and could be running by 2007. CKM plans to use a novel technique for $`K^+\pi ^+\nu \overline{\nu }`$ : a decay in flight experiment, with redundant kinematic constraints from a reasonably conventional momentum spectrometer and a novel velocity spectrometer based on RICH counters. CKM expects to observe 100 SM signal events in a two year run, using Main Injector pulses not needed by the Tevatron. CKM will require less than 20% of the flux from the Main Injector, but will require a slow extracted spill of $``$1 second. Figure 4 shows a projected measurement of B($`K^+\pi ^+\nu \overline{\nu }`$ ) from CKM, assuming the current branching ratio.

## 4 FUTURE PLANS

Other HEP experiments that are under consideration to run at the AGS include additional running for measurement of g-2 of the muon, a search for muon to electron conversion (MECO) and a search for and measurement of the $`K_L^{}\pi ^{}\nu \overline{\nu }`$ branching ratio (K0PI0). The National Science Board of the National Science Foundation (NSF) has has approved the construction of the two new large experiments: K0PI0 and MECO, as components of the Rare Symmetry Violation Proposal (RSVP). RSVP is planned to be one of the next Major Research Equipment construction projects at the NSF.

### 4.1 K0PI0

The KOPIO experiment is designed to discover the $`K_L^{}\pi ^{}\nu \overline{\nu }`$ decay and measure its branching ratio to $``$20%. KOPIO will make use of a time-of-flight technique to measure the momentum of the $`K_L`$ and will operate at a large targeting angle to improve the $`p_K`$ resolution and soften the neutron spectrum to reduce $`\pi ^{}`$ hadroproduction. KOPIO will have a substantial photon veto system and make the same sort of background measurements as E787, directly from the data, with two independent tools for attacking the major background, $`K_L^{}\pi ^{}\pi ^{}`$ , through both kinematics and photon veto. K0PI0 expects to observe 50 SM events, with a background of 50%. This will allow a determination of $`Im\lambda _t`$ to 10%. K0PI0 is expected to start data collection in $``$2006.

## 5 CONCLUSION

The next decade will be an exciting time for improved understanding of CP-violation and quark mixing. It is quite likely that precise measurements of all four “Golden Modes” will be made: $`B_d^{}\psi K_S^{}`$ at the B-factories, $`\mathrm{\Delta }M_{B_d}/\mathrm{\Delta }M_{B_s}`$ most likely at the Tevatron, and $`K\pi \nu \overline{\nu }`$ at BNL, KEK and FNAL. These measurements will allow a precise determination of CKM parameters and provide a critical test of the SM picture of CP-violation.