format: text

text:
# Generalized model for dynamic percolation

## I Introduction.

The percolation concept has turned out very useful for understanding transport and conduction processes in a wide range of disordered media, as exemplified by ionic conduction in polymeric, amorphous or glassy ceramic electrolytes, diffusion in biological tissues and permeability of disordered membranes .

Most of the situations discussed in Refs. pertain, however, to systems with ”frozen” disorder; that is, the random environment in which a given transport process takes place does not change in time. This is certainly the case in many instances, but it is not true in general. As a matter of fact, there are many experimental systems in which the static percolation picture does not apply since the structure of the host material undergoes essential structural reorganizations on a time scale comparable to that at which the transport itself occurs. A few stray examples of such systems include certain biomembranes , solid protonic conductors , oil-continuous microemulsions  and polymer electrolytes .

More specifically, ionic transport across a biomembrane, such as, e.g. gramicidin-$`A`$, occurs by the motion of ions through molecular channels along which they encounter potential barriers that fluctuate in time. The fluctuations of potential barriers may hinder significantly the transport and constitute an important transport-controlling factor . In the case of protonic conduction by the Grotthus mechanism , site-to-site motion of carriers occurs only between those neighboring $`H_2O`$ or $`NH_3`$ groups that have a favorable relative orientation; thermally activated rotation of these groups is the structural host-reorganization process interacting with the carrier motion. Similarly, within oil-continuous microemulsions, the charge transport proceeds by charge being transfered from one water globule to another, as globules approach each other in their Brownian motion . Lastly, in polymer electrolytes, such as, e.g. polyethylene oxide complexed nonstoichiometrically with the ionic salt $`NaSCN`$, the $`Na^+`$ ions are largely tetrahedrally coordinated by polyether oxygenes, but at the same time that $`Na^+`$ ions hop from one fourfold coordination site to another, the oxygens themselves, along with the polymeric backbone, undergo large-amplitude wagging and even diffusive motion .

Clearly, all the above mentioned examples involve two characteristic time scales, one which describes the typical time $`\tau `$ between two successive hops of the carrier, while the other is associated with a typical renewal time $`\tau ^{}`$ of the environment itself; namely, the time needed for the host medium to re-organize itself and thereby provide a new set of available pathways for transport. Consequently, the static percolation picture applies only when the characteristic time $`\tau ^{}`$ gets infinitely large. For a finite $`\tau ^{}`$ $`dynamicpercolation`$ has to be considered, and one encounters quite a different behavior when compared to the random environments with quenched disorder. As a result, one observes Ohmic-type or Stokes-type linear velocity-force relation for the carrier’s terminal velocities as a function of the applied field, in contrast to the threshold behavior predicted by the static percolation theory. The prefactor in the linear velocity-force relation may depend, however, in a non-trivial way on the system’s parameters and this dependence consitutes the main challenge for the theoretical analysis here. On the other hand, we note that in the above mentioned examples of the dynamic percolative environments quite different physical processes are responsible for the time evolution of the host medium. Consequently, one expects that the prefactor in the Stokes-type velocity-force relation should also be dependent on the precise mechanism which underlies the temporal re-organization of the environment.

Theoretical modelling of dynamic percolative environments has followed several avenues, which differ mostly in how the time evolution of the disorder is constrained; Namely, is it constrained (a) by conservation laws or (b) by spatial and temporal correlations in the renewal events? Early models of dynamic percolation  described the random environment within the framework of a standard bond-percolation model, in which the strength of each bond fluctuates in time between zero and a finite value. The dynamics of the host medium in these models  was accounted for by a series of instantaneous renewal events. These events were assumed to occur at random times, chosen from a renewal time distribution. In the renewal process the positions of all unblocked bonds are being reassigned, such that after each renewal event a carrier sees a newly defined network. This approach is thus characterized by a $`globaldynamicaldisorder`$ without global conservation laws and correlations, since the entire set of random hopping rates is simultaneously renewed independently of the previous history. Another model characterized by a $`localdynamicaldisorder`$ has been proposed in Refs. and , and subsequently generalized to the non-Markovian case in Ref.. This model appears to be similar to the previous one, except that here the hopping rates at different sites fluctuate $`independently`$ of each other. That is, individual bonds, rather than the whole lattice change in the renewal events. To describe the dynamical behavior in the $`localdynamicaldisorder`$ case, a dynamical mean-field theory has been proposed , based on the effective medium approximation introduced for the analysis of random walks on lattices with static disorder , and has been generalized to include the possibility of multistate transformations of the $`dynamically`$ random medium . More recently, several exactly solvable one-dimensional models with $`global`$ and $`local`$ dynamical disorder have been discussed .

In the second approach, which emerged within the context of the ionic conductivity in superionic solids, the dynamical percolative environment has been considered as a multicomponent mixture of mobile species in which one or several neutral components block the carrier component . In particular, such a situation can be observed in a superionic conductor $`\beta ^{\prime \prime }`$-alumina, doped with two different ionic species (e.g. $`Na^+`$ and $`Ba^{2+}`$), where small $`Na^+`$ ions are rather mobile, while the larger $`Ba^{2+}`$ ions move essentially slower and temporarily block the $`Na^+`$ ions. Contrary to the previous line of thought, the dynamics of such a percolative environment has essential correlations, generated by hard-core exclusion interactions between the species involved, and moreover, it obeys the conservation law \- the total number of the particles involved is conserved. In Ref., the frequency-dependent ionic conductivity of the light species has been analysed combining a continuous time random walk approach for the dynamical problem with an effective medium approximation describing the frozen environment of slow species. Next, as an explanation of the sharp increase of electrical conductivity transition in water-in-oil microemulsions when the volume fraction of water is increased towards a certain threshold value, in Refs. and  it has been proposed that the charge carriers are not trapped in the finite water clusters, but rather a charge on a water globule can propagate by either hopping to a neighboring globule, when they approach each other, or via the diffusion of the host globule itself. This picture has been interpreted in terms of a model similar to that employed in Ref., with the only difference being that here the ”blockers” of Ref. play the role of the transient charge carriers. In the model of Refs. and , in which the host dynamics is influenced by spatial correlations and conservation of the number of the water globules involved, the conductivity depends hence, on the rate of cluster rearrangement. Lastly, a similar problem of a carrier diffusion in an environment created by mobile hard-core lattice-gas particles has been analysed in Ref. by using the developed dynamic bond percolation theory of Refs. and .

In this paper we propose a generalized model of dynamic percolation which shares common features with both bond-fluctuating models of Refs. as well as models involving mobile blockers of Refs.. The system we consider consists of a host lattice, which here is a regular cubic lattice whose sites support at most a single occupancy, hard-core ”environment” particles, and a single hard-core carrier particle. The ”environment” particles move on the lattice by performing a random hopping between the neighboring lattice sites, which is constrained by the hard-core interactions, and may disappear from and re-appear (renewal processes) on the empty sites of the lattice with some prescribed rates***We hasten to remark that diffusive processes, of course, also result in a certain renewal of the environment; diffusive processes, as compared to the spontaneous creation and annihilation of particles, have however completelely different underlying physics and influence in a completely different fashion the evolution of the system, as we proceed to show. Following the terminology of Refs., we thus choose here to distinguish between diffusive and creation/annihilation processes, referring to the latter as the renewal ones.. In turn, the carrier particle is always present on the lattice, i.e., it can not disappear spontaneously, and is subject to a constant external force $`\stackrel{}{E}`$. Hence, the carrier performs a biased random walk, which is constrained by the hard-core interactions with the ”environment” particles, and probes the response of the percolative environment to the internal perturbancy or, in other words, the frictional properties of such a dynamical environment.

An important aspect of our model, which makes it different to the previously proposed models of dynamic percolation, is that we include the hard-core interaction between ”environment” particles and the carrier molecule, such that the latter may influence the dynamics of the environment. This results, as we proceed to show, in the emergence of complicated density profiles of the ”environment” particles around the carrier. These profiles, as well as the terminal velocity $`V_c`$ of the carrier, are determined here explicitly, in terms of an approximate approach of Ref., which is based on the decoupling of the carrier-particle-particle correlation functions into the product of pair-wise correlations. We show that the ”environment” particles tend to accumulate in front of the driven carrier creating a sort of a ”traffic jam”, which impedes its motion. Thus the density profiles around the carrier are highly asymmetric: the local density of the ”environment” particles in front of the carrier is higher than the average and approaches the average value as an exponential function of the distance from the carrier. The characteristic length and the amplitude of the density relaxation function are calculated explicitly. On the other hand, past the carrier the local density is lower than the average: We show that depending on the condition whether the number of particles in the percolative environment is explicitly conserved or not, the local density past the carrier may tend to the average value either as an exponential or even as an $`\mathrm{a}\mathrm{l}\mathrm{g}\mathrm{e}\mathrm{b}\mathrm{r}\mathrm{a}\mathrm{i}\mathrm{c}`$ function of the distance, revealing in the latter case especially strong memory effects and strong correlations between the particle distribution in the environment and the carrier position. Further on, we find that the terminal velocity of the carrier particle depends explicitly on the excess density in the ”jammed” region in front of the carrier, as well as on the ”environment” particles density past the carrier. Both, in turn, are dependent on the magnitude of the velocity, as well as on the rate of the renewal processes and the rate at which the ”environment” particles can diffuse away from the carrier. The interplay between the jamming effect of the environment, produced by the carrier particle, and the rate of its homogenization due to diffusive smoothening and renewal processes, manifests itself as a medium-induced frictional force exerted on the carrier, whose magnitude depends on the carrier velocity. As a consequence of such a non-linear coupling, in the general case, (i.e. for arbitrary rates of the renewal and diffusive processes), $`V_c`$ can be found only implicitly, as the solution of a non-linear equation relating $`V_c`$ to the system parameters. This equation simplifies considerably in the limit of small applied external fields $`\stackrel{}{E}`$ and we find that the force-velocity relation to the field becomes linear. This implies that the frictional force exerted on the carrier particle by the environment is $`\mathrm{v}\mathrm{i}\mathrm{s}\mathrm{c}\mathrm{o}\mathrm{u}\mathrm{s}`$. This linear force-velocity relation can be therefore interpreted as the analog of the Stokes formula for the dynamic percolative environment under study; in this case, the carrier velocity is calculated explicitly as well as the corresponding friction coefficient. In turn, this enable us to estimate the self-diffusion coefficient of the carrier in absence of external field; we show that when only diffusive re-arrangement of the percolative environment is allowed, while the renewal processes are suppressed, the general expression for the diffusion coefficient reduces to the one obtained previously in Refs. and . We note that the result of Refs. and  is known to serve as a very good approximation for the self-diffusion coefficient in hard-core lattice-gases .

We finally remark, that a qualitatively similar physical effect was predicted recently for a different model system involving a charged particle moving at a constant speed a small distance above the surface of an incompressible, infinitely deep liquid. It has been shown in Refs., that the interactions between the moving particle and the fluid molecules induce an effective frictional force exerted on the particle, producing a local distortion of the liquid interface, - a bump, which travels together with the particle and increases effectively its mass. The mass of the bump, which is analogous to the jammed region appearing in our model, depends itself on the particle’s velocity resulting in a non-linear coupling between the medium-induced frictional force exerted on the particle and its velocity .

