format: text

text:
# Matching Conditions on Capillary Ripples: Polarization Estudio de las condiciones de empalme para las oscilaciones de intercara entre dos fluidos. Polarización.

## I Introduction

The boundary conditions constitute the key feature in any theory of surface waves. It is through them that one introduces in the analysis the physical consequences of specific surface effects, such as local changes of mass or stresses, thus going beyond the simple case in which one merely matches two semi-infinite bulk media at an interface (which, in particular, can be a free surface).

In the last years, increased attention has been paid to the properties of capillary waves by physicists and chemists . Ripples represent one of the few cases in which the relation between the dynamical properties of a surface and liquid flows can be predicted completely. The study of capillary ripples has clarified which properties of a liquid surface determine the surface’s resistance against deformation.

The boundary conditions for the stress at the interface are derived from the principle that the forces acting upon such an “interfacial” element result not only from viscous stresses in the liquid but also from stresses existing in the deformed interface. The cause of the difference is that an interface, unlike a three dimensional liquid, can not enjoy the property of incompressibility . Work done on an element of liquid is partially degrade into heat by viscous friction and partially transformed into kinetic energy which is transmitted to adjoining elements. Work done on an interfacial element leads, at least partially, to an increment of the surface potential energy. It is this potential energy of the deformed interface which enables the whole system, including the interface, to carry out an oscillatory motion.

In this article, first, it will be shown a non-traditional way to obtain the matching conditions at the interface between two non-mixed fluid at rest, considering the equation of motion of the whole media directly. Usually, university courses do not use this approach and state the boundary conditions from outside the constitutive equations governing the studied problem. It is important from a pedagogical point of view to evidence that in fact, the matching conditions at the interfaces are contained, almost in all cases, in the equation of movement for the whole system taken as the composition of all media. The other cases arise when the fluid interfaces have intrinsic properties not included on the equations of the media in that cases mentioned before.

This point of view was already used in  to elaborate a formalism, based on Surface Green Functions, which establish an isomorphism between solids and fluids related to the interface normal modes. The aim of this paper is to emphasize the fact that this way to establish the matching conditions allows us to introduce the physical characteristics related with the boundary conditions in a firmer floating.

## II Equations for the bulk in a viscous incompressible fluid.

To achieve a system of equations which describes the oscillatory motion of a viscous incompressible fluid, the starting point is the linearized Navier-Stokes equation for a viscous fluid , which reads:

$$\rho \frac{V_i}{t}=\frac{\tau _{ij}}{x_j}+\rho F_{i,\text{ext}}$$
(1)
where $`\rho `$ is the density of equilibrium, $`V_i`$ the components of the fluid particle, $`F_{i,\text{ext}}`$ are the external forces and $`\tau _{ij}`$ is the stress tensor which, for a viscous incompressible fluid, has the following form

$$\tau _{ij}=p\delta _{ij}+\eta \left(\frac{V_i}{x_j}+\frac{V_j}{x_i}\right)$$
(2)
Here $`p`$ are the small variations of pressure, $`\eta `$ is the viscosity parameter and $`\delta _{ij}`$ is the unit matrix elements. The subscripts $`i`$ and $`j`$ take the values $`y`$ and $`z`$ identically. In (1) the sum over repeated subscripts is understood. It had been assumed that the system is symmetric with respect to the $`x`$ direction.

Putting (2) in (1) it is obtained

$`\rho {\displaystyle \frac{V_y}{t}}+{\displaystyle \frac{}{y}}\left(p+2\eta {\displaystyle \frac{V_y}{y}}\right)+`$ (3)
$`+{\displaystyle \frac{}{z}}\left[\eta \left({\displaystyle \frac{V_y}{z}}+{\displaystyle \frac{V_z}{y}}\right)\right]`$ $`=`$ $`0`$ (4)
$`\rho {\displaystyle \frac{V_z}{t}}+{\displaystyle \frac{}{y}}\left[\eta \left({\displaystyle \frac{V_z}{y}}+{\displaystyle \frac{V_y}{z}}\right)\right]+`$ (5)
$`+{\displaystyle \frac{}{z}}\left(p+2\eta {\displaystyle \frac{V_z}{z}}\right)`$ $`=`$ $`0`$ (6)
where we have neglected the external forces.

Also the continuity equation is needed, which expresses that the volume of an element of the incompressible fluid does not change during the motion. It has the form

$$\frac{V_y}{y}+\frac{V_z}{z}=0$$
(7)
The first step is to obtain the equation of motion and its solution of the each medium taken as infinite. So, the parameters as density and viscosity are considered constants in the whole medium and it leads to transform equations (4) and (6) to:

$`\rho {\displaystyle \frac{V_y}{t}}{\displaystyle \frac{p}{y}}+\eta ^2V_y`$ $`=`$ $`0`$ (8)
$`\rho {\displaystyle \frac{V_z}{t}}{\displaystyle \frac{p}{z}}+\eta ^2V_z`$ $`=`$ $`0`$ (9)
where $`^2=(^2/y^2+^2/z^2)`$.

The solution of the system (7), (8) and (9) written as a vector field velocity, can be putted as the sum of an irrotational field (related with the longitudinal mode) and a divergence free field (related with the transverse modes) , i.e.:

$$𝑽=𝑽_1+𝑽_2$$
(10)
which satisfy:

$`\times 𝑽_1`$ $`=`$ $`0`$ (11)
$`𝑽_2`$ $`=`$ $`0`$ (12)
Any irrotational field is characterized by a scalar function, the “potential” function $`\phi (y,z,t)`$, such that

$$𝑽_1=\phi $$
(13)
and the divergence free field can be described by a vector function , the “stream” or vorticity function $`\psi (y,z,t)`$, such that

$$𝑽_2=(\frac{\psi }{z},\frac{\psi }{y})$$
(14)
The velocity components can thus be written in terms of the potential and the stream functions, as:

$`V_y`$ $`=`$ $`{\displaystyle \frac{\phi }{y}}{\displaystyle \frac{\psi }{z}}`$ (15)
$`V_z`$ $`=`$ $`{\displaystyle \frac{\phi }{z}}+{\displaystyle \frac{\psi }{y}}`$ (16)
Substitution of eqns. (15) and (16) in the continuity condition (7) gives

$$^2\phi =0$$
(17)
and substitution into eqns. (8) and (9) leads to

$`{\displaystyle \frac{}{y}}\left\{\rho {\displaystyle \frac{\phi }{t}}p\right\}+{\displaystyle \frac{}{z}}\left\{\rho {\displaystyle \frac{\psi }{t}}\eta ^2\psi \right\}`$ $`=`$ $`0`$ (18)
$`{\displaystyle \frac{}{z}}\left\{\rho {\displaystyle \frac{\phi }{t}}p\right\}{\displaystyle \frac{}{y}}\left\{\rho {\displaystyle \frac{\psi }{t}}\eta ^2\psi \right\}`$ $`=`$ $`0`$ (19)
Equations (18) and (19) are simultaneously satisfied if one considers:

$`\rho {\displaystyle \frac{\phi }{t}}p`$ $`=`$ $`C_1`$ (20)
$`\rho {\displaystyle \frac{\psi }{t}}\eta ^2\psi `$ $`=`$ $`C_2`$ (21)
The constants $`C_1`$ and $`C_2`$ are obtained from the condition at zero flow, and it gives raise to $`C_1=p_o`$ and $`C_2=0`$ where $`p_o`$ is the reference atmospheric pressure.

The solution of equations (17), (20) and (21) are looking for as the following

$`\phi (y,z,t)`$ $`=`$ $`\mathrm{\Phi }(z)\text{e}^{i\left(\kappa y\omega t\right)}`$ (22)
$`\psi (y,z,t)`$ $`=`$ $`\mathrm{\Psi }(z)\text{e}^{i\left(\kappa y\omega t\right)}`$ (23)
where the parameters $`\kappa `$ and $`\omega `$ are the wavevector and frequency of the wave respectively.

Substituting eqns. (22) and (23) in (17) and (21) gives the $`z`$-dependence of the functions $`\phi `$ and $`\psi `$, which satisfy:

$`{\displaystyle \frac{\text{d}^2\mathrm{\Phi }(z)}{\text{d}z^2}}\kappa ^2\mathrm{\Phi }(z)`$ $`=`$ $`0`$ (24)
$`{\displaystyle \frac{\text{d}^2\mathrm{\Psi }(z)}{\text{d}z^2}}\text{q}_\text{t}^2\mathrm{\Psi }(z)`$ $`=`$ $`0`$ (25)
with $`\text{q}_\text{t}^2=\kappa ^2i\rho \omega /\eta `$

Equations (24) and (25) lead to the solutions:

$`\mathrm{\Phi }(z)`$ $`=`$ $`C_1\text{e}^{\kappa z}+C_2\text{e}^{\kappa z}`$ (26)
$`\mathrm{\Psi }(z)`$ $`=`$ $`C_3\text{e}^{\text{q}_\text{t}z}+C_4\text{e}^{\text{q}_\text{t}z}`$ (27)
and in combination with eqns. (22) and (23) they give a solution of the form:

$`\phi (y,z,t)`$ $`=`$ $`\left(C_1\text{e}^{\kappa z}+C_2\text{e}^{\kappa z}\right)\text{e}^{i\left(\kappa y\omega t\right)}`$ (28)
$`\psi (y,z,t)`$ $`=`$ $`\left(C_3\text{e}^{\text{q}_\text{t}z}+C_4\text{e}^{\text{q}_\text{t}z}\right)\text{e}^{i\left(\kappa y\omega t\right)}`$ (29)
where $`C_1,C_2,C_3`$ and $`C_4`$ are constants to be determined by boundary and matching conditions.

Then, the expressions for the varying velocity components and the pressure are finally obtained by substitution of eqns. (28) and (29) into (15), (16) and (20), respectively.

## III Interface problem: matching conditions.

Now we will match the two media. We consider a surface which at rest coincides with the plane $`z=0`$ and it separates medium $`M_1`$ at $`z<0`$ from medium $`M_2`$ at $`z>0`$. Each one is viscous and incompressible.

First of all, any solution has to fulfill the continuity of the velocity field across the surface according to

$`V_y^{(1)}`$ $`=`$ $`V_y^{(2)}\text{at}z=0`$ (30)
$`V_z^{(1)}`$ $`=`$ $`V_z^{(2)}\text{at}z=0`$ (31)
Superscript $`1`$ denotes medium $`M_1`$ and superscript $`2`$ denotes medium $`M_2`$. The rest of the matching conditions at the interface are derived from the system of equations which govern the whole system. This is no usually done in the normal program courses at the Universities and we consider that it is important to state that the matching conditions are, almost in all of the cases, contained in the equation of movement for the whole system taken as the composition of each one. These are eqns. (4) and (6) where the parameters $`\eta `$ and $`\rho `$ are constant, but taking different values on each medium. Integrating these equations through the surface about $`z=0`$ from $`ϵ`$ to $`+ϵ`$ and later taking $`ϵ0`$, it is obtained from eq. (4)

$`\left[\eta _1\left({\displaystyle \frac{V_y^{(1)}}{z}}+{\displaystyle \frac{V_z^{(1)}}{y}}\right)\right]_{z=0^{}}=`$ (32)
$`=\left[\eta _2\left({\displaystyle \frac{V_y^{(2)}}{z}}+{\displaystyle \frac{V_z^{(2)}}{y}}\right)\right]_{z=0^+}`$ (33)
and from eq. (6)

$`\left[p_1+2\eta _1{\displaystyle \frac{V_z^{(1)}}{z}}\right]_{z=0^{}}=`$ (34)
$`=\left[p_2+2\eta _2{\displaystyle \frac{V_z^{(2)}}{z}}\right]_{z=0^+}p_\gamma `$ (35)
where $`p_\gamma `$ is the jump due to the surface tension according with the Laplace Law . The other terms in eqns. (4) and (6) vanish when $`ϵ0`$ because they are continuous or have a finite jump at $`z=0`$.

It is easily seen according to eq. (2) that eqns. (33) and (35) are the conditions of the continuity of the stress tensor components through the interface, as expected.

On the other hand, the boundary conditions of the problem are regularity at $`z\pm \mathrm{}`$, which leads to eliminate 4 of the 8 constants appearing in (28) and (29) for the two media. Using the resulting functions $`\phi `$ and $`\psi `$ for each medium in conditions (30)-(35) the following system is met:

$$\begin{array}{cccc}i\kappa & i\kappa & \text{q}_{\text{t1}}& \text{q}_{\text{t2}}\\ 1& 1& i& i\\ 2i\kappa ^2\eta _1& 2i\kappa ^2\eta _2& \eta _1Q_{t1}& \eta _2Q_{t2}\\ d_{41}& \eta _2Q_{t2}& d_{43}& 2i\kappa \eta _2\text{q}_{\text{t2}}\end{array}\begin{array}{c}C_1\\ C_2\\ C_3\\ C_4\end{array}=0$$
(36)
to determine the remaining constants, with the definitions $`d_{41}=\eta _1Q_{t1}i\mathrm{{\rm Y}}`$, $`d_{43}=2i\kappa \eta _1\text{q}_{\text{t1}}\mathrm{{\rm Y}}`$, $`\mathrm{{\rm Y}}=\gamma \kappa ^3/\omega `$, $`Q_{t1}=2\kappa ^2i\rho _1\omega /\eta _1`$ and $`Q_{t2}=2\kappa ^2i\rho _2\omega /\eta _2`$.

In deducing the expression in system (36) which comes from (35) it was considered that

$$p_\gamma =\gamma \frac{^2z}{y^2}$$
(37)
on the surface, but as we are dealing with the velocity field, then it is obtained that:

$$\frac{p_\gamma }{t}=\gamma \frac{^2V_z}{y^2}$$
(38)
which finally leads to:

$$p_\gamma =\frac{\gamma \kappa ^2}{i\omega }\left(\kappa C_1+i\kappa C_3\right)\text{e}^{i\left(\kappa y\omega t\right)}$$
(39)
according to the substitution $`/ti\omega `$ and $`/y\kappa ^2`$

The vanishing of the determinant of (36) leads to the equation for the dispersion relation of the existing surface modes. It can be written as:

$`\omega ^2\left[(\rho _1+\rho _2)(\rho _1\text{q}_{\text{t2}}+\rho _2\text{q}_{\text{t1}})\kappa (\rho _1\rho _2)^2\right]`$ $`+`$ (40)
$`+\gamma \kappa ^3\left[\rho _1(\kappa \text{q}_{\text{t2}})+\rho _2(\kappa \text{q}_{\text{t1}})\right]`$ $`+`$ (41)
$`+\mathrm{\hspace{0.33em}4}\kappa ^3(\eta _2\eta _1)^2(\kappa \text{q}_{\text{t1}})(\kappa \text{q}_{\text{t2}})`$ $`+`$ (42)
$`+\mathrm{\hspace{0.33em}4}i\kappa ^2\omega (\eta _2\eta _1)(\rho _1\kappa \rho _2\kappa \rho _1\text{q}_{\text{t2}}+\rho _2\text{q}_{\text{t1}})`$ $`=`$ $`0`$ (43)
This is the dispersion relation for capillary waves for the interface of two viscous incompressible fluids.

In order to obtain the constants $`C_1,C_2,C_3`$ and $`C_4`$, an initial stimulus is needed according to an initial value problem  but this method, simple at the beginning, becomes rather complicated later and it is not good for a quite general study.

As was said on the introduction, the aim of this paper is not to get inside the dispersion relation of the normal modes at the interface of two viscous fluids at rest. For a better study of this subject we recommend paper . We get here to show a way of solution also different from the usually taken as only the velocity vector as a function of an ”stream function”. From now on, we will put our attention on the matching conditions and we will show that more that a mathematical information of the matching can be found on it, but also the physics of the polarization can be deduced and how the interface moves in its oscillation.

From eqns. (33) and (35) it can be obtained more information about the polarization of the modes on the interface. This will be done in the next section.

## IV Matching conditions and polarization.

There are two possible modes on the fluid: one in which the fluid particle moves in the direction of the wave propagation called longitudinal with notation L($`V_y`$,$`0`$) and another transverse to the direction of the wave propagation and normal to the interface with notation TN($`0`$,$`V_z`$).

The longitudinal mode is related to the fluid compressibility because this motion of the fluid particles is only possible when its volume changes . This analysis also holds when the mode is on the interface, but this does not mean that there are no longitudinal modes of oscillation on the interface when the fluids involved are incompressible. It had been shown in  that the interface, when oscillating, must be considered as a compressible one, because precisely its change in area is responsable for the increasing of its potential energy and therefore, for its oscillation. It is important to state that this is a fundamental argument to understand the movement of any interface in hydrodynamics.

Nevertheless, now it can be shown that in spite of the compressibility of the interface, there does not exist pure longitudinal mode if we are dealing with incompressible media. Let us demonstrate this.

If a point $`y_o`$ on the interface is considered moving with velocity, say $`V_{yo}^S`$ in the $`y`$-axis direction, then, according to the continuity of velocity, the point ($`y_o`$,$`ϵ`$) in $`M_1`$ and the point ($`y_o`$,$`+ϵ`$) in $`M_2`$ must have the same velocity $`V_{yo}^S`$ if $`ϵ`$ is small enough. As the interface is compressible, at the point $`y_1`$ near enough $`y_o`$ the velocity can be, for instance, $`V_{y1}^S`$ different in general from $`V_{yo}^S`$ but as the media are incompressible, at the point ($`y_1`$,$`ϵ`$) and ($`y_1`$,$`+ϵ`$) the velocity must be $`V_{yo}^S`$. See Fig. 1. This is not in conformity with the continuity of the velocity through the surface and hence the pure longitudinal mode is not possible and only the TN mode seems to be valid when the media are incompressible.

After these considerations during the above demonstration, the student can keep the idea that the interface oscillations can only occur in the $`z`$-axis. This is the accurate moment to show to the student that things not always are as they apparently seem to be, because that assumption does not take into account the different properties of each medium, whose response depends on its fundamental parameters, as density and viscosity, which are different for each medium. Hence, it is evident that it must be analyzed, precisely, the interface matching conditions.

It is useful to compare and to support the previous qualitative analysis with a quantitative and more profound one regarding the interface matching conditions.

Recalling carefully eqns. (33) and (35) and supposing that such a wave propagates in $`y`$ direction with movement only in $`z`$ direction (TN mode), then $`V_y^{(1)}=V_y^{(2)}=0`$ and eq. (33) becomes

$$\eta _1\frac{V_z^{(1)}}{y}|_{z=0^{}}=\eta _2\frac{V_z^{(2)}}{y}|_{z=0^+}$$
(44)
It is known that $`V_z`$ is continuous along the interface for all points. Then, the derivative with respect to $`y`$ is also the same in both hands of (44) and this expression only holds if $`\eta _1=\eta _2`$, i.e., if the media have the same viscosity. It does not mean for the interface to disappear because the density of each medium can be different. Then, if the viscosities of the media have not the same value, the velocity component $`V_y`$ along the interface must be different from zero to compensate the inequality (44) and to fulfill the continuity of the stress tensor in the $`y`$ direction given by eq. (33) yielding to a component of movement along $`y`$ direction. This mode will be called Sagittal mode or S$`(V_y,V_z)`$.

The above analysis was done for the general case. Now we are able to take the particular case in which one of the media is vacuum, for instance, $`M_2`$, with $`\eta _2=0`$. Then, condition (33) becomes

$$\left[\eta _1\left(\frac{V_y^{(1)}}{z}+\frac{V_z^{(1)}}{y}\right)\right]_{z=0^{}}=0$$
(45)
and it can be seen that also $`V_y`$ must be non zero on the surface to hold eq. (45) with the corresponding Sagittal polarization movement.

With this analysis on the conditions of stress component continuity in $`y`$ direction along the interface, it can be seen that if the two media are viscous (at least one of them), the fluid particle of the interface moves in a Sagittal mode which combines movement in both directions: along the wave propagation in $`y`$ direction, and normal to the interface in $`z`$ direction.

Then, it is qualitatively clear that the viscosity of each media plays a fundamental role in the coupling of modes even for incompressible fluids. In spite of that, it could be a mistake to say that the modes decouple if the viscosities are equal. It should not be forgotten that eq. (35) is also important in the characterization of the interface particle behaviour and it includes the pressure on each side of the surface. According to eq. (20), the pressure is associated with the longitudinal mode and the inertial effect of the fluid particle according with the density of the media. Then, it contains the information of each components of the velocity and also of the density and according to eq. (28) the pressure has a jump through the interface. This result, in combination with the analysis of eq. (33) make difficult to understand the role played by the densities of each media on the surface polarization movement, and it can not be reached from this only analysis. This point is still a matter of investigation.

## V Conclusions

The present work is an attempt to give an example, using the hydrodynamics, of how the study of the interface matching conditions allows us to make a plentiful and rich in details discussion. Moreover, of how the interface matching conditions content a sufficient information to conclude that the interface oscillation must be with a Sagittal mode and not with neither a pure longitudinal, nor a pure transversal one. This movement has been shown to be close related to the physical properties of the media such as viscosities and that fact allows us to establish rigorously that those are the parameter which characterize the interface movement and the response of each medium to an stimulus coming from the other one.

It was seen that viscosity is the main parameter in the coupling of the two modes to achieve a Sagittal one, nevertheless within the framework of this formalism it is difficult to determine the role of viscosity and of the density ratios in the coupling of modes. This aspect needs further investigation.


text:
# Switching and symmetry breaking behaviour of discrete breathers in Josephson ladders

(April 4, 2000)

We investigate the roto-breathers recently observed in experiments on Josephson ladders subjected to a uniform transverse bias current. We describe the switching mechanism in which the number of rotating junctions increases. In the region close to switching we find that frequency locking, period doubling, quasi-periodic behaviour and symmetry breaking all occur. This suggests that a chaotic dynamic occurs in the switching process. Close to switching the induced flux increases sharply and clearly plays an important role in the switching mechanism. We also find three critical frequencies which are independent of the dissipation constant $`\alpha `$, provided that $`\alpha `$ is not too large.

The recent discovery of discrete breathers in Josephson ladder arrays (Fig. 1) driven by a transverse dc bias current has shown not only that these localised excitations exist but that they exhibit remarkable behaviour, which is the subject of this paper.

The Josephson ladder has been studied theoretically for many years. In the absence of a driving current the interaction between vortices has been found to be exponential and this leads to the vortex density exhibiting a devil’s staircase as the magnetic field is increased. Quantum fluctuations, meta-stable states and inductance effects have been studied as has the interaction of vortices with a transverse dc bias current.

In the case of a Josephson ladder, a roto-breather is a stable group of vertical junctions rotating together ($`\theta _j^{}\theta _j\omega t`$, see Fig. 1). Until recently, the interest has focussed on the case of a sinusoidal bias current, where it has been shown that roto-breather solutions exist and that a discrete breather may even include a chaotic trajectory.