The paper is structured as follows: In Section II we formulate the model and introduce basic notations. In Section III we write down the dynamical equations which govern the time evolution of the ”environment” particles and of the carrier. Section IV is devoted to the analytical solution of these evolution equations in the limit $`t\mathrm{}`$ ; here we also present some general results on the shape of the density profiles around stationary moving carrier and on the carrier terminal velocity, which is given implicitly, as the solution of a transcendental equation defining the general force-velocity relation for the dynamic percolative environment under study. In Section V we derive explicit asymptotic results for the carrier terminal velocity in the limit of small applied external fields $`\stackrel{}{E}`$ and obtain the analog of the Stokes formula for such a percolative environment; as well, we present here explicit results for the friction coefficient of the host medium and for the self-diffusion coefficient of the carrier in the absence of external field. Asymptotic behavior of the density profiles of the ”environment” particles around the carrier is discussed in Section VI. Finally, we conclude in Section VII with a brief summary and discussion of our results.

## II The model.

The model for dynamic percolation we study here consists of a three-dimensional simple cubic lattice of spacing $`\sigma `$, the sites of which are partially occupied by identical hard-core ”environment” particles and a single, hard-core, carrier particle (see Fig.1). For both types of particles the hard-core interactions prevent multiple occupancy of the lattice sites; that is, no two ”environment” particles or the ”carrier” and an ”environment” particle can occupy simultaneously the same site, and particles can not pass through each other.

The occupation of the lattice sites by the ”environment” particles is characterized by the time-dependent occupation variable $`\eta (\stackrel{}{r})`$, $`\stackrel{}{r}`$ being the lattice-vector of the site in question. This variable assumes two values:

$$\eta (\stackrel{}{r})=\{\begin{array}{cc}1,\text{ if the site }\stackrel{}{r}\text{ is occupied}\hfill & \\ 0,\text{ if the site }\stackrel{}{r}\text{ is empty}\hfill & \end{array}$$
(1)
Next, we assume the following dynamics of the ”environment” particles: The particles can spontaneously disappear from the lattice, and may re-appear at random positions and random time moments, which is reminiscent of the host medium dynamics stipulated in Refs.. We refer to these two processes generally as renewal processes. In addition, the environment particles move randomly within the lattice by performing nearest-neighbor random walks constrained by the hard-core interactions, which is the main feature of the approach in Refs.. We stipulate that any of the ”environment” particles waits a time $`\delta \tau `$, which has an exponential probability distribution with a mean $`\tau ^{}`$, and then chooses from a few possibilities: (a) disappearing from the lattice at rate $`g`$, which is realized instantaneously, or (b) attempting to hop, at rate $`l/6`$, onto one of $`6`$ neighboring sites. The hop is actually fulfilled if the target site is not occupied at this time moment by any other particle; otherwise, the particle attempting to hop remains at its initial position, and (c) particles may re-appear on any $`\mathrm{v}\mathrm{a}\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{t}`$ lattice site with rate $`f`$.

Note that, for simplicity, we assumed that the characteristic diffusion time and the renewal times of the ”environment” particles are equal to each other. These times, i.e. $`\tau _{dif}`$, mean creation time $`\tau _{cr}`$ and mean annihilation time $`\tau _{an}`$ may, however, be different, and can be restored in our final results by a mere replacement $`ll\tau ^{}/\tau _{dif}`$, $`ff\tau ^{}/\tau _{cr}`$ and $`gg\tau ^{}/\tau _{an}`$.

Note also that the number of particles is not explicitly conserved in such a dynamical model of the environment, which happens because of the presence of the renewal processes; the particles diffusion, on contrary, conserves the particles number. However, in the absence of attractive particle-particle interactions and external perturbances, the particles distribution on the lattice is uniform and the average occupation $`\rho (t)=\overline{\eta (\stackrel{}{r})}`$ of the lattice tends, as $`t\mathrm{}`$, to a constant value, $`\rho _s=f/(f+g)`$. This relation can be thought of as the Langmuir adsorption isotherm .

Hence, the limit $`\tau _{dif}\mathrm{}`$ (or, $`l0`$) corresponds to the ordinary site percolation model with immobile blocked sites. The limit $`f,g0`$, ($`\tau _{cr},\tau _{an}\mathrm{}`$), while keeping the ratio $`f/g`$ fixed, $`f/g=\rho _s/(1\rho _s)`$, corresponds to the usual hard-core lattice-gas with the number of particles conserved.

At time $`t=0`$ we introduce at the origin of the lattice an extra particle, the carrier, whose role is to probe the response of the environment modeled by dynamic percolation to an external perturbance. We stipulate that only the carrier out of all participating particles can not disappear from the system, and moreover, its motion is biased by some external constant force. As a physical realization, we envisage that the carrier is charged, while all other particles are neutral, and the system is exposed to constant external electric field $`\stackrel{}{E}`$. The dynamics of the carrier particle is defined as follows: We suppose that the waiting time between successive jumps of the carrier has also an exponential distribution with a mean value $`\tau `$, which may in general be different from the corresponding waiting time of the environment particles. Attempting to hop, the carrier first chooses a hop direction with probabilities

$$p_\mu =exp\left[\frac{\beta }{2}\left(\stackrel{}{E}\stackrel{}{e}_\mu \right)\right]/\underset{\nu }{}exp\left[\frac{\beta }{2}\left(\stackrel{}{E}\stackrel{}{e}_\nu \right)\right],$$
(2)
where $`\beta `$ is the reciprocal temperature, $`\stackrel{}{e}_\nu `$ (or $`\stackrel{}{e}_\mu `$) stand for six unit lattice vectors, $`\nu ,\mu =\{\pm 1,\pm 2,\pm 3\}`$, connecting the carrier position with $`6`$ neighboring lattice sites, and $`(\stackrel{}{E}\stackrel{}{e}_\nu )`$ denotes the scalar product. We adopt the convention that $`\pm 1`$ corresponds to $`\pm X`$, $`\pm 2`$ corresponds to $`\pm Y`$ while $`\pm 3`$ stands for $`\pm Z`$. The jump is actually fulfilled when the target lattice site is vacant. Otherwise, as mentioned the carrier remains at its position. For simplicity we assume in what follows that the external field is oriented along the $`X`$-axis in the positive direction, such that $`\stackrel{}{E}=(E,0,0)`$. Note also that for the choice of the transition probabilities as in Eq.(2), the detailed balance is naturally preserved.

## III Evolution Equations

Let $`P(\stackrel{}{R}_c,\eta ;t)`$ denote the joint probability that at time moment $`t`$ the carrier occupies position $`\stackrel{}{R}_c`$ and all ”environment” particles are in configuration $`\eta \{\eta (\stackrel{}{r})\}`$. Next, let $`\eta ^{\stackrel{}{r},\mu }`$ denote particles’ configuration obtained from $`\eta `$ by exchanging the occupation variables of the sites $`\stackrel{}{r}`$ and $`\stackrel{}{r}+\stackrel{}{e}_\stackrel{}{\mu }`$, i.e. $`\eta (\stackrel{}{r})\eta (\stackrel{}{r}+\stackrel{}{e}_\stackrel{}{\mu })`$, and $`\widehat{\eta }^\stackrel{}{r}`$ be the configuration obtained from $`\eta `$ by changing the occupation of the site $`\stackrel{}{r}`$ as $`\eta (\stackrel{}{r})1\eta (\stackrel{}{r})`$. Clearly, the first type of process appears due to random hops of the ”environment” particles, while the second one stems from the renewal processes, i.e. random creation and annihilation of the ”environment” particles. Then, summing up all possible events which can result in the configuration $`(\stackrel{}{R}_c,\eta )`$ or change this configuration for any other, we find that the temporal evolution of the system under study is governed by the following master equation:

$`_tP(\stackrel{}{R}_c,\eta ;t)={\displaystyle \frac{l}{6\tau ^{}}}{\displaystyle \underset{\mu }{}}{\displaystyle \underset{\stackrel{}{r}\stackrel{}{R}_c\stackrel{}{e}_\stackrel{}{\mu },\stackrel{}{R}_c}{}}\left\{P(\stackrel{}{R}_c,\eta ^{\stackrel{}{r},\mu };t)P(\stackrel{}{R}_c,\eta ;t)\right\}+`$ (3)
$`+`$ $`{\displaystyle \frac{1}{\tau }}{\displaystyle \underset{\mu }{}}p_\mu \left\{\left(1\eta (\stackrel{}{R}_c)\right)P(\stackrel{}{R}_c\stackrel{}{e}_\stackrel{}{\mu },\eta ;t)\left(1\eta (\stackrel{}{R}_c+\stackrel{}{e}_\stackrel{}{\mu })\right)P(\stackrel{}{R}_c,\eta ;t)\right\}+`$ (4)
$`+`$ $`{\displaystyle \frac{g}{\tau ^{}}}{\displaystyle \underset{\stackrel{}{r}\stackrel{}{R}_c}{}}\left\{\left(1\eta (\stackrel{}{r})\right)P(\stackrel{}{R}_c,\widehat{\eta }^\stackrel{}{r};t)\eta (\stackrel{}{r})P(\stackrel{}{R}_c,\eta ;t)\right\}+`$ (5)
$`+`$ $`{\displaystyle \frac{f}{\tau ^{}}}{\displaystyle \underset{\stackrel{}{r}\stackrel{}{R}_c}{}}\left\{\eta (\stackrel{}{r})P(\stackrel{}{R}_c,\widehat{\eta }^\stackrel{}{r};t)\left(1\eta (\stackrel{}{r})\right)P(\stackrel{}{R}_c,\eta ;t)\right\}.`$ (6)
Note that the terms in the first (resp. second) line of Eq.(6) describe random hopping motion of the ”environment” particles (resp. biased motion of the carrier) in terms of the Kawasaki-type particle-vacancy exchanges, while the terms in the third and the fourth lines account for the Glauber-type decay and creation of the ”environment” particles.

### A Mean velocity of the carrier and correlation functions.

From Eq.(6) we can readily compute the velocity of the carrier. Multiplying both sides of Eq.(6) by $`(\stackrel{}{R}_c\stackrel{}{e}_1)`$ and summing over all possible configurations $`(\stackrel{}{R}_c,\eta )`$ we find that the carrier’s mean velocity $`V_c(t)`$, defined as

$$V_c(t)\frac{\mathrm{d}}{\mathrm{dt}}\left(\overline{\stackrel{}{R}_c\stackrel{}{e}_1}\right),$$
(7)
obeys:

$$V_c(t)=\frac{\sigma }{\tau }\left\{p_1(1k(\stackrel{}{e}_1;t))p_1(1k(\stackrel{}{e}_1;t))\right\},$$
(8)
where $`k(\stackrel{}{\lambda };t)`$ stands for the carrier-”environment” particles pair correlation function

$$k(\stackrel{}{\lambda };t)\underset{\stackrel{}{R}_c,\eta }{}\eta (\stackrel{}{R}_c+\stackrel{}{\lambda })P(\stackrel{}{R}_c,\eta ;t).$$
(9)
In other words, $`k(\stackrel{}{\lambda };t)`$ can be thought of as the density distribution of the ”environment” particles as seen from the carrier which moves with velocity $`V_c(t)`$.

Hence, $`V_c(t)`$ depends explicitly on the local density of the ”environment” particles in the immediate vicinity of the carrier. Note that if the ”environment” is perfectly homogeneous, i.e. if for any $`\stackrel{}{\lambda }`$ the density profile is constant, $`k(\stackrel{}{\lambda };t)=\rho _s`$, which immediately implies decoupling between $`\eta (\stackrel{}{R}_c+\stackrel{}{\lambda })`$ and $`P(\stackrel{}{R}_c,\eta ;t)`$ in Eq.(9), then we obtain from Eq.(8) a trivial mean-field-type result

$$V_c^{(0)}=(p_1p_1)(1\rho _s)\frac{\sigma }{\tau },$$
(10)
which states that the frequency of jumps of the carrier particles ($`\tau ^1`$) only gets renormalized by a factor $`1\rho _s`$, which gives the fraction of successful jumps.

The salient feature of our model is that there are essential backflow effects. The carrier effectively perturbs the spatial distribution of the ”environment” particles so that stationary density profiles emerge. This can be contrasted to the earlier dynamic percolation models  in which the carrier had no impact on the embedding medium and hence there was no re-arrengement of the host medium around the carrier particle. As a consequence $`k(\stackrel{}{\lambda };t)\rho _s`$, and $`k(\stackrel{}{\lambda };t)`$ approaches $`\rho _s`$ only at infinite separations from the carrier, i.e. when $`|\stackrel{}{\lambda }|\mathrm{}`$. Therefore, we rewrite Eq.(8) in the form

$$V_c(t)=V_c^{(0)}\frac{\sigma }{\tau }\{p_1(k(\stackrel{}{e}_1;t)\rho _s)p_1(\rho _sk(\stackrel{}{e}_1;t))\},$$
(11)
which shows explicitly the deviation of the mean velocity of the carrier from the mean-field-type result in Eq.(10) due to the formation of the density profiles.

### B Evolution equations of the pair correlation functions.

From Eq.(8) it follows that in order to obtain $`V_c(t)`$, it suffices to compute $`k(\stackrel{}{e}_{\pm 1};t)`$. Consequently, we have to evaluate the equation governing the time evolution of the pair correlation functions. Multiplying both sides of Eq.(6) by $`\eta (\stackrel{}{R}_c)`$ and summing over all configurations $`(\stackrel{}{R}_c,\eta )`$, we find that $`k(\stackrel{}{\lambda };t)`$ obeys

$`_tk(\stackrel{}{\lambda };t)`$ $`=`$ $`{\displaystyle \frac{l}{6\tau ^{}}}{\displaystyle \underset{\mu }{}}(_\mu \delta _{\stackrel{}{\lambda },\stackrel{}{e}_\mu }_\mu )k(\stackrel{}{\lambda };t){\displaystyle \frac{(f+g)}{\tau ^{}}}k(\stackrel{}{\lambda };t)+{\displaystyle \frac{f}{\tau ^{}}}+`$ (12)
$`+`$ $`{\displaystyle \frac{1}{\tau }}{\displaystyle \underset{\mu }{}}{\displaystyle \underset{\stackrel{}{R}_c,\eta }{}}p_\mu \left(1\eta (\stackrel{}{R}_c+e_\mu )\right)_\mu \eta (\stackrel{}{R}_c+\stackrel{}{\lambda })P(\stackrel{}{R}_c,\eta ;t),`$ (13)
where $`_\mu `$ denotes the ascending finite difference operator of the form

$$_\mu f(\stackrel{}{\lambda })f(\stackrel{}{\lambda }+\stackrel{}{e}_\mu )f(\stackrel{}{\lambda }),$$
(14)
and

$$\delta _{\stackrel{}{r},\stackrel{}{r}^{}}=\{\begin{array}{cc}1,\text{ if the site }\stackrel{}{r}=\stackrel{}{r}^{}\hfill & \\ 0,\text{ otherwise.}\hfill & \end{array}$$
(15)
The Kroneker-delta term $`\delta _{\stackrel{}{\lambda },\stackrel{}{e}_\mu }`$ signifies that the evolution of the pair correlations, Eq.(13), proceeds differently at large separations and at the immediate vicinity of the carrier. This stems from the asymmetric hopping rules of the carrier particle defined by Eq.(2).

Note next that the contribution in the second line in Eq.(13), which is associated with the biased diffusion of the carrier, appears to be non-linear with respect to the occupation numbers, such that the pair correlation function gets effectively coupled to the evolution of the third-order correlations of the form

$$T(\stackrel{}{\lambda },\stackrel{}{e}_\nu ;t)\underset{\stackrel{}{R}_c,\eta }{}\eta (\stackrel{}{R}_c+\stackrel{}{\lambda })\eta (\stackrel{}{R}_c+\stackrel{}{e}_\nu )P(\stackrel{}{R}_c,\eta ;t).$$
(16)
That is, Eq.(13) is not closed with respect to the pair correlations but rather represents a first equation in the infinite hierachy of coupled equations for higher-order correlation functions. One faces, therefore, the problem of solving an infinite hierarchy of coupled differential equations and needs to resort to an approximate closure scheme.

### C Decoupling Approximation

Here we resort to the simplest non-trivial closure approximation, based on the decoupling of the third-order correlation functions into the product of pair correlations. More precisely, we assume that for $`\stackrel{}{\lambda }\stackrel{}{e}_\nu `$, the third-order correlation fulfils

$`{\displaystyle \underset{\stackrel{}{R}_c,\eta }{}}\eta (\stackrel{}{R}_c+\stackrel{}{\lambda })\eta (\stackrel{}{R}_c+\stackrel{}{e}_\nu )P(\stackrel{}{R}_c,\eta ;t)`$ (17)
$``$ $`\left({\displaystyle \underset{\stackrel{}{R}_c,\eta }{}}\eta (\stackrel{}{R}_c+\stackrel{}{\lambda })P(\stackrel{}{R}_c,\eta ;t)\right)\left({\displaystyle \underset{\stackrel{}{R}_c,\eta }{}}\eta (\stackrel{}{R}_c+\stackrel{}{e}_\nu )P(\stackrel{}{R}_c,\eta ;t)\right),`$ (18)
or, in other words,

$$\underset{\stackrel{}{R}_c,\eta }{}\eta (\stackrel{}{R}_c+\stackrel{}{\lambda })\eta (\stackrel{}{R}_c+\stackrel{}{e}_\nu )P(\stackrel{}{R}_c,\eta ;t)k(\stackrel{}{\lambda };t)k(\stackrel{}{e}_\nu ;t)$$
(19)
The approximate closure in Eq.(19) has been already employed for studying related models of biased carrier diffusion in hard-core lattice gases and has been shown to provide quite an accurate description of both the dynamical and stationary-state behavior. The decoupling in Eq.(19) was first introduced in Ref. to determine the properties of a driven carrier diffusion in a one-dimensional hard-core lattice gas with a conserved number of particles, i.e. without an exchange of particles with the reservoir. Extensive numerical simulations performed in Ref. have demonstrated that such a decoupling is quite a plausible approximation for the model under study. Moreover, rigorous probabilistic analysis of Ref. has shown that for this model the results based on the decoupling scheme in Eq.(19) are exact. Furthermore, the same closure procedure has been recently applied to study spreading of a hard-core lattice gas from a reservoir attached to one of the lattice sites . Again, a very good agreement between the analytical results and the numerical data has been found. Next, the decoupling in Eq.(19) has been used in a recent analysis of a biased carrier dynamics in a one-dimensional model of an adsorbed monolayer in contact with a vapour phase , i.e. a one-dimensional version of the model to be studied here. Also in this case an excellent agreement has been observed between the analytical predictions and the Monte Carlo simulations data . We now show that the approximate closure of the hierarchy of the evolution equations in Eq.(19) allows us to reproduce in the limit $`f,g=0`$ and $`f/g=const`$ the results of Refs. and , which are known (see e.g. Ref.) to provide a very good approximation for the carrier diffusion coefficient in three-dimensional hard-core lattice gases with arbitrary particle density. We expect therefore that such a closure scheme will render a plausible description of the carrier dynamics in a three-dimensional generalized dynamic percolation model. We base our further analysis on this approximation.

Making use of Eq.(19), we find from Eq.(13) that the pair correlations obey the following equations:

$$_tk(\stackrel{}{\lambda };t)=\frac{l}{6\tau ^{}}\stackrel{~}{L}k(\stackrel{}{\lambda };t)+\frac{f}{\tau ^{}},$$
(20)
which hold for all lattice sites except for those at the immediate vicinity of the carrier, i.e. for all $`\stackrel{}{\lambda }`$ except for $`\stackrel{}{\lambda }=\{\mathrm{𝟎},\stackrel{}{e}_{\pm 1},\stackrel{}{e}_{\pm 2},\stackrel{}{e}_{\pm 3}\}`$, while at the sites adjacent to the carrier one has

$$_tk(\stackrel{}{e}_\nu ;t)=\frac{l}{6\tau ^{}}\left(\stackrel{~}{L}+A_\nu (t)\right)k(\stackrel{}{e}_\nu ;t)+\frac{f}{\tau ^{}},$$
(21)
where $`\nu =\{\pm 1,\pm 2,\pm 3\}`$. The operators $`\stackrel{~}{L}`$ and coefficients $`A_\nu (t)`$ are given explicitly by

$$\stackrel{~}{L}\underset{\mu }{}A_\mu (t)_\mu \frac{6(f+g)}{l},$$
(22)
and

$$A_\mu (t)1+\frac{6\tau ^{}}{l\tau }p_\mu (1k(\stackrel{}{e}_\mu ;t)),$$
(23)
where $`_\mu `$ has been defined previously in Eq.(14), $`\mu =\{\pm 1,\pm 2,\pm 3\}`$. It is important to emphasize that all coefficients $`A_\mu (t)=A_\mu (E,V_c;t)`$, i.e. are functions of both the applied field and the carrier velocity.

Now, several comments about equations (20) and (21) are in order. First of all, let us note that Eq.(21) represents, from the mathematical point of view, the boundary conditions for the general evolution equation (20), imposed on the sites in the immediate vicinity of the carrier. Equations (20) and (21) have a different form since in the immediate vicinity of the carrier its asymmetric hopping rules perturb essentially the ”environment” particles dynamics. Equations (20) and (21) possess some intrinsic symmetries and hence the number of independent parameters can be reduced. Namely, reversing the field, i.e. changing $`EE`$, leads to the mere replacement of $`k(\stackrel{}{e}_1;t)`$ by $`k(\stackrel{}{e}_1;t)`$ but does not affect $`k(\stackrel{}{e}_\nu ;t)`$ with $`\nu =\{\pm 2,\pm 3\}`$, which implies that

$$k(\stackrel{}{e}_1;t)(E)=k(\stackrel{}{e}_1;t)(E),\text{and}k(\stackrel{}{e}_\nu ;t)(E)=k(\stackrel{}{e}_\nu ;t)(E)\text{for}\nu =\{\pm 2,\pm 3\},$$
(24)
Besides, since the transition probabilities in Eq.(2) obey

$$p_2=p_2=p_3=p_3$$
(25)
one evidently has that

$$k(\stackrel{}{e}_2;t)=k(\stackrel{}{e}_2;t)=k(\stackrel{}{e}_3;t)=k(\stackrel{}{e}_3;t),$$
(26)
and, by symmetry,

$$A_2(t)=A_2(t)=A_3(t)=A_3(t)$$
(27)
which somewhat simplifies equations (20) and (21). Lastly, we note that despite the fact that using the decoupling scheme in Eq.(19) we effectively close the system of equations on the level of the pair correlations, the solution of Eqs.(20) and (21) still poses serious technical difficulties. Namely, these equations are strongly non-linear with respect to the carrier velocity, which introduces the gradient term on the rhs of the evolution equations for the pair correlation, and depends by itself on the values of the ”environment” particles densities in the immediate vicinity of the carrier. Below we discuss a solution to this non-linear problem, focusing on the limit $`t\mathrm{}`$.

## IV Solution of the decoupled evolution equations in the stationary state.

Consider the limit $`t\mathrm{}`$ and suppose that the density profiles and the stationary velocity of the carrier have non-trivial stationary values

$$k(\stackrel{}{\lambda })\underset{t\mathrm{}}{lim}k(\stackrel{}{\lambda };t),V_c\underset{t\mathrm{}}{lim}V_c(t),\text{and}A_\mu =\underset{t\mathrm{}}{lim}A_\mu (t).$$
(28)
Define next the local deviations of $`k(\stackrel{}{\lambda })`$ from the unperturbed density as

$$h(\stackrel{}{\lambda })k(\stackrel{}{\lambda })\rho _s.$$
(29)
Choosing $`h(\mathrm{𝟎})=0`$, we obtain the following fundamental system of equations:

$$\stackrel{~}{L}h(\stackrel{}{\lambda })=0,$$
(30)
which holds for $`\stackrel{}{\lambda }\{\mathrm{𝟎},\stackrel{}{e}_{\pm 1},\stackrel{}{e}_{\pm 2},\stackrel{}{e}_{\pm 3}\}`$, while for the special sites adjacent to the carrier, i.e. for $`\stackrel{}{\lambda }=\{\mathrm{𝟎},\stackrel{}{e}_{\pm 1},\stackrel{}{e}_{\pm 2},\stackrel{}{e}_{\pm 3}\}`$, one has

$$(\stackrel{~}{L}+A_\nu )h(\stackrel{}{e}_\nu )+\rho _s(A_\nu A_\nu )=0,$$
(31)
Equations (30) and (31) determine the spatial distribution of the deviation from the unperturbed density $`\rho _s`$ in the stationary state. Note also that in virtue of the symmetry relations in Eqs.(26) and (27), $`h(\stackrel{}{e}_{\pm 2})=h(\stackrel{}{e}_{\pm 3})`$ and $`A_2=A_2=A_3=A_3`$.

The method for solving the coupled non-linear Eqs.(8),(30) and (31) is as follows: We first solve these equations supposing that the carrier stationary velocity is a given parameter, or, in other words, assuming that $`A_\nu `$ entering Eqs.(30) and (31) are known. In doing so, we obtain $`h(\lambda )`$ in the parametrized form

$$h(\stackrel{}{\lambda })=h(\stackrel{}{\lambda };A_{\pm 1},A_2).$$
(32)
Then, substituting into Eq.(32) particular values $`\stackrel{}{\lambda }=\{\stackrel{}{e}_{\pm 1},\stackrel{}{e}_{\pm 2},\stackrel{}{e}_{\pm 3}\}`$ and making use of the definition of $`A_\mu `$ in Eq.(23), we find a system of three linear equations with three unknowns of the form

$$A_\nu =1+\frac{6\tau ^{}}{l\tau }p_\nu \left(1\rho _sh(\stackrel{}{e}_\nu ;A_{\pm 1},A_2)\right),$$
(33)
where $`\nu =\{\pm 1,2\}`$, which will allow us to define all $`A_\nu `$ explicitly (and hence, all $`h(\stackrel{}{e}_\nu ))`$. Finally, substituting the results into Eq.(8), which can be written down in terms of $`A_\nu `$ as

$$V_c=\frac{l\sigma }{6\tau ^{}}(A_1A_1),$$
(34)
we arrive at a closed-form equation determining implicitly the stationary velocity.

### A Formal expression for the density profiles in the dynamic percolative environment as seen from the stationary moving carrier.

The general solution of Eqs.(30) and (31) can be most conveniently obtained by introducing the generating function

$$H(w_1,w_2,w_3)\underset{n_1,n_2,n_3}{}h(\stackrel{}{\lambda })w_1^{n_1}w_2^{n_2}w_3^{n_3},$$
(35)
where $`n_1`$,$`n_2`$ and $`n_3`$ are the components of the vector $`\stackrel{}{\lambda }`$, $`\stackrel{}{\lambda }=\stackrel{}{e}_1n_1+\stackrel{}{e}_2n_2+\stackrel{}{e}_3n_3`$. Multiplying both sides of Eqs.(30) and (31) by $`w_1^{n_1}w_2^{n_2}w_3^{n_3}`$ and performing summation, we find then that $`H(w_1,w_2,w_3)`$ is given explicitly by

$$H(w_1,w_2,w_3)=l\frac{_\nu \left(A_\nu (w_{|\nu |}^{\nu /|\nu |}1)h(\stackrel{}{e}_\nu )+\rho _s(A_\nu A_\nu )w_{|\nu |}^\nu \right)}{l_\nu A_\nu (w_{|\nu |}^{\nu /|\nu |}1)6(f+g)},$$
(36)
an expression which allows us to determine the stationary density profiles as seen from the carrier which moves with a constant velocity $`V_c`$.

Inversion of the generating function defined by Eq.(36) yields then, after rather lenghty but straightforward calculations, the following explicit result for the local deviation from the unperturbed density:

$`h(\stackrel{}{\lambda })`$ $`=`$ $`\alpha ^1\{{\displaystyle \underset{\nu }{}}A_\nu h(\stackrel{}{e}_\nu )_\nu `$ (37)
$``$ $`\rho _s(A_1A_1)(_1_1)\}F(\stackrel{}{\lambda }),`$ (38)
where $`F(\stackrel{}{\lambda })`$ is given by

$`F(\stackrel{}{\lambda })=`$ $`\left({\displaystyle \frac{A_1}{A_1}}\right)^{n_1/2}{\displaystyle _0^{\mathrm{}}}e^x\mathrm{I}_{n_1}\left(2{\displaystyle \frac{\sqrt{A_1A_1}}{\alpha }}x\right)\times `$ (39)
$`\times `$ $`\mathrm{I}_{n_2}\left(2{\displaystyle \frac{A_2}{\alpha }}x\right)\mathrm{I}_{n_3}\left(2{\displaystyle \frac{A_2}{\alpha }}x\right)\mathrm{d}x,`$ (40)
and

$`\alpha `$ $`=`$ $`{\displaystyle \underset{\nu }{}}A_\nu +{\displaystyle \frac{6(f+g)}{l}}=`$ (41)
$`=`$ $`A_1+A_1+4A_2+{\displaystyle \frac{6(f+g)}{l}}`$ (42)
Consequently, the particles density distribution as seen from the carrier moving with a constant velocity $`V_c`$ obeys

$`k(\stackrel{}{\lambda })`$ $`=`$ $`\rho _s+\alpha ^1\{{\displaystyle \underset{\nu }{}}A_\nu h(\stackrel{}{e}_\nu _\nu `$ (43)
$``$ $`\rho _s(A_1A_1)(_1_1)\}F(\stackrel{}{\lambda }),`$ (44)
where we have to determine three yet unknown parameters $`A_1`$, $`A_1`$ and $`A_2`$.

To determine these parameters, we set in Eq.(37) $`\stackrel{}{\lambda }=\stackrel{}{e}_1`$, $`\stackrel{}{\lambda }=\stackrel{}{e}_1`$ and $`\stackrel{}{\lambda }=\stackrel{}{e}_2`$, which results in the system of three closed-form equations determining the unknown functions $`A_\nu `$, $`\nu =\{\pm 1,2\}`$,

$$A_\nu =1+\frac{6\tau ^{}}{l\tau }p_\nu \left\{1\rho _s\rho _s(A_1A_1)\frac{det\stackrel{~}{C}_\nu }{det\stackrel{~}{C}}\right\},$$
(45)
where $`\stackrel{~}{C}`$ is a square matrix of the third order defined as

$$\left(\begin{array}{ccc}A_1_1F(\stackrel{}{e}_1)\alpha & A_1_1F(\stackrel{}{e}_1)& A_2_2F(\stackrel{}{e}_1)\\ A_1_1F(\stackrel{}{e}_1)& A_1_1F(\stackrel{}{e}_1)\alpha & A_2_2F(\stackrel{}{e}_1)\\ A_1_1F(\stackrel{}{e}_2)& A_1_1F(\stackrel{}{e}_2)& A_2_2F(\stackrel{}{e}_2)\alpha \end{array}\right)$$
(46)
while $`\stackrel{~}{C}_\nu `$ stands for the matrix obtained from $`\stackrel{~}{C}`$ by replacing the $`\nu `$-th column by a column vector $`\left((_1_1)F(\stackrel{}{e}_\nu )\right)_\nu `$. Equation (43), together with the definition of the coefficients $`A_\nu `$, constitutes the first general result of our analysis defining the density distribution in the percolative environment under study.

### B General force-velocity relation.

Substituting Eqs.(37) and (46) into (34), we find that the stationary velocity of the carrier particle is defined implicitly as the solution of equation:

$$V_c=\frac{\sigma }{\tau }(p_1p_1)(1\rho _s)\left\{1+\rho _s\frac{6\tau ^{}}{l\tau }\frac{p_1det\stackrel{~}{C}_1p_1det\stackrel{~}{C}_1}{det\stackrel{~}{C}}\right\}^1,$$
(47)
where $`\stackrel{~}{C}_1`$ and $`\stackrel{~}{C}_1`$ are the following square matrices of the third order:

$$\stackrel{~}{C}_1=\left(\begin{array}{ccc}(_1_1)F(\stackrel{}{e}_1)& A_1_1F(\stackrel{}{e}_1)& A_2_2F(\stackrel{}{e}_1)\\ (_1_1)F(\stackrel{}{e}_1)& A_1_1F(\stackrel{}{e}_1)\alpha & A_2_2F(\stackrel{}{e}_1)\\ (_1_1)F(\stackrel{}{e}_2)& A_1_1F(\stackrel{}{e}_2)& A_2_2F(\stackrel{}{e}_2)\alpha \end{array}\right)$$
(48)
and

$$\stackrel{~}{C}_1=\left(\begin{array}{ccc}A_1_1F(\stackrel{}{e}_1)\alpha & (_1_1)F(\stackrel{}{e}_1)& A_2_2F(\stackrel{}{e}_1)\\ A_1_1F(\stackrel{}{e}_1)& (_1_1)F(\stackrel{}{e}_1)& A_2_2F(\stackrel{}{e}_1)\\ A_1_1F(\stackrel{}{e}_2)& (_1_1)F(\stackrel{}{e}_2)& A_2_2F(\stackrel{}{e}_2)\alpha \end{array}\right).$$
(49)
Equation (47) represents our second principal result defining the force-velocity relation in the dynamic percolative environment for an arbitrary field and arbitrary rates of the diffusive and renewal processes.

## V Carrier velocity in the limit of small applied field $`E`$, friction coefficient and carrier diffusivity in dynamic percolative environment.

We consider now the case when the applied external field $`E`$ is small. Expanding the transition probabilities $`p_1`$ and $`p_1`$ in the Taylor series up to the first order in powers of the external field, i.e.

$$p_{\pm 1}=\frac{1}{6}\pm \frac{\sigma \beta E}{12}+𝒪\left(E^2\right),$$
(50)
we find that $`V_c`$ defined by Eq.(34) follows

$$V_c\frac{\sigma }{6\tau }\left\{\sigma \beta E(1\rho _s)(h(\stackrel{}{e}_1)h(\stackrel{}{e}_1))\right\}.$$
(51)
On the other hand, Eq.(37) entails that

$$h(\stackrel{}{e}_1)h(\stackrel{}{e}_1)=\frac{2\sigma \rho _s(1\rho _s)\tau ^{}}{l\tau \left(\alpha _0(2A_0/\alpha _0)A_0\right)+2\rho _s\tau }\beta E+𝒪\left(E^2\right),$$
(52)
where

$$A_0=lim_{E0}A_\nu =1+\frac{\tau ^{}}{l\tau }(1\rho _s),$$
(53)
and

$$\alpha _0=lim_{E0}\alpha =6\left(1+\frac{\tau ^{}(1\rho _s)}{l\tau }+\frac{f+g}{l}\right),$$
(54)
while

$`(x)`$ $``$ $`\left\{{\displaystyle _0^{\mathrm{}}}e^t\mathrm{I}_0^2(xt)\left(\mathrm{I}_0(xt)\mathrm{I}_2(xt)\right)dt\right\}^1`$ (55)
$`=`$ $`\left\{P(\mathrm{𝟎};3x)P(2\stackrel{}{e}_1;3x)\right\}^1,`$ (56)
$`P(\stackrel{}{r};\xi )`$ being the generating function,

$$P(\stackrel{}{r};\xi )\underset{j=0}{\overset{+\mathrm{}}{}}P_j(\stackrel{}{r})\xi ^j,$$
(57)
of the probability $`P_j(\stackrel{}{r})`$ that a walker starting at the origin and performing a Polya random walk on the sites of a three-dimensional cubic lattice will arrive on the $`j`$-th step to the site with the lattice vector $`\stackrel{}{r}`$ .

Consequently, we find that in the limit of a small applied field $`E`$ the force-velocity relation in Eq.(47) attains the physically meaningful form of the Stokes formula $`E=\zeta V_c`$, which signifies that the frictional force exerted on the carrier by the environment particles is $`\mathrm{v}\mathrm{i}\mathrm{s}\mathrm{c}\mathrm{o}\mathrm{u}\mathrm{s}`$. The effective friction coefficient $`\zeta `$ is the sum of two terms,

$$\zeta =\zeta _0+\zeta _{coop}$$
(58)
where the first term represents a mean-field-type result $`\zeta _0=6\tau /\beta \sigma ^2(1\rho _s)`$ (see Eq.(7)), while the second one, $`\zeta _{coop}`$, obeys

$$\zeta _{coop}=\frac{12\rho _s\tau ^{}}{\beta \sigma ^2l(1\rho _s).\left(\alpha _0(2A_0/\alpha _0)A_0\right)}$$
(59)
The second contribution has a more complicated origin and is associated with the cooperative behavior - formation of a inhomogeneous stationary particle distribution around the carrier moving with constant velocity $`V_c`$. Needless to say, such an effect can not be observed within the framework of previous models of dynamic percolation, since there the carrier does not influence the host medium dynamics .

Let us now compare the relative importance of two contributions, i.e. $`\zeta _0`$ and $`\zeta _{coop}`$, to the overall friction. In Fig.2 we plot the ratio $`\zeta /\zeta _0`$ versus the creation rate $`f`$ for three different values of the density $`\rho _s`$, $`\rho _s=0.9,0.7`$ and $`0.5`$, while the annihilation rate is prescribed by the relation $`g=f(1\rho _s)/\rho _s`$. This figure shows that the cooperative behavior clearly dominates at small and moderate $`f`$ (which entails also small values of $`g`$), while for larger $`f`$, when $`\zeta /\zeta _0`$ tends to $`1`$, the mean-field behavior becomes most important. The cooperative behavior also appears to be more pronounced at larger densities $`\rho _s`$.

Consider next some analytical estimates. We start with the situation, in which diffusion of the environment particles is suppressed, i.e. when $`l=0`$. In this case, we get

$$\frac{\zeta _{coop}}{\zeta _0}=\frac{2\rho _s}{(1\rho _s)\left(\frac{2}{y}(y)1\right)},$$
(60)
where

$$y=\frac{1}{3}\left(1+\frac{\tau }{\tau ^{}}\frac{(f+g)}{(1\rho _s)}\right)^1.$$
(61)
Suppose first that $`\rho _s`$ is small, $`\rho _s1`$. Then, $`y1/3(1+\tau /\tau ^{}(f+g))`$ and we can distinguish between two situations: when $`\tau (f+g)/\tau ^{}`$, i.e. when the carrier moves faster than the environment re-organizes itself, and and the opposite limit, $`\tau (f+g)/\tau ^{}`$, when the environment changes very rapidly compared to the motion of the carrier. In the former case we find that $`y1/3`$, which yields $`\zeta _{coop}/\zeta _02\rho _s/(6(1/3)1)`$, $`(1/3)0.7942`$, while in the latter case we have $`y\tau ^{}/3\tau (f+g)`$ and $`\zeta _{coop}/\zeta _0\rho _s\tau ^{}/3\tau (f+g)`$. Note, that in both cases the ratio $`\zeta _{coop}/\zeta _0`$ appears to be small, which signifies that at small densities $`\rho _s`$ the mean-field friction dominates. Such a result is consistent with the behavior depicted in Fig.2 and is not counterintuitive, of course, since in the absence of the particles’ diffusion, which couples effectively the density evolution at different lattice sites, no significant cooperative behavior can emerge at small densities. On the other hand, at relatively high densities $`\rho _s1`$ and $`\tau /(1\rho _s)\tau ^{}/(f+g)\tau `$, when the carrier moves at much faster rate than the host medium reorganizes itself, we find that $`\zeta _{coop}/\zeta _0\tau ^{}/3\tau (f+g)1`$. This result stems from the circumstance that in sufficiently dense environments modeled by dynamic percolation a highly inhomogeneous density profile emerges even in the absence of particles diffusion; Here, on the one hand, the carrier perturbs significantly the particle density in its immediate vicinity. On the other hand, the density perturbance created by the carrier does not shift the global balance between creation and annihilation events, i.e. the mean particle density still equals $`\rho _s`$, as we set out to show in what follows. The latter constraint induces then appearence of essential correlations in particles distribution and hence, appearence of cooperative behavior.

Let us consider the opposite case when the renewal processes are not allowed, which means that the particles number is conserved and local density in the percolative environment evolves only due to particles diffusion. In this case we find

$$\frac{\zeta _{coop}}{\zeta _0}=\frac{2\tau ^{}\rho _s}{(l\tau +\tau ^{}(1\rho _s))\left(6(1/3)1\right)}$$
(62)
Here, the ratio $`\zeta _{coop}/\zeta _0`$ can be large and the ”cooperative” friction dominates the mean-field one when $`l\tau \tau ^{}(3\rho _s1)`$, which happens, namely, at sufficiently high densities and in the limit when the carrier moves at a much faster rate than the environment reorganizes itself. Otherwise, the mean-field friction prevails.

To estimate the carrier particle diffusion coefficient $`D_c`$ we assume the validity of the Einstein relation, i.e. $`\beta D_c=\zeta ^1`$. We find that, in the general case, the carrier diffusion coefficient $`D_c`$ reads

$$D_c=\frac{\sigma ^2(1\rho _s)}{6\tau }\left\{1\frac{2\rho _s\tau ^{}}{l\tau }\left(\alpha _0(2A_0/\alpha _0)1+\frac{\tau ^{}(3\rho _s1)}{l\tau }\right)^1\right\}$$
(63)
In the particular case of conserved particles number, when $`f,g0`$ but their ratio $`f/g`$ is kept fixed, $`f/g=\rho _s/(1\rho _s)`$, the latter equation reduces to the classical result

$$D_c^{NK}=\frac{\sigma ^2(1\rho _s)}{6\tau }\left\{1\frac{2\rho _s\tau ^{}}{l\tau }\left(6A_0(1/3)1+\frac{\tau ^{}(3\rho _s1)}{l\tau }\right)^1\right\},$$
(64)
obtained earlier in Refs. and  by different analytical techniques. The result in Eq.(64) is known to be exact in the limits $`\rho _s1`$ and $`\rho _s1`$, and serves as a very good approximation for the self-diffusion coefficient in hard-core lattice gases of arbitrary density .

It seems also interesting to analyse how random annihilation and creation of particles can modify the self-diffusion coefficient compared to the situation when the particles number is conserved. In Figure 3 we plot the ratio $`D_c^{NK}/D_c`$ ($`=\zeta /\zeta _{NK}`$) versus the creation rate $`f`$ for three different values of the density $`\rho _s`$, $`\rho _s=0.9,0.7`$ and $`0.5`$. Again, the value of the annihilation rate $`g`$ is prescribed by the relation $`g=f(1\rho _s)/\rho _s`$. Figure 3 shows that the renewal processes affect considerably the friction coefficient and the ratio $`\zeta /\zeta _{NK}`$ deviates strongly from the unity with the growth of the creation rate. The overall friction also falls off when the density increases.

Finally, in the absence of particle diffusion (fluctuating-site percolation), our result for the carrier particle diffusion coefficient reduces to

$$D_c^{per}=\frac{\sigma ^2(1\rho _s)}{6\tau }\left\{12\rho _s\left(4[(1\rho _s)+(f+g)\tau /\tau ^{}](y)+3\rho _s1\right)^1\right\}$$
(65)
Note, however, that this result only applies when both $`f`$ and $`g`$ are larger than zero, such that the renewal processes take place. In fact, the underlying decoupling scheme is only plausible in this case. Similarly to the approximate theories in Refs. and , our approach predicts that in the absence of the renewal processes $`D_c^{per}`$ vanishes only when $`\rho _s1`$, which is an incorrect behavior.

## VI Asymptotic behavior of the density profiles at large distances in front of and past the carrier.

The density profiles at large separations in front of and past the carrier can be readily deduced from the asymptotical behavior of the following generating function

$$N(w_1)\underset{n_1=\mathrm{}}{\overset{+\mathrm{}}{}}h(n_1,n_2=0,n_3=0)w_1^{n_1}.$$
(66)
Inversion of Eq.(36) with respect to the symmetric coordinates $`n_2`$ and $`n_3`$ yields then

$`N(w_1)={\displaystyle \frac{\left(A_1h(\stackrel{}{e}_1)+\rho _s(A_1A_1)\right)\left(w_11\right)+\left(A_1h(\stackrel{}{e}_1)\rho _s(A_1A_1)\right)\left(w_1^11\right)}{\alpha A_1w_1^1A_1w_1}}\times `$ (67)
$`\times `$ $`{\displaystyle _0^{\mathrm{}}}\mathrm{exp}[x]\mathrm{I}_0^2({\displaystyle \frac{2A_2}{\alpha A_1w_1^1A_1w_1}}x)\mathrm{d}x+{\displaystyle \frac{4A_2h(\stackrel{}{e}_2)}{\alpha A_1w_1^1A_1w_1}}\times `$ (68)
$`\times `$ $`{\displaystyle _0^{\mathrm{}}}\mathrm{exp}[x]\mathrm{I}_0({\displaystyle \frac{2A_2}{\alpha A_1w_1^1A_1w_1}}x)(\mathrm{I}_1({\displaystyle \frac{2A_2}{\alpha A_1w_1^1A_1w_1}}x)`$ (69)
$``$ $`\mathrm{I}_0({\displaystyle \frac{2A_2}{\alpha A_1w_1^1A_1w_1}}x))\mathrm{d}x`$ (70)
We notice now that $`N(w_1)`$ is a holomorphic function in the region $`𝒲_1<w_1<𝒲_2`$, where

$`𝒲_1={\displaystyle \frac{\alpha 4A_2}{2A_1}}\sqrt{\left({\displaystyle \frac{\alpha 4A_2}{2A_1}}\right)^2{\displaystyle \frac{A_1}{A_1}}}`$ (71)
and

$`𝒲_2={\displaystyle \frac{\alpha 4A_2}{2A_1}}+\sqrt{\left({\displaystyle \frac{\alpha 4A_2}{2A_1}}\right)^2{\displaystyle \frac{A_1}{A_1}}}`$ (72)
As a consequence, the asymptotic behavior of $`h(n_1,n_2=0,n_3=0)`$ in the limit $`n_1\mathrm{}`$ (resp. $`n_1\mathrm{}`$) is controlled by the behavior of $`N(w_1)`$ in the vicinity of $`w_1=𝒲_2`$ (resp. $`w_1=𝒲_1`$) (see, for example, the analysis of the generating function singularities developed in Ref.).

#### 1 Asymptotics of the density profiles at large separations in front of the carrier.

Consider first the asymptotic behavior of the density distribution of the ”environment” particles at large separations in front of the carrier. Using the fact that

$$_0^{\mathrm{}}exp\left[x\right]\mathrm{I}_0(yx)\left(\mathrm{I}_1(x)\mathrm{I}_0(x)\right)dx$$
(73)
is a regular function when $`y1/2`$, while

$$_0^{\mathrm{}}exp\left[x\right]\mathrm{I}_0^2(yx)dx\frac{1}{\pi }\mathrm{ln}\left(\frac{1}{12y}\right),$$
(74)
we find that

$`N(w_1)`$ $``$ $`{}_{w_1𝒲_2}{}^{}[{\displaystyle \frac{\left(A_1h(\stackrel{}{e}_1)+\rho _s(A_1A_1)\right)\left(𝒲_21\right)}{4\pi A_2}}+`$ (75)
$`+`$ $`{\displaystyle \frac{\left(A_1h(\stackrel{}{e}_1)\rho _s(A_1A_1)\right)\left(𝒲_2^11\right)}{4\pi A_2}}\left]\mathrm{ln}\right(𝒲_2w_1)`$ (76)
Then, (cf, Flajolet et al., Ref.), we obtain the following asymptotical result

$$h(n_1,0,0)_{n_1\mathrm{}}\frac{K^+}{n_1}e^{n_1/\lambda _+},$$
(77)
where

$$\lambda _+\mathrm{ln}^1\left(\frac{\alpha /22A_2}{A_1}+\sqrt{\left(\frac{\alpha /22A_2}{A_1}\right)^2\frac{A_1}{A_1}}\right),$$
(78)
and

$`K^+`$ $`=`$ $`[{\displaystyle \frac{\left(A_1h(\stackrel{}{e}_1)+\rho _s(A_1A_1)\right)\left(𝒲_21\right)}{4\pi A_2}}+`$ (79)
$`+`$ $`{\displaystyle \frac{(A_1h(\stackrel{}{e}_1)\rho _s(A_1A_1)(𝒲_2^11)}{4\pi A_2}}]>0,`$ (80)
which signifies that the density of the ”environment” particles in front of the carrier is higher than the average value $`\rho _s`$ and approaches $`\rho _s`$ at large separations from the carrier as an exponential function of the distance.

#### 2 Asymptotics of the density profiles at large separations past the carrier.

We consider next the asymptotic behavior of the ”environment” particles density profiles past the carrier particle, which turns out to be very different depending on whether the dynamics of the percolative environment obeys the strict conservation of the ”environment” particles number or not (the renewal processes are suppressed or allowed). The sketch of this behavior is presented in Fig.4.

##### a Non-conserved particles number.

In the case when partciles may disappear and re-appear on the lattice, one has that the root $`𝒲_1<1`$. We find then, following essentially the same lines as in the previous subsection, that

$`N(w_1)`$ $``$ $`{}_{w_1𝒲_1}{}^{}[{\displaystyle \frac{\left(A_1h(\stackrel{}{e}_1)+\rho _s(A_1A_1)\right)\left(𝒲_11\right)}{4\pi A_2}}+`$ (81)
$`+`$ $`{\displaystyle \frac{\left(A_1h(\stackrel{}{e}_1)\rho _s(A_1A_1)\right)\left(𝒲_1^11\right)}{4\pi A_2}}\left]\mathrm{ln}\right({\displaystyle \frac{1}{w_1𝒲_1}}).`$ (82)
Hence, in the non-conserved case the approach to the unperturbed value $`\rho _s`$ is also exponential when $`n_1\mathrm{}`$, and follows

$$h_{n_1,0,0}_{n_1\mathrm{}}\frac{K^{}}{|n_1|}e^{|n_1|/\lambda _{}},$$
(83)
where

$$\lambda _{}\mathrm{ln}^1\left(\frac{\alpha /22A_2}{A_1}\sqrt{\left(\frac{\alpha /22A_2}{A_1}\right)^2\frac{A_1}{A_1}}\right)$$
(84)
and

$`K^{}`$ $`=`$ $`[{\displaystyle \frac{\left(A_1h(\stackrel{}{e}_1)+\rho _s(A_1A_1)\right)\left(𝒲_11\right)}{4\pi A_2}}+`$ (85)
$`+`$ $`{\displaystyle \frac{\left(A_1h(\stackrel{}{e}_1)\rho _s(A_1A_1)\right)\left(𝒲_1^11\right)}{4\pi A_2}}]<0`$ (86)
which implies that the particles density past the carrier is lower than the average. Note that, in the general case, $`\lambda _+<\lambda _{}`$, which means that the depleted region past the carrier is more extended in space than the traffic-jam-like region in front of the carrier. The density profiles are therefore asymmetric with respect to the origin, $`n_1=0`$. Since creation of particles is favored (suppressed) in depleted (jammed) regions, while annihilation is suppressed (favored), one might expect that this will shift the overall density in the system, i.e. the average density of the ”environment” particles will differ from $`\rho _s`$. Interestingly, the overall deviation, i.e. the sum of local deviations over the volume of the system, of the density of the ”environment” particles from the average value $`\rho _s`$, appears to be equal exactly to zero,

$$H(w_1=1,w_2=1,w_3=1)0,$$
(87)
and hence, the driven carrier does not perturb the global balance between creation and annihilation of the ”environment” particles. This is not, however, an $`\mathrm{a}\mathrm{p}\mathrm{r}\mathrm{i}\mathrm{o}\mathrm{r}\mathrm{i}`$ evident result in view of the asymmetry of the density profiles.

##### b Conserved particles number.

Finally, we turn to the analysis of the shape of the density profiles of the percolative environment past the carrier in the particular limit when the host medium evolves only due to diffusion, while creation and annihilation of particles are completely suppressed. In this case, in which the particles number is explicitly conserved, one has that for arbitrary value of the field and particles’ average density, the root $`𝒲_11`$ and, consequently, the form of the generating function is qualitatively different from that in Eqs.(75) and (81),

$`N(w_1)`$ $``$ $`{}_{w_11^+}{}^{}[{\displaystyle \frac{(A_1h(\stackrel{}{e}_1)A_1h(\stackrel{}{e}_1)}{4\pi A_2}}+`$ (88)
$`+`$ $`{\displaystyle \frac{2\rho _s(A_1A_1))}{4\pi A_2}}](w_11)\mathrm{ln}({\displaystyle \frac{1}{w_11}}).`$ (89)
Equation (88) implies that in the limit when the particle number is conserved the large-$`n_1`$ asymptotic behavior of $`h_{n_1,0,0}`$ is described by an algebraic function of $`n_1`$ with a logarithmic correction; that is,

$$h_{n_1,0,0}\frac{K_{}\mathrm{ln}(|n_1|)}{n_1^2},$$
(90)
where $`K_{}`$ is an $`n_1`$-independent constant. Remarkably, the power-law decay of correlations implies existence of a quasi-long-range order in the percolative environment past the carrier. In the conserved case the mixing of the three-dimensional percolative environment is not very efficient and there are considerable memory effects - the host medium remembers the passage of the carrier on large space and time scales.

## VII Conclusions

To conclude, we have studied analytically the dynamics of a carrier driven by an external field $`\stackrel{}{E}`$ in a three-dimensional environment modeled by dynamic percolation on cubic lattice partially filled with mobile, hard-core ”environment” particles which can spontaneously disappear and reappear (renewal processes) in the system with some prescribed rates. Our analytical approach has been based on the master equation, describing the time evolution of the system, which has allowed us to evaluate a system of coupled dynamical equations for the carrier velocity and a hierarchy of correlation functions. To solve these coupled equations, we have invoked an approximate closure scheme based on the decomposition of the third-order correlation functions into a product of pairwise correlations, which has been first introduced in Ref. for a related model of a driven carrier dynamics in a one-dimensional lattice gas with conserved particles number. Within the framework of this approximation, we have derived a system of coupled, discrete-space equations describing evolution of the density profiles of the environment, as seen from the moving carrier, and its velocity $`V_c`$. We have shown that $`V_c`$ depends on the density of the ”environment” particles in front of and past the carrier. Both densities depend on the magnitude of the velocity, as well as on the rate of the renewal and diffusive processes. As a consequence of such a non-linear coupling, in the general case, (i.e. for an arbitrary driving field and arbitrary rates of renewal and diffusive processes), $`V_c`$ has been found only implicitly, as the solution of a non-linear equation relating its value to the system parameters. This equation, which defines the force-velocity relation for the dynamic percolation under study, simplifies considerably in the limit of small applied field $`\stackrel{}{E}`$. We find that in this limit it attains the physically meaningful form of the Stokes formula, which implies, in particular, that the frictional force exerted on the carrier by the environment modeled by dynamic percolation is $`viscous`$. In this limit, the carrier velocity and the friction coefficient are calculated explicitly. In addition, we determine the self-diffusion coefficient of the carrier in the absence of the field and show that it reduces to the well-know result of Refs. and  in the limit when the particles number is conserved. Further more, we have found that the density profile around the carrier becomes strongly inhomogeneous: the local density of the ”environment” particles in front of the carrier is higher than the average and approaches the average value as an exponential function of the distance from the carrier. On the other hand, past the carrier the local density is lower than the average, and depending on whether the number of particles is explicitly conserved or not, the local density past the carrier may tend to the average value either as an exponential or even as an $`\mathrm{a}\mathrm{l}\mathrm{g}\mathrm{e}\mathrm{b}\mathrm{r}\mathrm{a}\mathrm{i}\mathrm{c}`$ function of the distance. The latter reveals especially strong memory effects and strong correlations between the particle distribution in the environment and the carrier position.

Figure Captions.

Fig.1. A generalized model of dynamic percolation. Grey spheres denote the hard-core ”environment” particles, which perform symmetric random hopping among the sites of a simple cubic lattice, and can be spontaneously annihilated and created. The lighter sphere is the carrier, which performs a biased random walk due to an external field $`\stackrel{}{E}`$, consrained by hard-core exclusion with the ”environment” particles.

Fig.2. The ratio of the overall friction coefficient and the mean-field friction versus the creation rate for three different values of the mean density $`\rho _s`$. The upper curve corresponds to $`\rho _s=0.9`$, the intermediate - to $`\rho _s=0.7`$, and the lower - to $`\rho _s=0.5`$.

Fig.3. The ratio of the overall friction coefficient and the friction coefficient in the conserved particles number case versus the creation rate for three different values of the mean density $`\rho _s`$. The upper curve corresponds to $`\rho _s=0.5`$, the intermediate - to $`\rho _s=0.7`$, and the lower - to $`\rho _s=0.9`$.

Fig.4. A sketch of the asymptotic density profiles in front of and past the stationary moving carrier. The abbreviation CPN stands for the ”conserved particles number”. The two solid lines in the $`\lambda >0`$ and $`\lambda <0`$ domains denote exponential profiles, Eqs.(63) and (67). The dashed line in the domain $`\lambda <0`$ stands for the algebraic law, Eq.(72).


text:
# Flow equations in the light-front QCD: mass gap and confinement

## 1 Introduction

The perturbative aspects of non-abelian gauge theories were underestood many years ago, and the perturbative calculations provided convincing proof that QCD is the theory of strong interactions. However nonperturbative QCD phenomena have been difficult to analyze mainly because calculational techniques are still lacking, even though the qualitative features have been more or less understood.

In particular, it is widely believed that pure Yang-Mills theory, with no dynamical quarks, posseses the features like asymptotic freedom, mass generation through the transmutation of dimensions and confinement: linear potential between static (probe) quarks. In the past few years there were several studies addressing the issue of confinement and generation of mass gap in the framework of the Schrödinger picture $`^\mathrm{?}`$, $`^\mathrm{?}`$. To mention a few early works, see refs. $`^\mathrm{?}`$. We have tried to understand these nonperturbative mechanisms also in a Hamiltonian framework, solving flow equations for canonical QCD Hamiltonian in the light-front (LF) gauge selfconsistently for the few lowest sectors. Dynamics of quarks has been excluded to disentangle the complexity of chiral symmetry breaking. Early studies of confinement in the LF frame were conducted using similarity renormalization $`^\mathrm{?}`$ and transverse lattice $`^\mathrm{?}`$, based on the fact that QCD in $`3+1`$ already has a confining interaction in the form of the instantaneous four-fermion interaction, which is the confining interaction in $`1+1`$. However, in our study the instantaneous interaction is canceled by the corresponding term in the dynamical interaction, generated by flow equations. The result is a three dimensional linear rising potential in the infrared region.

Our basic strategy has been to generate dynamically an effective gluon mass (energy) through interaction with the LF ’vacuum’ (exact zero modes on the light-front, $`x=0`$), and construct an effective interaction between static quarks as an exchange with this dynamical gluon, which forms a flux between quarks. One could construct a LF field theory with massive photons and gluons, as was done in $`^\mathrm{?}`$, but this theory is in a Coulomb phase rather than confining because of the lack of a dynamical mass generation mechanism. Dynamical mass generation in the light-front frame dates back to the works of Cornwall $`^\mathrm{?}`$, where a kind of pinch technique was used to define a gluon propagator and hence a gluon effective energy, extracted from a Green’s function for physical observables (for example two-body scattering amplitude). It was suggested in $`^\mathrm{?}`$, that in physical terms, a gluon ’mass’ may lead to a vortex condensation aaaIn the massive gauge theory one needs to introduce a scalar field, which condenses, in order to preserve gauge invariance. Topologically this field represents vortices $`^\mathrm{?}`$.  with $`<TrG_{\mu \nu }^2>0`$, and conversely $`<TrG_{\mu \nu }^2>0`$ generates a gluon mass. This provides a link with the existing vortex picture of confinement.

A subtle point in the LF field theory is that the LF vacuum is just empty space. Therefore it seems a problem how confinement can occur in the LF frame, and what quantity sets up a scale for a dynamical mass gap and the string tension. We adopt a model, suggested by Susskind and Burkardt in the context of chiral symmetry breaking in the LF frame $`^\mathrm{?}`$. In the parton model one pictures a fast moving hadron as being some collection of constituents with relatively large momentum, such that when one boosts the system, doubles its momentum, all these partons double their momenta and so forth. Therefore one can formulate a field theory on the axis of the light-front momentum $`x`$. Partons that form a hadron are at positive, finite $`x`$ and according to Feynman and Bjorken fill the $`x`$-axis in a way which gets denser and denser as one goes to smaller $`x`$; and the vacuum is at $`x=0`$. The fundamental property of LF Hamiltonians, that under a rescaling of the LF momentum, $`x\lambda x`$, the LF Hamiltonian scales like $`HH/\lambda `$, can be interpreted as a dilatation symmetry along the $`x`$-axis, if one thinks of the $`x`$-axis as a spacial axis. This symmetry holds on a classical level and it is brocken on a quantum level by anomali. As one approaches small $`x`$, interaction between partons gets stronger, contributing divergent matrix elements. A natural cutoff is provided by $`\delta x=\epsilon x`$, a minimal spacing between constituents, which plays the role of UV-regulator bbbSmall $`x`$ correspond to the large light-front energies.. As long as $`\delta x`$ is finite, i.e. as long as the density of partons on the $`x`$-axis is not infinite, one obtains finite matrix elements. Cutoff $`\delta x`$ breaks the dilatation symmetry along the $`x`$-axis and gauge symmetry, and sets up a scale for quantities of dimension of energy (mass) and higher powers in energy. In particular, one generates dynamically an effective gluon mass, which in turn will define string tension between quark and antiquark. Formation of the $`q\overline{q}`$ bound state through breaking an internal symmetry can be viewed analogously to the creation of Cooper pairs in superconducter.

In terms of effective theory, the Hamiltonian below the LF cutoff, $`H_\epsilon `$, which desribes high energy (UV) strongly correlated states, can be substituted by the v.e.v., since strongly coupled configurations are frozen. Hamiltonian above the cutoff $`H_{>\epsilon }`$ is treated by flow equations (renormalization group transformation). One way of looking at the physics behind this v.e.v. and mass generation in gluodynamics is that the composite field $`\varphi `$ which creates $`0^+`$ glueballs has a finite v.e.v. $`^\mathrm{?}`$, $`^\mathrm{?}`$.

The dilatation symmetry, discussed above, reflects some underlying scale invariance of the LF field theory formulated on $`x`$-axis. Introducing the cutoff, breaks this symmetry. Because physics should remain independent on the cutoff, one must be looking for a fixed point of the renormalization group. Therefore the right tool for studing such a system is the renormalization group, which is provided by the method of flow equations $`^\mathrm{?}`$.

In the method of flow equations an effective interaction between quarks arises through elimination of quark-gluon coupling. Procedure converges in the UV for large gluon momenta transfer $`\stackrel{}{q}^2`$, but in the IR one is not able to eliminate vanishing gluon momenta, because the similarity factor which governes the elimination does not decay for vanishing arguments. In physical terms, one can not integrate soft gluons in the same fashion. As soon as gluon aquires a dynamical mass, which vanishes at large gluon momenta and not equal zero only at small momenta, the similarity factor decays even for small momenta transfer, because the argument contains an effective gluon energy instead of only gluon mometum. Therefore we can eliminate by flow equations even soft gluon modes, that are responsible for the long-range part of an effective quark potential. In the IR this gives $`1/q^4`$ behavior for an effective $`q\overline{q}`$ interaction.

In the next section we solve flow equations for an effective gluon mass and quark-antiquark effective interaction, based on the QCD Hamiltonian in the light-front gauge.

## 2 Gluon mass gap and an effective quark interaction

Let $`Q`$ being a projection operator on a one-gluon state, and $`P`$ on a $`q\overline{q}`$ state, $`Q|\psi =|g`$ and $`P|\psi =|q\overline{q}`$. Flow equations for matrix elements of the QCD Hamiltonian between these states read

$`{\displaystyle \frac{dE_q(l)}{dl}}`$ $`=`$ $`{\displaystyle \underset{p}{}}{\displaystyle \frac{1}{E_q(l)E_p(l)}}{\displaystyle \frac{d}{dl}}\left(h_{qp}(l)h_{pq}(l)\right)`$ (1)
$`{\displaystyle \frac{dh_{pp^{}}(l)}{dl}}`$ $`=`$ $`{\displaystyle \underset{q}{}}\left({\displaystyle \frac{dh_{pq}(l)}{dl}}{\displaystyle \frac{1}{E_p(l)E_q(l)}}h_{qp^{}}(l)+h_{pq}(l){\displaystyle \frac{1}{E_p^{}(l)E_q(l)}}{\displaystyle \frac{dh_{qp^{}}(l)}{dl}}\right),`$
where $`p`$ ($`p^{}`$) runs through all states in the $`P`$-subspace, and $`q`$ in the $`Q`$-subspace. Flow equations for the quark-gluon coupling $`h_{pq}`$ (coupling between $`P`$ and $`Q`$ sectors) are

$`{\displaystyle \frac{dh_{pq}(l)}{dl}}`$ $`=`$ $`\left(E_p(l)E_q(l)\right)\eta _{pq}(l)`$
$`\eta _{pq}(l)`$ $`=`$ $`{\displaystyle \frac{h_{pq}(l)}{E_p(l)E_q(l)}}{\displaystyle \frac{d}{dl}}\left(\mathrm{ln}f(z_{pq}(l))\right)`$
$`z_{pq}(l)`$ $`=`$ $`l\left(E_p(l)E_q(l)\right)^2,`$ (2)
where $`f`$ is the similarity function with properties $`f(0)=1`$ and $`f(z\mathrm{})=0`$. From Eq. (2), elimination of quark-gluon coupling

$`h_{pq}(l)=h_{pq}(0)f(z_{pq}(l)),`$ (3)
is governed by the similarity function $`f(l(E_p(l)E_q(l))^2)`$, which vanishes for the matrix elements with energy differences exceeding the cutoff $`\lambda `$, $`|E_p(l)E_q(l)|1/\sqrt{l}=\lambda `$, $`h_{pq}(0)`$ is the initial value, and $`E_q(l)`$ ($`E_p(l)`$) is a solution of the flow (renormalization group) equation, Eq. (1). As long as the argument of the similarity function is not equal zero, one can eliminate quark-gluon coupling and solve flow equations for an effective gluon energy and effective $`q\overline{q}`$-interaction, Eq. (1). However, if one neglects the cutoff dependence of energies, $`E_p(l)E_p(0)`$, the argument of the similarity function may become zero in the degenerate case, $`E_p(0)=E_q(0)`$, and then the energy denominator blows up in effective interactions in Eq. (1).

In physical terms, the leading behavior of the argument in the quark-gluon similarity function is given by an effective gluon energy transfer between quark and antiquark, $`\left(E_p(l)E_q(l)\right)E_{gluon}^{eff}(l)`$, which reduces at the initial value of flow parameter, $`l=0`$, to a gluon momentum transfer, $`E_{gluon}^{eff}(0)\stackrel{}{q}^2`$. Neglecting the cutoff dependence of the effective gluon energy, one eliminates quark-gluon coupling only at large gluon momenta, and fails at small gluon momenta, because the argument of the similarity function tends to zero. For small gluon momenta the dependence of gluon effective energy on the cutoff becomes important, given by a dynamically generated gluon mass, that makes elimination of quark-gluon coupling possible even for small gluon momenta. We therefore solve flow equations, Eq. (1), consistently for an effective gluon mass and an effective $`q\overline{q}`$-interaction below. We show, that elimination of quark-gluon coupling for high and low gluon momenta contribute correspondingly to the UV and IR parts of quark-antiquark potential, that behave in momentum space like $`1/\stackrel{}{q}^2`$ and $`1/\stackrel{}{q}^4`$, respectively, where $`\stackrel{}{q}`$ is a gluon momentum transfer.

### 2.1 Gluon gap equation

In the light-front frame flow equation for an effective gluon mass (the first equation in Eq. (1)), with connection to the light-front energy $`q^{}=(q_{}^2+\mu ^2(\lambda ))/q^+`$, reads

$`{\displaystyle \frac{d\mu ^2(\lambda )}{d\lambda }}`$ $`=`$ $`2T_fN_f{\displaystyle _0^1}{\displaystyle \frac{dx}{x(1x)}}{\displaystyle _0^{\mathrm{}}}{\displaystyle \frac{d^2k_{}}{16\pi ^3}}g_q^2(\lambda ){\displaystyle \frac{1}{Q_2^2(\lambda )}}{\displaystyle \frac{df^2(Q_2^2(\lambda );\lambda )}{d\lambda }}`$ (4)
$`\times `$ $`\left({\displaystyle \frac{k_{}^2+m^2}{x(1x)}}2k_{}^2\right)`$
$`+`$ $`2C_a{\displaystyle _0^1}{\displaystyle \frac{dx}{x(1x)}}{\displaystyle _0^{\mathrm{}}}{\displaystyle \frac{d^2k_{}}{16\pi ^3}}g_g^2(\lambda ){\displaystyle \frac{1}{Q_1^2(\lambda )}}{\displaystyle \frac{df^2(Q_1^2(\lambda );\lambda )}{d\lambda }}`$
$`\times `$ $`\left(k_{}^2(1+{\displaystyle \frac{1}{x^2}}+{\displaystyle \frac{1}{(1x)^2}})\right),`$
with

$`Q_1^2(\lambda )`$ $`=`$ $`{\displaystyle \frac{k_{}^2+\mu ^2(\lambda )}{x(1x)}}\mu ^2(\lambda ),Q_2^2(\lambda )={\displaystyle \frac{k_{}^2+m^2}{x(1x)}}\mu ^2(\lambda ),`$ (5)
where gluon couples to the quark-anti-quark pairs, with the quark-gluon coupling $`g_q(\lambda )`$, and pairs of gluons, with the three-gluon coupling $`g_g(\lambda )`$; for finite $`\lambda `$ these couplings are different from each other functions of momenta; connection between flow parameter and the cutoff, $`l=1/\lambda ^2`$, is used; similarity function $`f`$ plays the role of UV regulator in loop integrals; current quark mass is $`m`$; for $`SU(N_c)`$ $`T_f\delta _{ab}=\mathrm{Tr}(T^aT^b)=\frac{1}{2}\delta _{ab}`$ and the adjoint Casimir is $`C_a\delta _{ab}=f^{acd}f^{bcd}=N_c\delta _{ab}`$, $`N_c`$ is the number of colors (i.e., $`N_c=3`$).

Coupled system of equations for the cutoff dependent, effective gluon and quark masses was first derived in $`^\mathrm{?}`$. Here, neglecting the cutoff dependence of quark mass this system is reduced to Eq. (4).

Generally, it is difficult to solve Eq. (4). One of the reasons is that this equation contains (unknown) running coulings; therefore the cutoff dependence of couplings is neglected below. Second, initial condition for Eq. (4) is not known. The following renormalization condition to define the running gluon mass $`\mu (\lambda )`$ through the ’physical’ mass is assumed $`^\mathrm{?}`$: the effective Hamiltonian at the scale $`\lambda `$ has bosonic eigenstates with eigenvalues of the form $`q^{}=(q_{}^2+\stackrel{~}{\mu }^2)/q^+`$ $`^\mathrm{?}`$, i.e.

$`{\displaystyle \frac{q_{}^2+\stackrel{~}{\mu }^2}{q^+}}q^{},q`$ $`=`$ $`{\displaystyle \frac{q_{}^2+\mu ^2(\lambda )}{q^+}}q^{},q`$ (6)
$``$ $`{\displaystyle ^\lambda }𝑑\lambda ^{}{\displaystyle \underset{p}{}}\left(\eta _{q^{}p}(\lambda ^{})h_{pq}(\lambda ^{})h_{q^{}p}(\lambda ^{})\eta _{pq}(\lambda ^{})\right),`$
where $`\stackrel{~}{\mu }`$ denotes the ’physical’ gluon mass; the generator $`\eta _{pq}`$, given by Eq. (2), eliminates the quark-gluon (three-gluon) coupling $`h_{pq}`$; $`|q`$ denotes a single effective gluon state with momentum $`q^+`$ and $`q_{}`$, $`q^{}|q=16\pi ^3q^+\delta ^{(3)}(q^{}q)`$. Neglecting the cutoff dependence of the gluon mass on the r.h.s. of Eq. (4), the renormalization point is choisen at $`q^2=\stackrel{~}{\mu }^2`$ (for in- and out-going and intermediate state gluons). Third, severe divergences arise at small light-front momenta $`x`$ and must be regularized. For this purpose, as suggested by Zhang and Harindranath in $`^\mathrm{?}`$, the minimum cutoff $`u`$ for transverse momentum $`k_{}`$ is introduced. Therefore, flow equation, Eq. (4), is integrated in the finite limits, $`[u;\lambda ]`$. Explicitly, from the similarity function in gluon loop (the second term in Eq. (4)) one has $`u\stackrel{~}{Q}_1^2=(k_{}^2+\stackrel{~}{\mu }^2)/x(1x)\stackrel{~}{\mu }^2\lambda `$, that restricts the transverse momentum to $`k_{min}=(u^2+\stackrel{~}{\mu }^2)x(1x)\stackrel{~}{\mu }^2k_{}k_{max}=(\lambda ^2+\stackrel{~}{\mu }^2)x(1x)\stackrel{~}{\mu }^2`$, and the light-front momentum to $`\stackrel{~}{\mu }^2/(u^2+\stackrel{~}{\mu }^2)x1\stackrel{~}{\mu }^2/(u^2+\stackrel{~}{\mu }^2)`$. In the same way one finds restrictions on mometa in the quark-gluon loop (the first term in Eq. (4)). Integrating Eq. (4) with all constraints discussed above and assuming the condition of Eq. (6), gives $`^\mathrm{?}`$

$`\mu ^2(\lambda )`$ $`=`$ $`\stackrel{~}{\mu }^2+\delta \mu _{PT}^2(\lambda )+\delta \mu _{dyn}^2(\lambda ,\stackrel{~}{\mu },u).`$ (7)
where the perturbative term

$`\delta \mu _{PT}^2(\lambda )={\displaystyle \frac{g^2}{4\pi ^2}}\lambda ^2\left(C_a(\mathrm{ln}{\displaystyle \frac{u^2}{\stackrel{~}{\mu }^2}}{\displaystyle \frac{11}{12}})+T_fN_f{\displaystyle \frac{1}{3}}\right).`$ (8)
reproduces the result of the LF perturbation theory, with renormalization point $`q^2=0`$, $`^\mathrm{?}`$. When renormalization is performed in the second oder, the perturbative mass correction is absorbed by the mass counterterm, $`m_{CT}^2=\delta \mu _{PT}^2(\lambda =\mathrm{\Lambda }\mathrm{})`$. After the perturbative renormalization is completed, the dynamical mass

$`\mu _{dyn}^2(\lambda )`$ $`=`$ $`\stackrel{~}{\mu }^2+\delta \mu _{dyn}^2(\lambda ,\stackrel{~}{\mu },u)=\stackrel{~}{\mu }^2+\sigma (\stackrel{~}{\mu },u)\mathrm{ln}{\displaystyle \frac{\lambda ^2}{\stackrel{~}{\mu }^2}}`$
$`\sigma (\stackrel{~}{\mu },u)`$ $`=`$ $`{\displaystyle \frac{g^2}{4\pi ^2}}\stackrel{~}{\mu }^2\left(C_a({\displaystyle \frac{u^2}{\stackrel{~}{\mu }^2}}+\mathrm{ln}{\displaystyle \frac{u^2}{\stackrel{~}{\mu }^2}}{\displaystyle \frac{5}{12}})+T_fN_f({\displaystyle \frac{1}{3}}+{\displaystyle \frac{m^2}{\stackrel{~}{\mu }^2}})\right),`$ (9)
is left. Note, that the gluon mass $`\mu _{dyn}`$ is generated dynamically by flow equations. In the limit $`\stackrel{~}{\mu }0`$, one has $`\sigma =lim_{\stackrel{~}{\mu }0}\sigma (\stackrel{~}{\mu },u)=u^2g^2C_a/2\pi ^2`$, and, as shown below, $`\sigma `$ plays the role of the string tension between quark and antiquark. Scale $`u`$ is introduced by regulating the divergences at small LF momenta, $`x0`$, (small LF $`x`$ regularization). The value $`u`$ sets up a scale for the dynamical gluon mass and the string tension.

### 2.2 Effective $`q\overline{q}`$-interaction

An effective interaction, generated between quark and antiquark by flow equations in the LF frame (the second equation in Eq. (1)), reads

$`V_{q\overline{q}}=Const\gamma ^\mu \gamma ^\nu \underset{\stackrel{~}{\mu }0}{lim}B_{\mu \nu },`$ (10)
where instead of the strong coupling $`\alpha _s`$ some constant, $`Const`$, is introduced; $`\gamma ^\mu \gamma ^\nu `$ is a current-current term in the exchange channel. At the end we set the gluon mass renormalization point (’physical’ gluon mass) to zero, $`\stackrel{~}{\mu }0`$; thus eliminating the sensitivity to the renormalization point. The kernel includes dynamical interaction and instantaneous exchange with nonzero dynamical gluon mass, given by

$`B_{\mu \nu }=g_{\mu \nu }\left(I_1+I_2\right)+\eta _\mu \eta _\nu {\displaystyle \frac{\delta Q^2}{q^{+2}}}\left(I_1I_2\right),`$ (11)
where $`g_{\mu \nu }`$ is the LF metric tensor, and the LF vector $`\eta _\mu `$ is defined as $`\eta k=k^+`$. The cutoff dependence of four-momentum transfers along quark and antiquark lines is accumulated in the factor

$`I_1={\displaystyle _0^{\mathrm{}}}𝑑\lambda {\displaystyle \frac{1}{Q_1^2(\lambda )}}{\displaystyle \frac{df(Q_1^2(\lambda );\lambda )}{d\lambda }}f(Q_2^2(\lambda );\lambda ),`$ (12)
with mometum transfers defined in the light-front frame as

$`Q_1^2(\lambda )`$ $`=`$ $`Q_1^2+\mu _{dyn}^2(\lambda ),Q_1^2={\displaystyle \frac{(x^{}k_{}xk_{}^{})^2+m^2(xx^{})^2}{xx^{}}}`$
$`Q_2^2(\lambda )`$ $`=`$ $`Q_2^2+\mu _{dyn}^2(\lambda ),Q_2^2=Q_1^2|_{x(1x);x^{}(1x^{})},`$ (13)
and the dynamical gluon mass $`\mu _{dyn}`$ is given by Eq. (9); also the following momenta are used $`Q^2=(Q_1^2+Q_2^2)/2`$ and $`\delta Q^2=(Q_1^2Q_2^2)/2`$. Calculating Eq. (10) with an explicit form of similarity function, gives in the leading order $`\delta Q^2Q^2`$ the following $`q\overline{q}`$-interaction $`^\mathrm{?}`$

$`V_{q\overline{q}}=\gamma ^\mu \gamma _\mu \left(C_f\alpha _s{\displaystyle \frac{4\pi }{Q^2}}+\sigma {\displaystyle \frac{8\pi }{Q^4}}\right),`$ (14)
where to this order $`Q^2`$ reduces in the instant frame to the square of gluon momentum transfer, $`Q^2\stackrel{}{q}^2`$, with the gluon momentum $`\stackrel{}{q}=(q_z,q_{})`$. In order to reproduce the standard Coulomb and linear confining potentials, the correct prefactors are restored, using the freedom to fit $`Const`$ and $`\sigma `$ terms. Confining term in Eq. (14), with singular behavior like $`1/\stackrel{}{q}^4`$, arise from the elimination of the quark-gluon coupling at small gluon momenta, that is governed by the cutoff dependent, dynamical gluon mass.

## Acknowledgments

The author would like to thank the organizers of the workshop for hospitality and support. The author is thanful to Stan Brodsky and Lev Lipatov for helpful discussions. This work was supported by DOE grants DE-FG02-96ER40944, DE-FG02-97ER41048 and DE-FG02-96ER40947.

## References