Roto-breathers should also be stable in the experimentally easier case of a ladder array subjected to a uniform dc transverse driving current. The case of rotation occurring at a single vertical junction has been studied numerically. Two experimental groups have now independently confirmed the existence of such solutions and discovered further unpredicted behaviour. In an annular ladder (Fig. 1a) and in a linear ladder (Fig. 1b) breathers were observed with various numbers of “vertical” (i.e. radial) junctions rotating. However the most interesting features of the observations were the switching between breathers and their symmetry breaking behaviour. Once a breather has been initialised it is stable, but if the bias current is slowly decreased then, at some critical value of the current, the number of rotating junctions switches to a larger number and a new breather is formed. Furthermore, the new breather may have a different centre of symmetry from the old breather despite the symmetry of the ladder itself and the symmetry of the driving currents about a single vertical junction; this is particularly remarkable in the case of the annular ladder where translational invariance should be exact. The production of breathers with an even number of rotating junctions is, by itself, a demonstration of symmetry breaking. The symmetry breaking effect was not mentioned by the experimenters, presumably because it could be due to experimental imperfections. However we show that it should also arise in a perfect experiment, the apparent mechanism being the occurrence of chaotic dynamics in the switching region.

First we construct a new model, similar to but significantly different from those proposed previously. The current through a junction is determined by the RCSJ model

$$\frac{I}{I_c}=\frac{d^2\phi }{dt^2}+\alpha \frac{d\phi }{dt}+\mathrm{sin}\phi $$
(1)
where $`I_c`$ is the critical current and

$$\phi =\mathrm{\Delta }\theta \frac{2\pi }{\mathrm{\Phi }_0}𝐀.\mathrm{d}𝐥$$
(2)
where $`\mathrm{\Delta }\theta `$ is the change in superconducting order parameter $`\theta `$ across the junction and $`𝐀`$ is the vector potential. The “vertical” (i.e. radial) junctions may differ in area from the “horizontal” junctions by an anisotropy parameter $`\eta =I_{ch}/I_{cv}=C_h/C_v=R_v/R_h`$ where $`I_{ch}`$ ($`I_{cv}`$), $`C_h`$ ($`C_v`$) and $`R_h`$ ($`R_v`$) are, respectively, the critical current, capacitance and resistance of a horizontal (vertical) junction. From Fig. 1a we see that there are three unknowns per plaquette: $`\theta _j`$, $`\theta _j^{}`$ and $`f_j`$, where $`f_j`$ is the total flux threading the $`j`$th plaquette. To solve for these unknowns we construct three equations per plaquette as follows. The first two equations are obtained from current conservation at the top (inner) and bottom (outer) rails. The third equation is obtained by making the approximation that the induced flux $`f_jf_a`$ (where $`f_a`$ is the applied flux) is produced solely by the currents flowing around the immediate perimeter of the $`j`$th plaquette (Fig. 1c):

$$\frac{f_jf_a}{\mathrm{\Phi }_0}=\frac{\beta _L}{8\pi }\left(I_j^v+I_j^hI_{j+1}^vI_j^h\right)$$
(3)
where $`\beta _L`$ is an inductance parameter. This last equation makes the plausible assumption that the plaquettes are more or less square and that it is only the currents flowing around the plaquette that produce the induced field. This differs from, and for square plaquettes is more accurate than, the common assumption that the induced field is proportional to the loop (or “mesh”) current circulating the plaquette.

We impose the initial conditions, mimicing the experiments, that at $`t=0`$, $`\theta _j=\theta _j^{}=0`$, $`f_j=0`$ and $`I_j=0`$ for all $`j`$. This means that $`\theta _j+\theta _j^{}=0`$ for all time, i.e. there is a symmetry between the inner and outer rails. Using Landau gauge we then have the following pair of coupled differential equations for each plaquette:

$`\left\{{\displaystyle \frac{d}{dt^2}}+\alpha {\displaystyle \frac{d}{dt}}\right\}\left(\theta _{j1}^{}+4\theta _j^{}\theta _{j+1}^{}\right)=2I_j+\mathrm{sin}\theta _{j1}^{}`$ (4)
$`4\mathrm{sin}\theta _j^{}+\mathrm{sin}\theta _{j+1}^{}+{\displaystyle \frac{8\pi }{\beta _L}}(f_{j1}f_j)`$
$`\left\{{\displaystyle \frac{d}{dt^2}}+\alpha {\displaystyle \frac{d}{dt}}\right\}\left(2\pi \eta f_j+(1\eta )(\theta _{j+1}^{}\theta _j^{})\right)=\mathrm{sin}\theta _j^{}`$ (5)
$`\mathrm{sin}\theta _{j+1}^{}+{\displaystyle \frac{8\pi }{\beta _L}}(f_jf_a)+2\eta \mathrm{sin}\chi _j^{}`$
where $`\theta _j^{}=\theta _j^{}\theta _j`$ and $`\chi _j^{}=\frac{1}{2}(\theta _{j+1}^{}\theta _j^{}+2\pi f_j)`$.

We now focus on determining whether or not our model exhibits the interesting switching and symmetry breaking behaviour observed in the data of Binder $`et`$ $`al`$. All parameter values are chosen to mimic the experimental setup i.e. $`\eta =0.44`$, $`\beta _L=2.7`$, $`\alpha =0.07`$ and $`f_a=0`$. Let

$$I_j=\{\begin{array}{cc}I_B+I_\mathrm{\Delta }\hfill & \text{ if }j=0\hfill \\ I_B\hfill & \text{ otherwise}\hfill \end{array}$$
(6)
where $`I_B`$ is called the bias current. Again following the experiment (Fig. 1a), we slowly increase $`I_\mathrm{\Delta }`$ while keeping $`I_B=0`$ until rotation starts at site $`j=0`$. $`I_\mathrm{\Delta }`$ is then slowly decreased while at the same time increasing $`I_B`$ to keep $`I_0`$ constant. When $`I_B`$ has reached the desired value it is then held fixed while $`I_\mathrm{\Delta }`$ is slowly reduced to zero. Finally $`I_B`$ is slowly reduced while keeping $`I_\mathrm{\Delta }=0`$. Note that the dynamical equations, initial conditions and injected currents are exactly symmetrical about site $`j=0`$. One would expect only solutions which are also symmetrical about $`j=0`$. However, like the experiments, the simulations also produce breathers which are not symmetrical about $`j=0`$. Fig. 2 shows how the number $`N_R`$ of rotating junctions changes as the bias current $`I_B`$ is slowly reduced from various starting values. In agreement with experiment, breathers with even $`N_R`$ are commonly produced (these cannot be symmetrical about $`j=0`$), and $`N_R`$ switches to larger and larger values until eventually all junctions are rotating. While we have found that similar behaviour occurs in a previously published model, our model gives considerably better agreement with the experimentally observed switching currents, thus indicating the importance of inductance effects in the switching mechanism. We believe the origin of the symmetry breaking is the occurrence of chaotic dynamics in the switching region.

The chaotic nature of the switching is further suggested by the fact that when the same computer code which produced Fig. 2 is performed on a different computer (different floating point processor) we see significant changes in the switching behaviour of nearly all trajectories although the overall pattern remains identical. The minute differences in the handling of floating point numbers are amplified in the chaotic regime to produce significantly altered trajectories.

Fig. 3 shows the voltage-current characteristics obtained from the same simulations used in producing Fig. 2. At high frequencies ($`V=d\theta ^{}/dt>4.3`$) most of the breathers show purely resistive behaviour (i.e. $`V=d\theta ^{}/dt=I_BR`$), the value of the resistance $`R`$ increasing with $`N_R`$ according to an expression deduced from experiment:

$$\alpha R1/(1+\eta /N_R)$$
(7)
This is the relationship expected if the current through each rotating vertical junction is $`\alpha V`$ and the four horizontal junctions surrounding the rotating region each carry current $`\frac{1}{2}\eta \alpha V`$ (to satisfy Kirchoff’s rule). For $`N_R=1`$ this expression was derived in ref. . We find that no such resistive behaviour occurs for $`V<4.3`$.

A typical example of a resistive breather (far away from any switching region) is shown in Fig. 4. Although the breather is not symmetric about $`j=0`$ it shows exact symmetry about the midpoint of the rotating region and also appears to be exactly periodic, the period being two revolutions of a vertical junction (we call this “period 2”). The non-rotating junctions oscillate, the amplitude decreasing exponentially with distance from the rotating region. Far away from switching the magnitude of the flux $`f_j`$ is everywhere small ($`<0.1`$) and is limited to the edges of the rotating region.

Fig. 3 shows that there are at least two critical frequencies, $`d\theta ^{}/dt=6.5`$ and $`d\theta ^{}/dt=5.0`$, at which breathers become unstable and switching occurs. As the current is reduced the frequency falls until it reaches the critical value, at which point the breather becomes unstable and switches to a larger breather with a larger resistance (Eq. (7)) and therefore a higher frequency. The process then repeats as the current continues to be ramped down. We identify the larger of these two critical frequencies with that observed experimentally. Although, at fixed $`I_B`$, the frequency of rotation depends on the dissipation constant $`\alpha `$, we find that these two critical frequencies are more or less independent of $`\alpha `$, as is the critical frequency below which breathers cease to show constant resistance. This of course breaks down when the bias current required to achieve the upper critical frequency exceeds unity, i.e. when $`\alpha >0.15`$.

Note also that in the above simulations: (i) no breathers are seen for $`I_B>0.85`$, (ii) only the single site breather is seen for $`0.68<I_B<0.85`$ and (iii) all rotating junctions normally stop together when the current $`I_B`$ is reduced below 0.17, although sometimes a new breather containing only a few rotating junctions is produced which is then destroyed when $`I_B`$ falls below 0.14.

As well as resistive period 2 breathers we have found frequency-locked breathers. In fact the usual $`N_R=1`$ resistive breather makes a transition to a frequency-locked breather as the bias current is reduced below $`I_B=0.699`$. At this point the period doubles, and the breather becomes asymmetric and jumps to a higher rotation frequency $`d\theta ^{}/dt=6.32`$. At $`I_B=0.6895`$ the period doubles again and at $`I_B=0.6892`$ it becomes quasi-periodic (Fig. 5). Switching to multi-site breathers occurs at $`I_B0.6891`$. Such frequency locked breathers have not been reported experimentally but have probably been overlooked as they are only stable in narrow current intervals.

While the maximum flux $`f^{max}`$ is normally small ($`f^{max}<0.1`$), it sharply increase at switching to $`f^{max}1`$. It would appear that an important part of the switching mechanism is the instability caused by having a large flux concentrated in a small area. The importance of the flux in the switching process is further confirmed by the fact that when $`I_B`$ has finally been reduced all the way to zero we find that the ladder may contain one or more stable vortex-antivortex pairs.

The occurrence of frequency locking, period doubling, quasi-periodic behaviour and symmetry breaking as switching is approached strongly suggests that the origin of the symmetry breaking is the occurrence of chaotic dynamics in the switching region. From the theory of chaotic dynamics it is known that two or more coupled Josephson junctions (or, equivalently, two or more coupled pendula) may exhibit three types of motion: oscillatory, rotational and chaotic. Chaotic motion arises for particular initial conditions. An analogous situation arises in non-linear coupled lattices which also display both breathers and chaotic dynamics. In a large array of coupled Josephson junctions (or, analagously, an array of coupled pendula) we must also expect these three types of motion. While roto-breathers are the most characteristic stable solutions it appears that the switching between breathers occurs only in the chaotic regime. Of course, the chaotic dynamics has an influence on the selection of breather type. Any small difference in the initial conditions leads to a completely different breather. Thus, in this chaotic regime, any small experimental imperfections or small inaccuracies in a computer simulation may lead to different chaotic trajectories and hence different roto-breathers. Symmetry breaking arises if the perturbation itself breaks symmetry.

We conclude that our model exhibits most of the main features  of the roto-breathers recently observed in Josephson ladder arrays, sheds light on the switching mechanism and predicts further observable behaviour. The occurrence of frequency locking, period doubling, quasi-periodic behaviour and symmetry breaking in the switching region suggests that switching occurs in the chaotic regime. We find that the maximum flux increases sharply in the switching region and that the flux clearly plays an important role in the switching mechanism. We also find two critical switching frequencies and a critical frequency below which no resistive behaviour is observed. All three frequencies are approximately independent of the dissipation constant $`\alpha `$ (for $`\alpha <0.15`$),

We are grateful to A.V. Ustinov and P. Binder for kindly giving us their data and explaining it prior to publication. We are also grateful for discussions with M.V. Fistul and S. Flach and for the hospitality of the Max-Planck-Institut, Dresden and the Universität Erlangen-Nürnberg.


text:
# Fast fluctuating fields as the source of low-frequency conductance fluctuations in many-electron systems and failure of quantum kinetics

## I Introduction

The electrical 1/f-noise (flicker noise) was discovered in 1925. Up to now low frequency 1/f-noise is known in different physical, chemical, biological, cosmic, geophysical and social systems \[1-5\], but there is no conventional explanation of this phenomenon . In electronics frequently just the 1/f conductance fluctuations limit qualities of real devices. What is important, 1/f-noise is much higher sensitive to mesoscopic structure of conductors as well as to external influences than conductivity itself, hence, it may bring a rich and delicate information on microscopic mechanisms of conductivity.

Traditionally 1/f-noise is thought be reducable to some ”slow” random processes with a broad distribution of large life-times (relaxation times). The standard idea is that conductivity determined by ”fast” mechanisms is modulated by ”slow” varying parameters (number of charge carriers, local disorder, occupation of electron traps, Coulomb potentials, lattice magnetic moments, etc., for some hypotheses see \[1-4,6-14,17-18\]). But concrete origin of long life-times remains misterious even in best experimentally investigated situations such as the electron mobility 1/f fluctuations in semiconductors  (besides, this theory could scarcely likely say why 1/f-noises in liquid and solid metals are comparable , or why hypothetical activation energies required to fit observed 1/f spectra much exceed actual energies determined from the response to charge injections ).

The problem looks sharply in case of ”bad” (narrow-band, variable-range hopping, tunnel, magnetically disturbed) conductivity like in doped semiconductors , oxides , amorphous silicon \[11-14\], materials with colossal magnetoresistance \[15-16,19,20,32\], cermets , etc. In these systems long-range Coulombian (or magnetic) interactions are important stimulating the attempts to reduce 1/f-noise to slow random charge redistributions \[4,11-14,17-18\], but it is hard to explain why real 1/f spectra do not saturate below 100 Hz . The interesting concept of self-organized criticality  also connects low-frequency noise with large spatial and temporal scales.

The alternative way to understanding 1/f-noise is developed since 1982 \[5,21-25,28-30\]. It attributes 1/f-noise not to ”slow” processes but, in opposite, to ”fast” microscopic kinetic processes responsible for resistivity. The logics is simple: if producing stochastic behaviour any dynamical system constantly forgets its own history, hence, it is indifferent to a number of kinetic events happened in the past, hence, it has neither stimulus nor possibility to follow some absolutely certain ”probability of event per unit time”. Therefore, a dynamical (Hamiltonian) system what produces relaxation (irreversibility) and noise makes this, generally speaking, without keeping definite (well time-averagable) ”probabilities per unit time” and consequently produces fluctuations of the relaxation (dissipation) and noise rates. These fluctuations do not destroy detailed balance and do not cause any compensating reaction, hence, have no characteristic time scale (see  for more details and examples). In this theory the long-living statistical correlations associated with 1/f spectrum are the manifestation of pure freedom of random flow of kinetic events, not of long memory \[21-25\]. In this respect, as N.S.Krylov argued  in 1950, we all should get rid of the erroneous opinion that statistical correlations always reflect some actual causality.

The kinetic theory misses this 1/f-noise just because assumes quite certain ”probabilities per unit time” although this ansatz is masked by ansatzes like ”molecular chaos”, ”random phases”, ”thermodynamical limit”, ”continuous spectra”, etc. But the correctly performed derivation of gas kinetics from Hamiltonian dynamics  demonstrates violation of ”molecular chaos” and 1/f fluctuations of diffusivities and mobilities of gas particles. Thus 1/f-noise arises even in the system where wittingly nothing like slow processes and giant life-times do exist.

While a gas is characterized by strong but short (well time-space separated) kinetic events (collisions), many systems possess weak but long-lasting (bad separated) interactions. For example, phononic systems (dielectric crystals). Nevertheless, as was shown in , 1/f type fluctuations of dissipation (inner friction) and of Raman light scattering observed in these systems (quartz) also can be deduced just from the basic (phonon-phonon) kinetics, if only take into account that any particular interaction interplays parametrically with other simultaneous interactions. Thus, it is wrong to suppose ”elementary” kinetic events be statistically independent, in contrary to what we see in kinetics.

The purpose of the present paper is to show that similar situation may realize in quantum conducting many-electron systems. Concretely, that the very evolution of quantum amplitude and probability concerning one particular electron transfer may depend on fast fluctuating potentials and fields (electric, magnetic, etc.) induced by other simultaneously occuring transfers, i.e. many elementary kinetic events (quantum jumps) are essentially interplaying. As far as we know, this effect still was not under consideration. It is expected be especially strong if characteristic transfer time (the time of evolution of quantum transfer probability right up to certainty) noticably exceeds correlation time of fluctuating fields. Such a condition seems natural if any conducting electron feels displacements of many other electrons (or changes of many magnetic moments) by means of long-range interactions.

Our first principal statement is that short-correlated random fields if influencing the formation of quantum transfer amplitudes lead to strong uncertainty of transfer probabilities, thus well defined ”probabilities per unit time” do not exist. As a consequence, long-correlated conductance fluctuations arise whose principal scale properties coinside with thats of 1/f-noise. Hence, to explain 1/f-noise we have no need in extremely slow charge redistributions, instead the very fast ones are sufficient. The second statement is that these effects can not be catched if neglect the actual quantum discreteness of energy spectra.

For simplicity, here to comprehend the main idea we concentrate on the case of tunnel conduction (some aspects of the discreteness in tunnel junctions were touched in  but only relating to low temperatures without accounting for Coulombian effects). More extended argumentation including rigorous analysis of Hamiltonian models will be published separately.

## II Characteristic times of tunnel conductivity

Let us consider characteristic time scales relating to electron tunneling between metallic sides. Under a small voltage $`U`$ applied to tunnel junction the mean charge transported during time $`\mathrm{\Delta }t`$ can be phenomenologically expressed as

$`\mathrm{\Delta }Q=e{\displaystyle \frac{Ue}{\delta E}}{\displaystyle \frac{\mathrm{\Delta }t}{\tau _{trans}}}`$

Here $`\delta E`$ is the mean separation of electron energy levels, $`Ue/\delta E`$ is the number of levels effectivelly contributing to electric current, and $`\tau _{trans}`$ is the mean transmission time required for one electron jump from a given level at one side to wherever at opposite side. Though quantum jumps can realize in a moment, that moment is random and may lie approximately equally anywhere in the interval $`0<\mathrm{\Delta }t`$ $`<\tau _{trans}`$ , with $`\tau _{trans}`$ being the time necessary to accumulate the quantum transmission probability up to a value $`1`$ . Clearly, $`e/\tau _{trans}`$ serves as the mean current per level. Here from one gets the tunnel conductance as

$$G=\mathrm{\Delta }Q/U\mathrm{\Delta }t=e^2\nu \gamma ,\nu =1/\delta E,\gamma =1/\tau _{trans}$$
(1)
with $`\nu `$ being electron density of states and $`\gamma `$ mean jump probability per unit time.

Of course, any real junction possesses a finite capacity $`C`$ and thus finite characteristic time

$`\tau _{rel}RCC/G`$

Its physical role can be identified (as usually in RC-circuits) as the relaxation time of junction charging and the correlation time of thermal voltage fluctuations (if only suppose that Coulombian interaction between sides manifests itself in stochastic form). Compare the above defined times:

$$\tau _{trans}/\tau _{rel}e^2\nu /C=E_C/\delta E$$
(2)
with $`E_C=e^2/C`$ being the Coulomb energy. Clearly, (2) is merely number of levels promoting charge relaxation.

Now let us demonstrate that this ratio may essentially exceed unit,

$$\tau _{trans}/\tau _{rel}>>1$$
(3)
even if Coulombian effects are weak in the trivial sense $`E_C<<T`$ . For certainty, consider a flat junction with area $`S`$ , side thicknesses $`w`$ , and with $`d`$ and $`ϵ`$ being the thickness and dielectricity of isolating barrier, respectively, and use formula $`CϵS/4\pi d`$ . The estimate for the metallic density of states is $`\nu \mu /N_e`$ $`Sw/(\mathrm{}v_Fa^2)`$ , where $`\mu `$ is Fermi energy, $`N_e`$ number of metallic electrons, $`v_F`$ Fermi velocity and $`a`$ is atomic size ( $`Sw/N_e=a^3`$ is volume per one electron), thus

$`e^2\nu /C{\displaystyle \frac{e^2}{\mathrm{}c}}{\displaystyle \frac{c}{v_F}}{\displaystyle \frac{4\pi dw}{ϵa^2}}{\displaystyle \frac{dw}{a^2}}`$

where $`c`$ is the speed of light (introduced to see that $`c/v_F`$ overcompensates the fine structure constant) and typical value $`ϵ20`$ is substituted.

Obviously, the inequality (3) is satisfied if either $`d`$ or $`w`$ (moreover if both) a few times exceeds $`a`$ , i.e. practically always (in this sense Coulombian effects never could be neglected). Therefore the quantum probability of some particular electron jump grows inevitably under influence of relatively fast varying inter-side voltage $`u(t)`$ , $`u(t)\sqrt{T/C}`$ , produced by many other jumps at the same time factually happening in both directions between other levels. Not only a moment of jump realization is random but the jump probability itself turns be random, that is different kinetic events become entangled. This may be named quantum (Coulombian) interaction of electron transfers. It can not be completely described by one-electron language, but the insert of fluctuating voltage gives the quasi one-particle approximation which, to some extent, may substitute for a real picture.

## III Time ratios connected with energy discreteness

The standard kinetic scheme deals with the tunnel Hamiltonian

$`H=H_0+H_{tun},H_{tun}={\displaystyle \underset{kq}{}}g_{kq}(b_q^+a_k+a_k^+b_q)`$

and with three ansatzes attracted to avoid formal difficulties brought in by the discreteness of electron energy levels, namely: i) energy spectrum in sides is so dense that the continuous limit is possible, $`_k\mathrm{}`$ $``$ $`\mathrm{}\nu (E)𝑑E`$ ; ii) the Fermi’s golden rule $`p_{kq}(\mathrm{\Delta }t)`$ $``$ $`2\pi \mathrm{\Delta }t(g_{kq}^2/\mathrm{})`$ $`\delta (E_{kq})`$ , where $`p_{kq}(t)`$ is jump probability and $`E_{kq}`$ energy difference between states; iii) $`g_{kq}^2`$ is sufficiently smooth function of $`E_{kq}`$ .

Here $`\mathrm{\Delta }t`$ is a time interval necessary to adequately evaluate jump probabilities for kinetic equations. This scheme requires the restriction $`\mathrm{\Delta }t<<\tau _{gold}`$ , where $`\tau _{gold}=2\pi \mathrm{}/\delta E`$ . However, if we wanted to account for effects of the voltage fluctuations $`u(t)`$ we would need at least in $`\mathrm{\Delta }t`$ comparable with $`\tau _{trans}`$ , i.e. in the additional condition $`\tau _{trans}/\tau _{gold}<1`$ . But in an adequate model just the opposite relation must be expected,

$$\tau _{trans}/\tau _{gold}>1$$
(4)
It is easy seen if note that $`\tau _{trans}/\tau _{gold}`$ $`=\delta E/\mathrm{\Delta }E`$ , $`\mathrm{\Delta }E2\pi \mathrm{}/\tau _{trans}`$ , where $`\mathrm{\Delta }E`$ is the energy uncertainty associated with instability of intra-side electron states. Hence, if the desirable condition was valid, then close states would be undistinguishable, in other words, electron spectra in sides would undergoe mutual rebuilding because of too good transparency of tunnel barrier (see also  on the relation $`R/R_0>1`$ ).

To analyse the ratio (4) more carefully, estimate $`\tau _{trans}`$ . The mean transported charge $`\mathrm{\Delta }Q`$ is expressed by

$`\mathrm{\Delta }Q=e(\mathrm{\Delta }N_+\mathrm{\Delta }N_{})=e{\displaystyle \underset{kq}{}}[f(E_k^{})f(E_q^+)]p_{kq},\mathrm{\Delta }N_\pm ={\displaystyle \underset{kq}{}}f(E_k^{})[1f(E_q^\pm )]p_{kq}`$

where $`\mathrm{\Delta }N_\pm `$ is number of electrons tunneling in left (right) direction, $`E^\pm `$ are energy levels in the sides, $`f(E)=1/[1+\mathrm{exp}(E\mu )/T]`$ , and

$$p_{kq}=p_{kq}(\mathrm{\Delta }t,U)4g_{kq}^2\mathrm{sin}^2(E_{kq}\mathrm{\Delta }t/2\mathrm{})/E_{kq}^2,E_{kq}=E_q^+E_k^{}eU$$
(5)
is tunneling probability evaluated by ordinary perturbation theory. The corresponding low field conductance

$`G=[\mathrm{\Delta }Q/U\mathrm{\Delta }t]_{U0}=e^2{\displaystyle \underset{kq}{}}[f(E_k)/E_k]p_{kq}(\mathrm{\Delta }t,0)/\mathrm{\Delta }t`$

turns into Eq.1 with

$$\gamma =1/\tau _{trans}=\underset{q}{}p_{kq}(\mathrm{\Delta }t,0)/\mathrm{\Delta }t2\pi g^2\nu /\mathrm{},(E_k,E_q\mu )$$
(6)
and $`g`$ being characteristic magnitude of $`g_{kq}`$ . Combining this relation with Eq.1 one obtains

$$\frac{\tau _{trans}}{\tau _{gold}}=\frac{e^2}{2\pi \mathrm{}G}=\frac{R}{R_0}\left(\frac{\delta E}{2\pi g}\right)^2,R_0=\frac{2\pi \mathrm{}}{e^2}20\text{kOhm}$$
(7)
Hence, the condition for really weak interaction is the smallness of transfer matrix elements as compared with the energy level spacing.

## IV Discreteness, phase decoherence and fluctuations of quantum transfer probabilities

Since the wish $`\tau _{trans}/\tau _{gold}<1`$ was invoked by the continuous spectrum approximation, we may suspect that the discreteness must be involved in an evident form to describe the fluctuations of transfer probabilities, while the perturbation theory is still applicable. Consider the picture when quantum transfers from a fixed level ”k” at one side to any level ”q” at another side are influenced by a fast fluctuating field (FFF), here the voltage noise $`u(t)`$ , in its turn induced by electron jumps beyond our attention. Thus we use quasi one-electron picture, basing on the ansatz that in many-electron surroundings with a sufficiently rich energy spectrum $`u(t)`$ behaves as a random process. The modern theory of quantum chaos gives powerful support for this statement  (although taking in mind that in specific systems coherent collective charge oscillations are possible instead of stochasticity and relaxation).

In this section we may consider thermodynamical equilibrium taking $`U=0`$ . First introduce the randomly accumulating diffusive phase shift

$`\phi (t)={\displaystyle \frac{e}{\mathrm{}}}{\displaystyle _0^t}u(t^{})𝑑t^{}`$

Instead of (5), the standard Shroedinger equations for transfer amplitudes yield

$$p_{kq}\left|A_{kq}\right|^2,A_{kq}\frac{g_{kq}}{\mathrm{}}_0^{\mathrm{\Delta }t}\mathrm{exp}(iE_{kq}t/\mathrm{})Z(t)𝑑t,Z(t)=\mathrm{exp}[i\phi (t)]$$
(8)
Clearly, now $`p_{kq}`$ become random values governed (parametrically excited or damped) by the phase shift in its turn produced by fluctuating potential difference between ”in” and ”out” states. Further, introduce the phase correlation function, the phase correlation time and corresponding energetical measure of quantum coherence by

$$K(t_1t_2)=Z(t_1)Z^{}(t_2),\tau _{phase}=_0^{\mathrm{}}\left|K(\tau )\right|𝑑\tau ,\epsilon _{coh}=2\pi \mathrm{}/\tau _{phase}$$
(9)
where angle brackets denote averaging with respect to FFF.

As demonstrated below, expectedly $`\tau _{phase}`$ is rather short as compared with $`\tau _{trans}`$ and naturally limited by $`\tau _{rel}`$ . Therefore under the integral in (8) $`Z(t)`$ can be treated as complex shot noise. Consequently, at $`\mathrm{\Delta }t>>\tau _{rel}`$ the transfer amplitudes $`A_{kq}`$ behave approximately as complex Brownian paths, and regardless of details of $`Z(t)`$’s statistics we have reasons to write

$$p_{kq}^2=\left|A_{kq}\right|^42\left|A_{kq}\right|^2^2=2p_{kq}^2,p_{kq},p_{kq}p_{kq}^2$$
(10)
Here the Malakhov’s cumulant brackets are used,

$`x,yxyxy`$

This is our first finding: if considered at time intervals of order of the actual transition time the quantum amplitudes and probabilities may become 100% uncertain due to the phase decoherence (phase diffusion) forced by FFF. The second is that the mean probabilities grow linearly with time:

$`p_{kq}\mathrm{\Delta }t(g/\mathrm{})^2{\displaystyle K(\tau )\mathrm{exp}(iE_{kq}\tau /\mathrm{})𝑑\tau }\mathrm{\Delta }t,(E_{kq}=E_q^+E_k^{})`$

i.e. FFF remove a need in the golden rule to have uniform growth of probabilities even at $`\mathrm{\Delta }t>>\tau _{gold}`$ . Notice for the next that according to this formula now $`\epsilon _{coh}`$ (instead of $`2\pi \mathrm{}/\tau _{gold}`$ ) serves as the width of energy region available from a given level (”k” or ”q”).

But most interesting subject is the summary transition rate $`\gamma =p/\mathrm{\Delta }t`$ , with $`p=p_k`$ ,

$$p_k\underset{q}{}p_{kq}=_0^{\mathrm{\Delta }t}\mathrm{\Gamma }(t_1t_2)Z(t_1)Z^{}(t_2)𝑑t_1𝑑t_2,\mathrm{\Gamma }(\tau )=\underset{q}{}(\frac{g_{kq}}{\mathrm{}})^2\mathrm{exp}(i\tau E_{kq}/\mathrm{})$$
(11)
Here the kernel $`\mathrm{\Gamma }(\tau )`$ represents the discreteness. Its analytical properties are of much importance for all the theory . In the continuous approximation provided ansatzes i)-iii) one would have $`\mathrm{\Gamma }(\tau )=`$ $`\delta (\tau )/\tau _{trans}`$ , and $`p`$ would be definitely nonrandom, $`p=\mathrm{\Delta }t/\tau _{trans}`$ . In reality, because of the discreteness this kernel is quite non-local and does not decay completely at arbitrary long time. For visuality only, let us choose equidistant spectrum at the right-hand side, $`E_q^+E_k^{}=`$ $`n\delta E+\epsilon `$ ( $`n`$ is integer, $`\epsilon \delta E`$ ), then

$$\mathrm{\Gamma }(\tau )=\frac{1}{\tau _{trans}}\mathrm{exp}(i\epsilon \tau /\mathrm{})\underset{n}{}\delta (\tau n\tau _{gold})$$
(12)
Now the third important point is easy seen : if $`\epsilon _{coh}>\delta E`$ then the mean probability $`p`$ $`\mathrm{\Delta }t/\tau _{trans}`$ practically coinsides with what is given by usual kinetics, even in spite of formal violation of the golden rule (in this sense, FFF effectively expand applicability of usual scheme).

And the main fourth point is that the phase decoherence produced by FFF if combined with the discreteness results in randomness of the summary quantum jump probabilities. Indeed, due to the possibility to consider $`A_{kq}`$’s as Brownian walks we obtain

$$p,p\mathrm{}_0^{\mathrm{\Delta }t}\mathrm{\Gamma }(t_1t_2)\mathrm{\Gamma }(t_3t_4)K(t_1t_4)K(t_3t_2)𝑑t_1\mathrm{}𝑑t_4$$
(13)
Here under the same condition $`\epsilon _{coh}>\delta E`$ only the regions $`t_1t_4,`$ $`t_3t_2`$ are significant but many delta functions from (12) contribute. The resulting transfer probability variance is

$$p,p\frac{\mathrm{\Delta }t^2}{\tau _{trans}^2\tau _{gold}}\left|K(\tau )\right|^2𝑑\tau \frac{\tau _{phase}}{\tau _{gold}}p^2$$
(14)
(we took into account that ”width” of delta functions determined by the inverse width of whole energy band is wittingly smaller than $`\tau _{phase}`$ ).

If inequality (4) was inverted because of too small $`\delta E`$ the expression (13) would formally turn into zero, but in fact it always remains nonzero under any realistic slightly non-equidistant spectrum. Nevertheless, we may predict that violation of (4), as well as of (3), leads to supression of the probability fluctuations. Oppositely, large $`\delta E>\epsilon _{coh}`$ means most strong fluctuations with variance $`p^2`$ or greater. But thats are accompanied by decreasing correlations between transfers from different levels:

$$p_{k_1,}p_{k_2}p_{k_1}p_{k_2}\frac{1}{\tau _{gold}}\mathrm{exp}[i(E_{k_1}^{}E_{k_2}^{})\tau /\mathrm{}]\left|K(\tau )\right|^2𝑑\tau $$
(15)
Evidently, while at $`\epsilon _{coh}>\delta E`$ the currents from many levels fluctuate concordanly, at $`\epsilon _{coh}<\delta E`$ even close levels inject independently one on another. Besides, even $`p`$ become sensible to the energy shift $`\epsilon `$ in (12), therefore, this extreme case should be carefully analysed with accounting for realistic (uncommensurable) discrete energy spectra (and better out of the frame of the quasi one-particle picture).

To end this section, estimate the phase correlation time $`\tau _{phase}`$ (what may be named also phase decoherence time). Notice that $`K(\tau )`$ is nothing but characteristic function of the phase. In principle, it is determined just by the transfer statistics, in a complicated self-consistent picture. Since the latter now is out of our look, we confine ourselves by rough reasonings. For instance, at $`E_C`$ $`<<`$ $`T`$ it is likely natural to treat $`u(t)`$ as the Ornstein-Uhlenbeck random process, then

$`K(\tau )\mathrm{exp}\left[{\displaystyle \frac{T}{C}}\left({\displaystyle \frac{e}{\mathrm{}}}\tau _{rel}\right)^2\{1{\displaystyle \frac{\tau }{\tau _{rel}}}\mathrm{exp}({\displaystyle \frac{\tau }{\tau _{rel}}})\}\right]`$

The corresponding decoherence time is $`\tau _{phase}`$ $`(\mathrm{}/e)\sqrt{C/T}`$ $`(<<\tau _{gold})`$ . Perhaply, this is low bound for it (too rigid in the sense that it does not include the conductance $`G`$ ). At $`E_CT`$ (what qualifies small junctions) the charge quantization is essential and the better model for $`\phi (t)`$ is infinitely divisible random walk (on this subject see ) formed with rare increments by $`\mathrm{\Delta }\phi =`$ $`(e/\mathrm{})(\pm e/C)\theta `$ , where $`\theta `$ is random duration of charged stay distributed with some probability density $`W(\theta )`$ and total time fraction $`1`$ . The suitable expression is $`W(\tau )=`$ $`\tau _{rel}^1\mathrm{exp}(\tau /\tau _{rel})`$ what corresponds to the characteristic function

$`\mathrm{\Xi }(\eta ,\tau )\mathrm{exp}[i\eta {\displaystyle _0^\tau }u(t)𝑑t]\mathrm{exp}\left\{{\displaystyle \frac{|\tau |}{\tau _{rel}}}{\displaystyle _0^{\mathrm{}}}[\mathrm{cos}\left({\displaystyle \frac{\eta e}{C}}\theta \right)1]W(\theta )𝑑\theta \right\}`$

(here the Levy-Khinchin representation  was applied). Taking $`\eta =e/\mathrm{}`$ and assuming $`2\pi R/R_0`$ $`>1`$ , obtain $`K(\tau )=`$ $`\mathrm{\Xi }(e/\mathrm{},\tau )`$ $`\mathrm{exp}(|\tau |/\tau _{rel})`$ . Thus, in this oppositely extreme case $`\tau _{phase}`$ may be of order of $`\tau _{rel}`$ (likely representing upper bound for $`\tau _{phase}`$ ) and possibly larger than $`\tau _{gold}`$ .

## V Low-frequency conductance fluctuations

In this section we return to externally driven junction, so the total voltage will be $`U(t)=U+u(t)`$ with the applied voltage $`U`$ interpreted as its average value.

Of course, fluctuations of electron jump probabilities eventually result in more or less analogous conductance fluctuations. The conductance is implied as $`\mathrm{\Delta }Q/U\mathrm{\Delta }t`$ where $`\mathrm{\Delta }Q`$ is the conditional quantum average value of transported charge taken under a fixed realization of $`u(t)`$ . The charge transport depends also on occupancies of energy levels. If $`U=0`$ then $`\mathrm{\Delta }Q`$ , by its definition, must turn into zero regardless of $`u(t)`$ . This means the existence of some natural statistical correlations between voltage noise and occupancies whose accurate description would need in a self-consistent many-electron picture. To avoid this difficulty, let us work up the storage

$$\mathrm{\Delta }Q=e\underset{kq}{}[f(E_k^{})f(E_q^+)]\left|A_{kq}\right|^2$$
(16)
making time integration by parts in the amplitudes as if $`u(t)`$ was absent, to transform (16) into

$$\mathrm{\Delta }Q=e^2U_0^{\mathrm{\Delta }t}\mathrm{\Lambda }(t_1t_2)\mathrm{exp}[ieU(t_1t_2)/\mathrm{}]Z(t_1)Z^{}(t_2)𝑑t_1𝑑t_2$$
(17)
where new kernel is introduced,

$$\mathrm{\Lambda }(\tau )=\underset{kq}{}(\frac{g_{kq}}{\mathrm{}})^2\frac{f(E_k^{})f(E_q^+)}{E_q^+E_k^{}}\mathrm{exp}[i\tau (E_q^+E_k^{})/\mathrm{}]$$
(18)
This form is consistent with the detailed balance $`\mathrm{\Delta }Q(U=0)=0`$ and hence may serve for estimates.

In general, if electron spectra in sides are much wider than $`T`$ then at $`U<T/e`$ we believe that, approximately, i) conductivity obeys Ohmic law and ii) its relative fluctuations do not depend on $`U`$ . Then Eqs. 17-18 can be simplified by means of linearization into

$$\mathrm{\Delta }Qe^2U\underset{k}{}[f(E_k^{})/E_k^{}]p_k$$
(19)
with $`p_k`$ defined by (11) . Naturally, the conductance fluctuations essentially depend on relation between decoherence and level spacing. At sufficiently small phase decoherence time (”large” junction, widely correlated jumps from different levels), the Eq.15 helps easy obtain

$$G,G\frac{\delta E}{T}G^2,(\tau _{phase}<<\tau _{gold})$$
(20)
At large decoherence time (”small” junction, non-correlated jumps) conductance fluctuations are very sensible to concrete structure of energy spectra in sides, first of all to the degree of their relative commensurability. Omitting the details, the result is that the relative conductance variance may prove to be anywhere in the interval

$$\frac{\tau _{phase}}{\tau _{gold}}\frac{\delta E}{T}<\frac{G,G}{G^2}<1,(\tau _{phase}>\tau _{gold})$$
(21)
i.e can achieve as great magnitude as $`1`$ . In this case, accurate accounting for fluctuations of the occupancies may be especially important.

## VI Comparison with experiment

The good experimental illustration for permissibly exciting properties of quantum transfers was brought in  where 1/f-noise in the cermet (granular composite) $`Ni`$-nanoparticles $`(25\%)`$-$`Al_2O_3`$was investigated. In this system the parameters of a typical elementary tunnel junction formed by neighbouring metal particles are $`\delta E0.2`$meV, $`d2`$nm, $`C510^6`$cm, $`E_CT`$ (at room temperature), and $`R30`$MOhm, which mean that $`\tau _{gold}310^{11}`$s , $`\tau _{rel}1.510^{10}`$s , $`\tau _{trans}/\tau _{rel}=E_C/\delta E200`$ and $`\tau _{trans}310^8`$s . Both the inequalities (3) and (4) are well satisfied, thus giving us all grounds (see also Discussion) to suspect that 1/f-noise could be attributed to quantum Coulombian interactions, in the above mentioned sense (possibly with a contribution by electron-phonon processes). Moreover, $`\tau _{gold}`$ is even noticably smaller than $`\tau _{rel}`$, thus indicating, from our point of view, the possibility of highly maintained conductance fluctuations.

In fact, such cermet is characterized by giant 1/f conductance noise with relative spectrum density $`S_{\delta G}(f)\alpha /N_gf`$ where $`\alpha 610^3`$ and $`N_g`$ is the number of metal particles in a sample. In view of $`E_CT`$, $`N_g`$ approximately represents the number of active (simultaneously transported) electrons , hence that is almost standard noise with classical Hooge constant ($`\alpha =210^3`$) . It corresponds to the $`S_{\delta G}(f)\alpha /f`$ noise in separate elementary junction, demonstrating at least $`100\%\sqrt{\alpha \mathrm{ln}(1\text{s}/\tau _{rel})}`$ , i.e. $`40\%`$, uncertainty of its conductance.

But most beautiful observation  was the visible 1/f-noise sensibility to the discreteness of electron energy spectra in metal granules. When applied voltage per elementary junction exceeds $`\delta E/e`$ and thus corresponding current exceeds $`e/\tau _{trans}`$ (see Eq.1), then quadratic (low-field) dependence of variance of voltage 1/f noise on bias current transforms into ”non-Ohmic” linear one (although mean conductance obeys Ohm law up to $`T/\delta E100`$ times larger current).

The experimental appearance of the factor $`\delta E/e`$ is evident manifestation of the role of discreteness and it itself gives clear confirmation of our aproach to 1/f-noise in this system. At the same time, this means that the hypothesis ii) (in previous section) fails, i.e. linearized expression (19) becomes invalid if applied to fluctuations, and the estimates (20)-(21) should be multiplied by a decreasing factor $`D(U)`$ ($`D(0)=1`$ , $`D(U)<1`$ ). Indeed, as we underlined, any transfer from one side embraces $`\epsilon _{coh}/\delta E`$ levels at another side. If this number is smaller than the number of charge transporting levels, $`eU/\delta E`$ , then the latters act as uncorrelated quantum channels. Hence, $`\mathrm{\Delta }Q`$ variance may become linear function of $`U`$ being proportional to the number of channels, what corresponds in small junctions to $`D\delta E/eU`$ at $`eU>\delta E`$ .

Definite formal support for this euristic reasoning (literally appropriate at zero temperature only) could be found already in the frame of quasi one-electron picture, if consider full non-linear dependence (17) of $`\mathrm{\Delta }Q`$ on $`U`$ . The decreasing factor what follows from the Eqs.17-18 is

$$D(U)=S(eU)/S(0),S(E)K_\mathrm{\Lambda }(\tau )\mathrm{exp}(iE\tau /\mathrm{})\left[K(\theta )K(\tau \theta )𝑑\theta \right]𝑑\tau $$
(22)
where the correlation function $`K_\mathrm{\Lambda }(\tau )`$ is defined by

$$K_\mathrm{\Lambda }(\tau )=\frac{1}{2\mathrm{\Delta }t}_{\mathrm{\Delta }t}^{\mathrm{\Delta }t}\mathrm{\Lambda }(t+\frac{\tau }{2})\mathrm{\Lambda }^{}(t\frac{\tau }{2})𝑑t$$
(23)
In case of small junctions the averaging over an ensemble of probable electron spectra must be added what automatically takes place in cermets. In theory, we need in some statistics of energy levels (see  on the known variants). At present, notice only that the absense of a rigid measure of local spacing of energy levels (except the mean value $`\delta E`$ as related to whole spectra only) might results naturally in a weak dependence of $`K_\mathrm{\Lambda }(\tau )`$ on $`\tau `$ characterized by logarithm $`\mathrm{ln}(\tau _{gold}/\tau )`$ which implies just the inverse proportionality $`D(U)\epsilon _{coh}/eU`$ at $`eU>\epsilon _{coh}`$ .

## VII Discussion and resume

The principal peculiarity of above considered fluctuations of quantum transfer probabilities and corresponding conductance fluctuations is that their relative measure seems independent on the time under observation, $`\mathrm{\Delta }t`$ . It looks as if conductance undergoe fluctuations with non-decaying statistical correlations. Outwardly thats resemble static fluctuations investigated in the theory of disordered conductors (so-called universal conductance fluctuations) . But by essence these are different things: one starts from the decoherence just when other finishes at it.

Though our analysis was limited by times $`\mathrm{\Delta }t\tau _{trans}`$ , there is a feeling that similar statistics expands to arbitrarily longer time scales. Formal proof that it is really true, attributing to Hamiltonian models of quantum channels (with ”weak links” like tunnel barriers) subjected to FFF (Coulombian or magnetical), will be done separately. For the present, semi-formal arguments are: i) all what was obtained results from the trivial rule that in general one should manipulate with quantum amplitudes to find quantum probabilities as a final product; ii) all what was obtained results from ”fast” noise and phase decoherence, but not from some causal correlations which could not be continued to longer time. We believe that a self-consistent analysis connecting FFF and electron transport statistics will result in some slow (nearly logarithmic) dependence of $`G,G`$ on the observation time reflecting non-static character of the fluctuations and described by 1/f spectrum, instead of $`\delta (f)`$ having the same formal dimensionality (such dependence may begin already as a weak $`\mathrm{\Delta }t`$ dependence of the correlator (23)).

An accurate approach to a tunnel junction would deal with non-equilibrium steady state governed by the Hamiltonian

$`H=H_0+H_{tun},H_0=H_{}+H_++H_CU\mathrm{\Delta }Q`$

where $`H_C`$ describes inter-side Coulombian interaction, $`\mathrm{\Delta }Q`$ is the operator of transported charge, and $`H_\pm `$ describe two sides with their leads serving also as thermal baths ensuring relaxation to equilibrium occupancies. Like in practice, in theory it is rather hard to ”solder leads”, especially if wishing to avoid appeals to standard kinetic schemes. But, regardless of tehnical difficulties, we may state that none processes in leads could influence the formation of non-diagonal (inter-side) elements of the operator $`\mathrm{\Delta }Q`$ since thats are determined by the tunneling itself and charging of tunnel barrier only. Therefore, the better relaxation in leads the more grounds we have to extend the Eq.16 to arbitrary time intervals (now with $`\left|A_{kq}\right|^2`$ representing the number of passing electrons). The comparatively non-principal corrections to be performed are accounting for fluctuations of occupancies and including conductance fluctuations into higher-order statistics of FFF (voltage noise).

To resume, we demonstrated that if do not neglect the actual quantum discreteness when constructing kinetic models of transport processes then possible strong sensibility of quantum transfer amplitudes and probabilities to fast fluctuating fields becomes visible (in particular, created by the transfers themselves) which may result in fundamental 1/f type low-frequency fluctuations of transport rates. Hence, the now reigning quantum kinetic models ask for evident general comments. All come from the well known Pauli’s kinetic master equations. In  Van-Hove developed its formal groundation under so called $`\lambda ^2t`$ -limit. But, with no doubts, this theory does not foresee anything like 1/f noise. What is the matter, the question arises naturally (some comments on this issue were suggested in ). We hope that our above consideration highlights the possible principal answer: in fact the Van-Hove’s formalism supposes the limit $`\delta E0`$ (i.e. the continuous spectrum idealization) be performed before the limit $`g0`$ with $`g`$ ($`g=\lambda `$) representing (as above) the magnitude of weak interactions. Thus the influence of quantum discreteness onto statistics of quantum jumps (between eigenstates of unperturbed Hamiltonian) is lost. We think that it will be not hopeless to properly improve the present kinetics (see Introduction).

ACKNOWLEDGEMENTS

I acknowledge Dr. Yu.V.Medvedev and participants of his seminar at the Department of kinetic properties of disordered and nonlinear systems in DonPTI NAS of Ukraine for support, stimulating criticism and helpfull discussions.

REFERENCES

1. P.Dutta, P.Horn, Rev.Mod.Phys., 53(3), 497 (1981).

2. F.N.Hooge, T.G.M.Kleinpenning, L.K.J.Vandamme, Rep.Prog.Phys., 44, 481 (1981).

3. M.B.Weissman, Rev.Mod.Phys., 60, 537 (1988).

4. M.B.Weissman, Rev.Mod.Phys., 65, 829 (1993).

5. G.N.Bochkov and Yu.E.Kuzovlev. UFN, 141, 151 (1983), English translation: Sov.Phys.-Usp., 26, 829 (1983).

6. B.Raquet, J.M.D.Coey, S.Wirth and S. Von Molnar, Phys.Rev., B59, 12435 (1999).

7. J.G.Massey and M.Lee, Phys.Rev.Lett., 79, 3986 (1997).

8. M.J.C. van den Homberg, A.H.Verbruggen, P.F.A.Alkemade, S.Radelaar, E.Ochs, K.Armbruster-Dagge, A.Seeger and H.Stoll, Phys.Rev., B 57, 53 (1998).

9. A.Ghosh, A.K.Raychaudhuri, R.Streekala, M.Rajeswari and T.Venkatesan, Phys.Rev., B 58, R14666 (1998).

10 X.Y.Chen, P.M.Koenrad, F.N.Hooge, J.H.Wolter and V.Aninkevicius, Phys.Rev., B 55, 5290 (1997).

11. G.M.Khera and J.Kakalios, Phys.Rev., B 56, 1918 (1997).

12. M.Gunes, R.E.Johanson and S.O.Kasap, Phys.Rev., B 60, 1477 (1999).

13. K.M.Abkemeier, Phys.Rev., B 55, 7005 (1997).

14. G.Snyder, M.B.Weisman and H.T.Hardner, Phys.Rev., B 56, 9205 (1997).

15. A.Lisauskas, S.I.Khartsev and A.M.Grishin, Studies of 1/f-noise in $`La_{1x}M_xMnO_3`$ (M=Sr,Pb) epitaxial thin films, in J.Low.Temp.Phys. as MOS-99 Proceedings.

16. A.Lisauskas, S.I.Khartsev, A.M.Grishin and V.Palenskis, Electrical noise in ultra thin giant magnetoresistors, in Mat.Res.Soc.Proc. Spring-99 Meeting.

17 B.I.Shklovskii and A.L.Efros. Electronic properties of doped semiconductors, Springer-Verlag. Berlin, 1984.

18. V.I.Kozub, Solid State Commun., 97, 843 (1996).

19. V.Podzorov, M.Uehara, M.E.Gershenson and S.-W.Cheong, lanl arXiv cond-mat/9912064.

20. M.Viret, L.Ranno and J.M.D.Coey, Phys.Rev., B 55, 8067 (1997).

21. Yu.E.Kuzovlev and G.N.Bochkov. On the nature and statistics of 1/f-noise. Preprint No.157, NIRFI, Gorkii, USSR, 1982.

22. Yu.E.Kuzovlev and G.N.Bochkov.Izv.VUZov.-Radiofizika, 26, 310 (1983), transl. in Radiophysics and Quantum Electronics (RPQEAC, USA), No 3 (1983).

23. G.N.Bochkov and Yu.E.Kuzovlev. Izv.VUZov.-Radiofizika, 27, 1151 (1984), transl. in Radiophysics and Quantum Electronics (RPQEAC, USA), No 9 (1984).

24. G.N.Bochkov and Yu.E.Kuzovlev. On the theory of 1/f-noise. Preprint N 195, NIRFI, Gorkii, USSR, 1985.

25. Yu.E.Kuzovlev, lanl arXiv cond-mat/9903350 .

26. N.S.Krylov. Works on the foundations of statistical mechanics. Princeton U.P., Princeton, 1979 (Russian original book published in 1950).

27. Yu.A.Genenko and Yu.M.Ivanchenko, Teor.Mat.Fiz., 69, 142 (1986).

28. Yu.E.Kuzovlev, Zh. Eksp. Teor. Fiz, 94, No.12, 140 (1988), transl. in Sov.Phys.-JETP, 67 (12), 2469 (1988).

29. Yu.E.Kuzovlev, Phys.Lett., A 194, 285 (1994).

30. Yu.E.Kuzovlev, Zh. Eksp. Teor. Fiz, 111, No.6, 2086 (1997), transl. in JETP, 84(6), 1138 (1997).

31. P.Bak. Self-organized criticality: why nature is complex. Springer, N-Y, 1996.

32. J.M.D.Coey, M.Viret and S. von Molnar, Adv.Phys., 48, 167 (1999).

33. G.Casati and B.Chirikov. Fluctuations in quantum chaos. Preprint, Budker Inst. of Nuclear Physics SB RAS, 1993.

34. C.W.J.Beenakker, Rev.Mod.Phys., 69, N 3, 731 (1997).

35. L.Van Hove, Physica, 21, 517 (1955).

36. J.V.Mantese, W.I.Goldburg, D.H.Darling, H.G.Craighead, U.J.Gibson, R.A.Buhrman and W.W.Webb, Solid State Commun., 37, 353 (1981).

37. W.Feller, Introduction to probability theory and its applications, John Wiley, N-Y, 1966.

38. Yu.A.Genenko and Yu.M.Ivanchenko, Phys.Lett., 126, 201 (1987).


text:
# Small-World Rouse Networks as models of cross-linked polymers

## I Introduction

Topological properties of polymers can dramatically effect their dynamical properties, such as their collapse in bad solvents and their response to external forces . Such forces can be applied microscopically, either by having charged polymers (polyelectrolytes, polyampholytes) in electrical fields, or via optical tweezers or magnetic beads. In this communication we study the stretching of cross-linked objects, whose backbones are regular lattices (we will consider for simplicity a ring), a few elements of the backbone being chemically connected to each other via cross-links. An experimental realization may be a very dilute solution of linear chains which are then cross-linked by irradiation or chemically. Our model is a realization of the so-called small-world networks (SWN), and the disorder (the cross-links) is in statistical terms quenched. Considering all bonds to be equal we study the dynamics in the framework of the Rouse model, in which the monomers are connected by harmonic springs; we term this structure the small-world Rouse network (SWRN). Our SWRN is built out of an N-monomer ring with superimposed fixed links between randomly chosen monomers; no additional links are generated or broken in external fields. The SWRN is a new intermediate between linear chains and networks. Distinct from Cayley-trees, which model hyperbranched polymers without loops, loops are a fundamental ingredient here. As such the SWRN is interesting in its own right as a study of the interplay between dynamics and topology, and it belongs to the class of generalized Gaussian structures.

## II Dynamics of small-world Rouse networks

The construction of the SWN which we consider here (see also ), is slightly different from the original one by Watts and Strogatz, but preserves its main SWN-characteristics. Starting from a ring of $`N`$ monomers, we cross-link with probability $`p`$ each monomer randomly to any of the monomers of the network. Such cross-links thus connect monomers far apart along the chemical backbone (the ring), rendering them close in Eucledian space. The SWRN-monomers are, in accordance with the Rouse model, connected by harmonic springs of strength $`k`$. The position $`𝐑_n(t)`$ of the $`n`$th bead under the action of an external force $`𝐅_n(t)`$ and in the presence of thermal noise $`𝜼_n(t)`$ is governed by the Langevin equation

$$\gamma \frac{d𝐑_n(t)}{dt}=k\underset{j=1}{\overset{N}{}}A_{nj}𝐑_j(t)+𝐅_n(t)+𝜼_n(t).$$
(1)
Here $`\gamma `$ is the coefficient of friction, and the matrix $`A_{ij}`$ is the connectivity matrix of the network. It is defined as follows: Every connection between site $`i`$ and site $`j`$ contributes $`1`$ to $`A_{ij}`$, while $`A_{ii}`$ and $`A_{jj}`$ are determined from the condition that $`_jA_{ij}=_iA_{ij}=0`$. Monasson  recently published a detailed study of the spectrum of the small-world connectivity-matrix $`A`$. Among his findings was the existence of a “pseudo-gap” in the SWN density of states $`\rho (E)`$ which has the form

$$\rho (E)E^{1/2}\mathrm{exp}\left(\frac{C}{\sqrt{E}}\right),E0$$
(2)
This behavior could not be confirmed numerically through direct diagonlization, and appeared in the data as a real gap. There are, as we shall see, other ways of probing this behavior numerically, and our results support that $`\rho (E)`$ behaves as Eq. (2).

Writing

$$𝐑(t)(𝐑_1(t),𝐑_2(t),\mathrm{},𝐑_N(t))^T$$
(3)
and similary for the other quantities in Eq. (1), we can rewrite the equation of motion in a more condensed form as

$$\frac{d𝐑(t)}{dt}=\sigma A𝐑(t)+\frac{𝐅(t)}{\gamma }+𝐰(t).$$
(4)
Here we have introduced $`\sigma k/\gamma `$ and $`𝐰𝜼/\gamma `$. In the case of a spatially constant external force, the solution to Eq. (4) is obtained as

$$𝐑(t)=_{\mathrm{}}^t𝑑se^{\sigma A(ts)}\left(\frac{𝐅(s)}{\gamma }+𝐰(s)\right)$$
(5)
We now specialize to the following situation: the force is switched on at time $`t=0`$ and pulls only the $`m`$th bead in the $`y`$-direction, i.e. $`𝐅_i(t)=\theta (t)\delta _{i,m}F\widehat{𝐲}`$. Here $`\theta (t)`$ is the Heaviside step-function, $`\delta _{i,j}`$ Kronecker’s delta and $`\widehat{𝐲}`$ is a unit-vector pointing in the $`y`$-direction. We focus on the displacement of the $`m`$th bead along the $`y`$-axis, and average over thermal noise, using $`𝐰(t)=0`$. Finally we perform a structural average over $`m`$ and end up with (see e.g.  and references therin for details)

$$Y(t)\frac{1}{N}\underset{m}{}𝐑_{m,y}(t)=\frac{F}{N\gamma }_0^t𝑑s\underset{i}{}e^{\sigma \lambda _is}=\frac{Ft}{N\gamma }+\frac{F}{N\gamma \sigma }\underset{i=2}{\overset{N}{}}\frac{1e^{\sigma \lambda _it}}{\lambda _i}.$$
(6)
In this equation $`𝐑_{m,y}`$ is the $`y`$ component of $`𝐑_m`$ and $`\lambda _i`$ with $`i=1\mathrm{}N`$ are the eigenvalues of the connectivity matrix $`A`$. The last expression follows from the fact that for a connected structure, only one eigenvalue vanishes (say, $`\lambda _1`$). At times $`t`$ much smaller than the time scale set by the largest eigenvalue $`\lambda _{max}`$, i.e. when $`\sigma \lambda _{max}t1`$, $`Y(t)`$ increases linearly in time: $`Y(t)Ft/\gamma `$. That is, only the monomer being pulled moves with a constant speed, not yet feeling the influence of the other monomers. Likewise, at late times $`\sigma \lambda _{min}t1`$, where $`\lambda _{min}`$ is the lowest non-vanishing eigenvalue, the entire polymer is being pulled with a constant speed, $`Y(t)Ft/(N\gamma )`$. These observations are independent of the specific structure being pulled, and only in the intermediary regime $`\lambda _{max}^1\sigma t\lambda _{min}^1`$ does the particular topology of the polymer affect the dynamics, namely through the spectrum of the connectivity matrix.

The numerical computation of the quantity $`Y(t)`$ above proceeds as follows. From a specific realization of an $`N=1000`$ small-world network, we construct the corresponding connectivity matrix. Then we find the $`N`$ eigenvalues using standard routines, and implement Eq. (6). To get an idea of the importance of sample to sample fluctuations, we consider first $`10`$ different realisations of the SWRN for $`p=0.05`$. Plotted is in Fig. 1 on double logarithmic scales $`Y`$ as a function of $`t`$, where here and in the following we use the dimensionless variables $`Y^{}(t)\sigma \gamma Y(t)/F`$ and $`t^{}\sigma t`$. In Fig. 1 we display the envelope of all $`10`$ realizations, i.e. the two curves are the extremal two “worst” cases. We see that the difference in the results is quite small (and appears, as it should, only for intermediary times), and therefore we regard results from any specific realization as being typical. In Fig. 2 we analyze the dependence of $`Y(t)`$ on cross-linking, by varying $`p`$ from $`p=0`$ (standard Rouse-model of the ring) to $`p=0.01`$ and $`p=0.05`$. We note first that the differences are now considerably larger than in Fig. 1. Second, for $`p=0`$ i.e. for the Rouse chain, we have the standard picture: a subdiffusive $`\sqrt{t}`$ behavior at intermediary times is followed by a diffusive $`t`$ behavior at longer times. At very early times we also have a linear behavior, albeit not visible in the range of the figure. The initial and final dynamics are in accordance with the explanation given above. The intermediate behavior reflects the structure of the spectrum of $`A`$, and is also well understood in the Rouse case; it is a result of the rather slow propagation of disturbances through the chain (here the ring).

We infer from the other curves in Fig. 2, that even very small but nonvanishing $`p`$ affect the intermediate behavior of $`Y(t)`$ quite strongly. For increasing $`p`$ the curves bend downwards from the $`p=0`$ case, mirroring the increased stiffness of the polymer due to the additional links. Thus a ring with cross-links can be easily distinguished experimentally (say through NMR or electronic energy transfer) from one without cross-links, whose $`Y(t)`$ dynamics under $`𝐅`$ is never slower than $`\sqrt{t}`$. As expected, the very early and very late behavior in all three cases coincide, being independent of the specific structure under scrutiny.

In Fig. 3 we plot on logarithmic scales $`Y^{}(t)`$ as a function of time $`t^{}`$ for several values of $`p`$. Increasing $`p`$ increases the stiffness of the polymer, and this is reflected in the intermediary regime, which becomes almost flat for large $`p`$. Moreover the long time behavior of the polymer is reached much earlier for polymers with a higher number of cross-links. In line with the discussion above of the range of the intermediary regime ($`\lambda _{max}^1\sigma t\lambda _{min}^1`$), this feature means that the lowest non-vanishing eigenvalue $`\lambda _{min}`$ gets quite large, and it is hence related to the appearance of a (pseudo) gap in the spectrum of $`A`$.

## III Theoretical Analysis

As indicated earlier, the initial as well as the asymptotic behavior of $`Y(t)`$ are well understood. Thus we will concentrate here on the richer and much more complex intermediate behavior. We shall rewrite Eq. (6) using a continuous picture, based on the density of states $`\rho (\lambda )=lim_N\mathrm{}(1/N)_i\delta (\lambda \lambda _i)`$, but continue to separate out the vanishing eigenvalue from the rest. Hence, with $`ϵ`$ very small, $`ϵ0^+`$:

$$Y(t)=\frac{F}{N\gamma }_0^t𝑑s\underset{i}{}e^{\sigma \lambda _is}=\frac{Ft}{N\gamma }+\frac{F}{\gamma }_0^t𝑑s_ϵ^{\mathrm{}}𝑑\lambda \rho (\lambda )e^{\sigma \lambda s},$$
(7)
an expression which is a fortiori correct in the presence of a gap, where one can take $`0<ϵ\lambda _{min}`$. It will be convenient also to consider the stretching (relative motion) $`\mathrm{\Delta }(t)`$ separately:

$$\mathrm{\Delta }(t)Y(t)\frac{Ft}{N\gamma }=\frac{F}{\gamma }_0^t𝑑s_{0^+}^{\mathrm{}}𝑑\lambda \rho (\lambda )e^{\sigma \lambda s}.$$
(8)
We remark that the inner integral in Eq. (8) is related to the probability for a random walker to be present at the original site. We therefore first analyse the behavior of this quantity:

$$P_0(t)_0^{\mathrm{}}𝑑\lambda \rho (\lambda )e^{\lambda t}$$
(9)
The asymptotic temporal behavior is accessed through the behavior of $`\rho (\lambda )`$ at small $`\lambda `$. Inserting the expression of Monasson Eq. (2) into Eq. (9) we obtain:

$$P_0(t)_0^{\mathrm{}}𝑑\lambda \lambda ^{1/2}\mathrm{exp}\left(\frac{C}{\sqrt{\lambda }}\lambda t\right)=t^{2/3}\frac{d}{dC}_0^{\mathrm{}}𝑑y\mathrm{exp}\left(t^{1/3}(\frac{C}{\sqrt{y}}y)\right)$$
(10)
The asymptotic behavior of the integral follows readily from a saddle-point procedure , so that we end up with

$$P_0(t)t^{1/2}\mathrm{exp}(C^{}t^{1/3}),$$
(11)
whith $`C^{}=3\left(C/2\right)^{2/3}`$.This is by itself a quite interesting and novel result, and it compares favourably to our numerical simulations for $`P_0(t)`$. The result is close in form to that for Cayley trees, where $`P_0(t)t^{3/2}\mathrm{exp}(ct)`$. Inserting Eq. (11) in Eq. (8) and reintroducing $`\sigma `$, we get

$$\mathrm{\Delta }(t)\frac{3F}{C^{}\gamma \sigma }\left(\frac{1}{2}\sqrt{\frac{\pi }{C^{}}}(\sigma t)^{1/6}e^{C^{}(\sigma t)^{1/3}}\right).$$
(12)
Notice that for very large $`t`$ the stretching $`\mathrm{\Delta }(t)`$ of the SWRN tends to a constant $`\mathrm{\Delta }_{\mathrm{}}`$. In the units of our figure this constant depends mainly on $`C^{}`$ (since $`F`$, $`\gamma `$ and $`\sigma `$ drop out); theoretically one may obtain $`C^{}`$ and also $`C`$ out of $`\mathrm{\Delta }_{\mathrm{}}`$.

In Fig. 4 we plot the dimensionless stretching $`\mathrm{\Delta }^{}(t)Y^{}(t)\sigma t/N`$ for $`p=0.05`$ and compare it to the analytical form Eq. (12). We do this by fitting $`abt^{1/6}\mathrm{exp}(ct^{1/3})`$ to the data, and as can be inferred from Fig. 4, the agreement is very convincing. From the least-squares fit we obtain $`a=5.09`$, $`b=7.67`$ and $`c=0.54`$. We remark that the agreement at short times may be rendered even better by also keeping the next term in the expansion of the integral of $`P_0(t)`$, a term which is proportional to $`t^{1/6}\mathrm{exp}\left(C^{}t^{1/3}\right)`$. Furthermore we remark that for Cayley-trees the intermediate behavior of $`Y(t)`$ can also be determined in a similar manner: The saddle-point approach yields to leading order:

$$Y(t)\stackrel{~}{a}\stackrel{~}{b}t^{3/2}\mathrm{exp}(\stackrel{~}{c}t),$$
(13)
with $`\stackrel{~}{a}`$, $`\stackrel{~}{b}`$ and $`\stackrel{~}{c}`$ being constants.

## IV Conclusions

In this communication we have studied the behavior of a small-world network model (the SWRN) of a linked ring-polymer, focusing on its dynamics under external forces. Thus the motion of a monomer pulled by such a force is vastly different, depending on whether the monomer belongs to a SWRN or to a simple ring without cross links. This may enable via NMR or electronic energy transfer to distinguish clearly between cross-linked and non-cross-linked polymers. Our numerical results for the stretching of the SWRN under external forces are in excellent agreement with our analytical expressions, which used the pseudo-gap behavior of the SWN density of states, as postulated in former work. As discussed in the present communication, these results are directly connected to expressions for the return to the origin of a random walker on the SWN.

###### Acknowledgements.

The support of the DFG, of the GIF through grant I0423-061.14, and of the Fonds der Chemischen Industrie are gratefully acknowledged.

FIGURES

CAPTIONS

FIG. 1. Two different realizations give rise to similar behavior of $`Y(t)`$, here plotted on logarithmic-logarithmic scales for $`p=0.05`$.

FIG. 2. On double logarithmic scales we plot the position $`Y(t)`$ as a function of time. From upper to lower curve, $`p=0`$, $`p=0.01`$ and $`p=0.05`$.

FIG. 3. On double logarithmic scales we plot the position $`Y(t)`$ as a function of time. From upper to lower curve, $`p=0.01`$, $`p=0.05`$, $`p=0.1`$, $`p=0.2`$, $`p=0.5`$ and $`p=0.8`$.

FIG. 4. Comparison of the theoretical prediction (dashed) with the data (dash-dotted), for $`p=0.05`$.