text:
# Layer by layer epitaxy in limited mobility nonequilibrium models of surface growth

## I Introduction

Thin film growth, under solid-on-solid epitaxial conditions, from the vacuum vapor deposition of an atomic or a molecular beam (the so-called molecular beam epitaxy or MBE) is an important technological process used extensively to produce high quality thin films with smooth and flat surfaces and interfaces. It is also a growth process of considerable fundamental significance  in the statistical mechanics of nonequilibrium phenomena because MBE (at least in its ideal form , with no evaporation and vacancy or overhang formation at the growth front) in principle represents  a novel universality class of nonlinear surface growth outside the generic Kardar-Parisi-Zhang universality . A great deal of attention has therefore focused over the last ten years on the statistical properties of kinetically rough (and, in principle, generically scale invariant) surface growth in low temperature (room temperature or below) MBE, following the suggestions  of the possible importance of ideal MBE in defining novel growth universality classes in kinetic surface roughening . We note that the conserved surface current nature  of ideal MBE growth (i.e. solid-on-solid growth with no evaporation and vacancy or overhang formation) rules out a KPZ description of its growth dynamics.

It is interesting, perhaps even ironic, that MBE growth has played such a central role  in kinetic surface roughening phenomena because the primary materials science impetus for MBE growth is obviously to avoid kinetic roughening as completely as possible in order to produce smooth and flat thin films of high surface quality with minimal amount of surface roughness. From the materials science perspective of producing high quality smooth (i.e. manifestly non-rough) thin films, therefore, MBE is typically carried out at elevated temperatures ($``$ 500-1000 K), where fast surface diffusion enables one to produce smooth thin films with very little surface roughness. Smooth MBE growth, as opposed to kinetically rough (low temperature) growth, is characterized by layer-by-layer growth oscillations , where each layer of the growing thin film (on a singular high symmetry substrate) essentially fills up completely before the next layer deposition begins (on the other hand, kinetically rough growth is, by definition, multilayer as many layers are partially filled at the growth front producing increasing surface roughness with no layer-by-layer oscillations), consequently the surface morphology and associated properties oscillate as growth proceeds. In the ideal layer-by-layer growth mode, therefore, the interface width ($`W`$, the root mean square fluctuation in the interface height) of the growing film oscillates (nominally between 0 and 1, as measured in lattice units, indicating an empty or a filled layer) whereas in the kinetically rough growth mode $`W`$ increases monotonically as a power law in the average film thickness ($`h`$). Layer-by-layer epitaxy is, however, an initial transient growth regime which eventually crosses over to the asymptotic kinetically rough growth regime as the shot noise intrinsic in the incident deposition beam fluctuations always wins out to damp the layer-by-layer oscillations, and at long enough time scales (and for large enough lateral system sizes) statistically scale invariant kinetically rough growth would always emerge. (In fact, the noise associated with the stochastic diffusion process also contributes to kinetic roughening, but the shot noise associated with the incident beam fluctuations is the most relevant roughening mechanism.) This is also the experimental observation : layer-by-layer growth oscillations, as studied for example through RHEED intensity oscillations monitoring the dynamical surface evolution , eventually always damp out as the stochastic depositopn shot noise associated with incident particle beam fluctuations leads to kinetically rough multilayer growth after some characteristic time $`t_c`$. The damping time $`t_c`$, beyond which layer by layer growth dies out, depends on the growth temperature (which controls the surface diffusion rate), and is in general larger for higher temperatures because longer diffusion lengths at higher temperatures enhance layer by layer growth. (There is actually a complication, arising from the unaviodable vicinality in the starting substrate which can never really be precisely a high symmetry singular plane in real growth where layer by layer oscillations tend to disappear at both high and low growth temperatures — the low temperature behavior is from the multilayer kinetically rough growth as discussed above, but the high temperature disappearance arises from the so-called step flow growth mode which is caused by the very fast surface diffusion of deposited atoms at high temperatures leading to their moving directly to step edges, which must be present in any real substrate due to vicinality or miscut, without any layer by layer growth oscillations; we neglect considerations of such a step flow growth mode in this paper assuming all growth to be occuring strictly on singular high symmetry substrates.)

The ideal MBE growth (on flat singular high symmetry substrates) can be thought of as composed of two regimes — an early time ($`t<t_c`$) transient (fast diffusion driven) regime of layer by layer growth followed by the asymptotic ($`t>t_c`$) kinetically rough (deposition beam fluctuation driven) growth regime (with no layer by layer oscillations) characterized by power law evolution  in surface roughness. At low enough temperatures, when surface diffusion is extremely slow, $`t_c`$ could be less than the time it takes to grow one monolayer of deposit on the average, and in that situation layer by layer transient growth regime is invisible with the kinetically rough growth regime being dominant essentially right from the beginning. Conversely, sufficiently high temperature growth on a small substrate could continue in the layer by layer mode for a very long time although some damping of the growth oscillations is inevitable with time as distant spatial regions on the substrate must lose coherence due to the inherent shot noise fluctuations associated with the discrete deposition process in the incident beam. Thus layer by layer growth us purely a “finite size” (both spatially and temporally) transient phenomenon — if the substrate is made sufficiently large and/or if one waits for sufficiently long growth time, layer by layer epitaxy must necessarily cross over to kinetically rough growth. It is important to emphasize, however, that the surface diffusion length is typically an exponentially activated function of growth temperature, and therefore a small change in temperature could cause a sharp large change in the growth mode (layer by layer to rough and vice versa depending on whether the growth temperature is decreased or increased) for a given substrate, leading to the empirical concept  of an epitaxial growth temperature $`T_c`$ with growth being layer by layer (rough) for $`T>T_c`$ ($`T<T_c`$) — clearly $`T_c`$ is a loosely defined concept because it must be a weak (sub-logarithmic) function of the effective substrate (or the terrace) size for a given material . In general, “good” MBE growth aiming toward producing high quality smooth epitaxial thin flims is carried out at the highest possible growth temperature (within the constraint that evaporation or desorption from the growth front should be negligible so that the growth temperature cannot be arbitrarily high) so as to make atomic mobility at the growth front to be very high leading to large “surface diffusion length” $`l`$. Here $`l`$ is taken to be the linear size over which the surface is smooth due to atomic diffusion. Assuming the deposition process to be a random Poisson process it is then easy to see that the typical surface roughness over terraces of size $`l`$ would only grow as $`\sqrt{h}/l`$, where $`h`$ is the thickness of the grown film (and we measure all lengths in lattice units). Thus for large $`l`$, one would have to grow a very thick film of thickness $`l^2`$ before the surface roughness reaches even one monolayer fluctuation. One can, therefore, grow MBE thin films of very high smoothness and quality, without worrying at all about the kinetic surface roughening by properly adjusting the growth temperature  to make $`l`$ large.

Layer by layer MBE growth has been extensively studied  in the literature using computer simulations of MBE growth through the stochastic (or kinetic) Monte Carlo simulations, with the atomistic diffusion at the growth front assumed to be controlled by stochastic activated hopping process with the hopping rate determined by a local coordination number dependent activated Arrhenius hopping. Such activated diffusion Arrhenius hopping simulations (sometimes also referred to as “full diffusion” simulations to differentiate them from “limited mobility” growth models which are our main interest in this paper) involve continuous (possible) hopping of all surface atoms according to their local bonding configurations (which determine the activation energy for the hopping process). Such full diffusion simulations are obviously not well designed to study the kinetic roughening universality class of MBE growth because they are extremely time consuming and cannot really be carried out for large systems (particularly in the physically relevant 2+1 dimensions) for long times, an essential requirement for ascertaining the asymptotic universality class of a growth model. Although there are some notable exceptions , the full diffusion Arrhenius activated kinetic Monte Carlo simulations of MBE growth have not been used with particular success for understanding statistical scale invariance properties of kinetic surface roughening. Instead, important insights into the MBE universality class of kinetic surface roughening have come primarily from nonequilibrium limited mobility growth models — mainly the so-called Wolf-Villain  (WV) and the Das Sarma-Tamborenea  (DT) model — which were introduced specifically for the elucidation of the MBE growth universality.

In this paper we study the DT and the WV model (we emphasize that WV and DT models, in spite of their close similarity in growth rules, belong  to different asymptotic universality classes in both 1+1 and 2+1 dimensions although their pre-asymptotic scaling behavior is very similar which has led to considerable confusion in the literature) in the complementary layer by layer growth regime rather than the kinetically rough growth regime which motivated the introduction of these models. We mention in this context that some 1+1 dimensional studies of WV and DT models in the layer by layer growth regime have recently been reported in the literature . Our results, where applicable, agree with these earlier works , but our focus in this paper is 2+1 dimensional growth and the effect of long surface diffusion length in 1+1 dimensional growth, neither one of which has earlier been studied.

In limited mobility growth models (the models and growth rules used in this paper are described in section II of the paper — see, for example, Fig. 1), in sharp contrast to full diffusion MBE growth simulations, the goal is to suppress crossover and transient effects as much as possible (so as to efficiently reach the asymptotic kinetic surface roughening regime) and as such only the most recently deposited atom is allowed to diffuse or relax instantaneously to the appropriate incorporation site following the mobility rules of the specific model. This allows suppression of crossover effects invariably present in the full diffusion simulations arising from many different diffusion rates corresponding to many different possibilities for local bonding configurations — the only time scale in the limited mobility growth models being the deposition rate, which defines the time unit for the problem. From now on we take the time unit (sometimes referred to as a “second”) as the time to deposit one monolayer on the average. Thus the growth time in this paper also defines $`h`$ — the average thickness of the deposited film measured in units of monolayers or the lattice constant, which we take to be the unit of length throughout. (With no loss of generality we choose the lattice constant to be the same along the substrate and the growth directions.)

The limited mobility growth models  are by construction strongly dominated by the deposition shot noise because the goal is to study the scale invariant kinetic surface roughening behavior. This is particularly true in the original versions of the growth model where the surface diffusion length is choosen to be unity, $`l=1`$, i.e. the deposited atoms are allowed to move only to the nearest neighbor incorporation sites around the deposition site. The original DT and WV models therefore did not exhibit, by design, any layer by layer growth oscillations since the smoothing or the healing distance ($`l`$) was just one lattice unit. In order to manifest layer by layer epitaxy in limited mobility growth models one must therefore suppress the shot noise associated with the incident beam fluctuations.

In this paper we accomplish the noise suppression by two alternative techniques : The ‘multiple hit’ noise reduction technique  and the ‘long surface diffusion length ($`l>1`$)’ noise reduction technique . These techniques, described in section II of the paper, give rise to layer by layer growth (as monitored by an oscillatory surface roughness, i.e. $`W(t)`$ showing oscillations as a function of growth time $`t`$) in the limited mobility growth models as described in section III of the paper.

The rest of this paper is organized as follows. In section II we describe the limited mobility growth models and the noise reduction technique(s) employed by us. We also provide some theoretical background for our analysis of the simulation results. In section III we present and describe our numerical simulation results for layer by layer epitaxial growth in DT and WV models. We also discuss in section III various (approximate) scaling properties of our simulated layer by layer epitaxial growth results. We conclude in section IV with a general discussion of our results making connections with some of the existing results in the literature and pointing out possible future directions as well.

## II Models, Theory, and Background

The DT  and WV  models used in our simulations are shown in Fig. 1. We carry out growth simulations in both 1+1 dimensions and 2+1 dimensions (on 100 high symmetry substrates). Particles are dropped, obeying solid on solid constraint, sequentially (one by one) and randomly with an average rate of $`N`$ per second, where $`N=L`$ for d=1+1 and $`N=L^2`$ for d=2+1, on a substrate of lateral linear size of $`L`$ lattice units. We measure length in lattice units and time in inverse deposition rate (i.e. in units of monolayer filling since one average monolayer of deposition occures every “second”). Each deposited atom is allowed to “diffuse” instantaneously to its incorporation site following the mobility rules of the specific model. The diffusion rules in the DT model are that a deposited atom can move only if it has no lateral nearest neighbor in the same layer (if it does, then the atom is incorporated at the deposition site) — if the deposition site has no lateral nearest neighbor then the incident atom may move instantaneously to a neighboring empty site (within a lateral diffusion length of $`l`$, where $`l=1`$ in the original model and in most existing simulations) provided the final incorporation site has a higher lateral coordination number (i.e. one or higher) than the deposition site. If several neighboring sites satisfy the diffusion rule then the atom will move randomly to any one of them with equal probability. The rules for the WV model are superficially similar to the DT rules : In the WV model all deposited atoms (and not just the ones with no lateral bonds) can, in principle, move provided they can increase their local lateral bonding and the deposited atom always moves to the site with the maximum local bonding environment. In both models, the deposited atom is incorporated at the deposition site if it cannot satisfy the diffusion rules (i.e. no sites with higher coordination available) within the diffusion length.

Both DT and WV models have been extensively studied in the literature (mostly within nearest neighbor, $`l=1`$, diffusion rules) in the context of their kinetic surface roughening universality classes. Recently, layer by layer epitaxy in the WV  and the DT  model have been investigated in 1+1 dimensions using the multiple hit noise reduction technique. The very first simulational observation of layer by layer growth in a limited mobility growth model was reported in the DT model in Ref.  where it was studied in 1+1 dimensions using a long ($`l>1`$) surface diffusion length, but no details were investigated. We emphasize that the usual $`l=1`$ limited mobility growth model does not exhibit any layer by layer epitaxy by definition, and manifests kinetic roughening right from the beginning since for $`l=1`$ the layer by layer epitaxy regime is restricted to less than one monolayer coverage, i.e. in the standard limited mobility growth models  the layer by layer epitaxy regime does not exist.

To obtain a layer by layer growth regime in the DT and the WV model we use two distinct techniques to suppress noise and enhance diffusion, which enable our growth simulations to manifest strong layer by layer growth oscillations before crossing over to the kinetic surface roughening regime with pure multilayer growth. One technique, referred to as the noise reduction technique, has also been used by others  to produce layer by layer growth in various models, e.g. Eden model , single step model , and the WV and the DT model in 1+1 dimensions . We have earlier used this technique  to suppress corrections to scaling in the asymptotic kinetically rough surface growth regime in the DT and the WV model in order to accurately determine the dynamic scaling exponents and the associated growth universality class. In the noise reduction technique, characterized by an integer number $`m`$ (the usual growth model without any noise reduction is an $`m=1`$ model), a counter is put on each surface site and each discrete deposition event on a site advances the counter by unity. A deposition event at a particular site is accepted only when the counter reaches a pre-determined number $`m>1`$. Thus, this technique is the multiple hit noise reduction technique since $`m(>1)`$ deposition hits on a site are needed for a true deposition. After $`m`$ hits on a site (i.e. after the acceptance of a deposition event) the counter at that particular site is set back to zero, and the whole multiple hit process begins all over again. The multiple hit noise reduction technique is a coarse-graining procedure which suppresses the deposition shot noise, and the noise reduction is enhanced for larger values of $`m`$. The second technique applied by us for obtaining layer by layer growth oscillations in limited mobility growth models is to use long surface diffusion lengths ($`l>1`$) in the growth simulations. Obviously long diffusion lengths enhance the layer by layer growth regime, and in particular, for $`l>L`$ i.e. the surface diffusion length exceeding the system size, the layer by layer growth may persist essentially indefinitely since each deposited atom may always be able to seek out a desired epitaxial site for incorporation. In some sense the multiple hit noise reduction parameter $`m`$ is equivalent to the dimensionless diffusion length parameter $`l/L`$ because large values of both tend to enhance the layer by layer growth regime. In section III, where we present our simulation results, we will see the precise nature of this correspondence between $`m`$ and $`l/L`$ in our two methods of obtaining layer by layer epitaxy in the DT and the WV model.

The central quantity of interest in layer by layer epitaxial growth is the characteristic time $`t_c`$ at which layer by layer growth dies out, i.e. for the deposited average film thickness larger than $`t_c`$, measured in lattice units or in monolayers, there are no discernible layer by layer growth oscillations. It has been found in earlier numerical simulations of layer by layer growth in a variety of contexts that $`t_c`$ obeys an approximate scaling relation with the coarse graining parameter $`m`$ (or $`l/L`$ as the case may be), and we are interested in investigating whether the following scaling relations hold in the limited mobility growth models

$$t_c\{\begin{array}{cc}m^\mu \hfill & \text{for the noise reduction method}\hfill \\ (l/L)^\delta \hfill & \text{for the long diffusion length method.}\hfill \end{array}$$
(1)
If such scaling relations do hold in our simulations we are interested in obtaining the relationship, if any, between the exponents $`\delta `$ and $`\mu `$.

It is in fact fairly straightforward to obtain a relationship between the noise reduction parameter $`m`$ and the surface diffusion length parameter $`l/L`$ using only dimensional arguments. In particular, we note that for a $`d^{}`$-dimensional substrate ($`d^{}=2`$ for real surfaces and $`d^{}=1`$ for the 1+1 -dimensional growth) one obtains a simple relationship between $`\mu `$ and $`\delta `$ by noting that a surface diffusion length $`l`$ corresponds to $`ml^d^{}`$ since there are $`l^d^{}`$ available sites for a particular deposited particle to incorporate at. This immediately leads to

$$\delta =\mu d^{}.$$
(2)
Within the limited accuracy of our growth simulations (see section III for the results), we find the scaling relation defined by Eq. (2) to be valid.

One secondary theoretical goal of our work is to investigate the extent to which layer by layer growth oscillations obey scaling with respect to growth time (or equivalently, the average film thickness). Using Eq. (1) as the theoretical ansatz one can ask whether the surface roughness $`W`$, defined as the ensemble averaged (over many growth simulations) root mean square fluctuation in the interface height, which is oscillatory (and $`W<1`$ implying little roughness) in the layer by layer growth regime is a general scaling function of the layer by layer growth parameter $`m`$ or $`l/L`$ through a dependence of the form :

$$W(t)f_m(t/m^\mu )\text{ or }f_l(t/(l/L)^\delta ).$$
(3)
If a scaling form such as Eq. (3) holds in the layer by layer growth regime, then Eq. (1) for $`t_c`$ trivially follows from it — $`t_c`$ being the value of time where layer by layer oscillations cease to exist. We could go further in our scaling analyses and ask whether the scaling defined by Eq. (3) continues to hold (perhaps approximately) well beyond ($`t>t_c`$) the layer by layer growth regime establishing an approximate scaling relationship between the layer by layer growth regime ($`t<t_c`$) and the kinetically rough growth regime ($`t>t_c`$). Our results presented in section III indicate that such an approximate scaling relation does indeed exist between the layer by layer growth regime and the kinetically rough growth regime.

Finally, we note that there have been recent attempts  at developing a theory for layer by layer growth oscillations starting from continuum growth equations underlying the coarse grained long wavelength behavior of MBE growth. A simple dimensional argument , later followed up  by a renormalization group approach, leads to the following results for MBE growth :

$$\mu =\{\begin{array}{cc}4/3\hfill & \text{for }d^{}=1\hfill \\ 2\hfill & \text{for }d^{}=2\text{.}\hfill \end{array}$$
(4)
The exponents defined by Eq. (4) correspond to the so-called  conserved fourth order nonlinear growth equation or the nonlinear MBE growth equation (which is sometimes also referred to as the conserved KPZ equation with nonconserved noise or the Lai-Das Sarma-Villain , LDV, equation) which is given by

$$\frac{h}{t}=\nu _4^4h+\lambda _{22}^2(\mathbf{}h)^2+\eta ,$$
(5)
where $`h(𝐱,𝐭)`$ is the dynamical height fluctuation variable relative to the average interface $`h`$ at the substrate site $`𝐱`$ (with $`𝐱`$ is the lateral substrate coordinate), $`\mathbf{}\frac{}{𝐱}`$ is the gradient operator along the surface, $`\eta `$ is the deposition shot noise (which causes the kinetic surface roughening), and $`\nu _4`$, $`\lambda _{22}`$ are coefficients which in general depend on surface diffusion rate, deposition rate, etc. Since the continuum description of the DT and the WV model are actually quite complex , and are likely to be different in different dimensions , it is by no means clear that the exponents (given in Eq. (4)) corresponding to Eq. (5) apply without any qualifications to the DT and the WV model (as will be discussed in section III where we present our numerical simulations). We therefore also provide below the exponent $`\mu `$ for the linear second order Edwards-Wilkinson (EW) growth equation  which applies to the limited mobility Family (F) growth model  and may also have significant relevance to the DT and the WV model  :

$$\mu =\{\begin{array}{cc}2\hfill & \text{for }d^{}=1\hfill \\ \mathrm{}\hfill & \text{for }d^{}=2\text{.}\hfill \end{array}$$
(6)
The EW equation, whose layer by layer growth exponent in 1+1 and 2+1 dimensions is given in Eq. (6), is the following :

$$\frac{h}{t}=\nu _2^2h+\eta .$$
(7)
We note here that both sets of exponent values given by Eq. (4) and (6) will be relevant in our discussion of our simulation results to be presented in the next section.

## III Results and Discussions

We now present our 1+1 (i.e. $`d^{}=1`$) and 2+1 (i.e. $`d^{}=2`$) dimensional layer by layer growth simulation results for the discrete limited mobility DT and WV models in Figs. 2-6. In Figs. 2-4 we present 1+1 dimensional simulation results whereas Figs. 5 and 6 give 2+1 dimensional simulational results for the two growth models. In each figure (to be described below) the panel (a) gives the simulated $`W(t)`$ as a function of growth time $`t`$ for various values of the noise reduction parameter $`m`$ or the surface diffusion length $`l`$ with the layer by layer oscillations manifestly obvious for larger values of $`m`$ and $`l`$. (The original DT and WV models correspond to the simulation with $`m=1`$ and $`l=1`$, which has no layer by layer oscillations by construction.) In panel (b) of each figure we demonstrate our best computed scaling collapse of the $`Wt`$ plots (shown in panel (a)) for various values of $`m`$ or $`l`$ with a suitable scaling of time $`t`$ to $`t/m^\mu `$ or $`t/(l/L)^\delta `$ as the case may be. For each scaling collapse in (b) we try various different values of the exponent $`\mu `$ or $`\delta `$ ro obtain the best statistical scaling in the simulated data. The system size used in each simulation is indicated in the corresponding figures and the captions. Here we should emphasize that the results shown in this paper represent only a typical fraction of our extensive DT and WV layer by layer growth simulations. The representative results presented here are of course in complete agreement with the full set of our simulation data, and our conclusion is based on a very large set of simulation results and not just on the results presented in this paper.

In Fig. 2 we show our 1+1 dimensional DT layer by layer growth simulation results for finite surface diffusion length $`l`$ ($`m=1`$, $`l>1`$) for $`l=10,\mathrm{\hspace{0.17em}20},\mathrm{\hspace{0.17em}30},\mathrm{\hspace{0.17em}40},\mathrm{\hspace{0.17em}50}`$. The layer by layer oscillations are visually obvious in Fig. 2(a) — the damping time $`t_c`$ increases from $`t_c10`$ for $`l=10`$ to roughly $`t_c50`$ for $`l=50`$. In general, the magnitude of the oscillations decays exponentially with increasing time, and for $`t>t_c`$ we can only discern the power law increase of $`Wt^\beta `$ where $`\beta `$ is the growth exponent in the model. In Fig. 2(b) we show our scaling collapse of the $`W(t)`$ data from Fig. 2(a), leading to the exponent value $`\delta 1.5`$. Thus, $`t_c(l/L)^{1.5}`$ in $`d^{}=1`$ DT model. We note that the scaling collapse in the kinetically rough growth regime (i.e. for $`t>t_c`$) is not excellent, but it is remarkable that the exponent $`\delta `$ which is meaningfully defined only in the layer by layer growth regime continues to provide an approximate reasonable description of the kinetically rough growth regime.

In Fig. 3 we show our 1+1 dimensional DT model layer by layer growth oscillations (Fig. 3(a)) for the multiple hit noise reduction technique ($`l=1`$, $`m>1`$) for different values of the noise reduction factor $`m`$. In Fig. 3(b) we show the scaling collapse of the $`Wt`$ data in Fig. 3(a) for various $`m`$ values. The scaling is excellent with an exponent $`\mu =1.5`$. Thus, $`t_cm^{1.5}`$ in $`d^{}=1`$ DT model.

In Fig. 4 we depict our noise reduced layer by layer growth oscillations (Fig. 4(a)) and the corresponding scaling collapse of the data for various values of $`m`$ (with $`l=1`$) in the 1+1 dimensional WV model simulations. Again, the scaling exponent $`\mu `$ is found to be $`\mu =1.5`$ for the best scaling collapse, indicating $`t_cm^{1.5}`$ in both DT and WV noise reduced models in $`d^{}=1`$. The finding of the apparent same exponent value $`\mu =1.5`$ in both DT and WV models in 1+1 dimensional growth is consistent with the fact that the effective growth exponent $`\beta `$ (obtained by plotting $`\mathrm{ln}W`$ against $`lnt`$ in the simulated results) is almost identical in the two models : From the slope of the $`\mathrm{log}\mathrm{log}`$ plot in Fig. 2(b) we obtain $`\beta 0.338`$ for the $`d^{}=1`$ DT model whereas from the slope of the $`\mathrm{log}\mathrm{log}`$ plot in Fig. 3(b) we obtain $`\beta 0.339`$ for the $`d^{}=1`$ WV model. Thus within the effective time and length scales of our simulations the two models (DT and WV) have essentially the same effective dynamical universality class, which is consistent with the fact that they have the same effective exponent $`\mu =1.5`$ in $`d^{}=1`$. The fact that the asymptotic universality classes of the DT and WV models are different even in $`d^{}=1`$ dimension  does not seem to affect the effective values of $`\mu `$ we obtain in our simulations.

Before presenting our 2+1 dimensional simulation results in Fig. 5-7 we first discuss the exponent values $`\delta `$ and $`\mu `$, all of which have turned out to be approximately $`1.5`$ in the 1+1 dimensional DT and WV layer by layer epitaxial simulations. (We do not show here $`l>1`$, $`m=1`$ simulation results for the $`d^{}=1`$ WV model because they are very similar to those shown in Figs. 2-4 with the same exponent $`\mu 1.5`$.) First we note that the exponent value is very close to (but somewhat above) the theoretically “expected” exponent $`\mu =4/3=1.33`$ predicted in Refs.  assuming that the continuum growth equation for these discrete growth models is that given in Eq. (5). We also note that the expected relationship between $`\mu `$ and $`\delta `$, namely $`\delta =\mu d^{}`$ which becomes $`\delta =\mu `$ in $`d^{}=1`$, is obeyed by our simulation results. We also add that the measured exponent $`\mu 1.5`$ is consistent with other findings in the literature . The cause for our calculated $`\mu (1.5)`$ to be somewhat (by roughly $`10\%`$) higher than the theoretical value ($`\mu =4/3`$) is not very clear at this stage. We do not, however, believe this discrepancy to be particularly significant because of a number of reasons : (1) our scaling collapse are in fact not inconsistent with an exponent of $`1.33`$ although an exponent value of $`1.5`$ is definitely a statistically better fit for our scaling collapse; (2) it is quite conceivable that there are some systematic finite size and finite time corrections to scaling (which are known to be very important in DT and WV models); (3) finally, at least in the WV model which is definitely known  to asymptotically belong to the EW universality class, it is possible that our simulated exponent $`\mu 1.5`$ is showing some effects of the asymptotic universality class since the theoretically expected  $`\mu `$ for the linear EW equation (our Eq. (7)) is $`\mu =2`$ in $`d^{}=1`$ dimensions.

Our $`d^{}=2`$ dimensional noise reduced ($`m>1`$, $`l=1`$) results for the DT, the WV, and the F model are shown in Figs. 5-7 respectively. These 2+1 dimensional layer by layer growth results (as well as the results shown in Fig. 2) in limited mobility models are completely new and do not exist anywhere in the literature. We carried out our 2+1 layer by layer growth simulations using only the noise reduction technique since our 1+1 dimensional results (compare Figs. 2 and 3) indicate that the finite diffusion length and the multiple hit noise reduction techniques are essentially equivalent and the 2+1 dimensional simulations with $`l>1`$ are particularly cumbersome to carry out. Our $`d^{}=2`$ layer by layer epitaxy simulations seem to have produced a few surprising results as discussed below.

In Fig. 5 we present our $`d^{}=2`$ layer by layer growth simulations in the DT model using the noise reduction technique. The results are depicted in the same manner as in the $`d^{}=1`$ case shown in Fig. 3 — in particular, Fig. 5(a) shows the actual layer by layer growth oscillations for various values of $`m`$ whereas Fig. 5(b) shows the scaling collapse. It is obvious from comparing Fig. 3 and 5 that the layer by layer growth regime is substantially stronger in $`d^{}=2`$ case compared with the $`d^{}=1`$ case, which is consistent with a much larger value of the damping exponent $`\mu 2.5`$ (compared with $`1.5`$ in $`d^{}=1`$) in Fig. 5(b) in the $`d^{}=2`$ system. Our calculated damping exponent $`\mu 2.5`$ for the $`d^{}=2`$ dimensional DT model is substabtially (by $`25\%`$) higher than the corresponding theoretical prediction  of $`\mu =2`$ for the LDV equation (Eq.(5)), which is generally thought to be the continuum description for MBE growth. This large discrepancy between our simulated damping exponent ($`\mu 2.5`$) and the theoretical damping exponent ($`\mu =2`$) corresponding to Eq. (5) in $`d^{}=2`$ dimensions may be a real effect, arising from the recently discovered fact  that the DT model in 2+1 dimensions has actually a very small (but non-zero) EW $`^2h`$ term in its continuum description in contrast to the 1+1 dimensional DT model where the EW $`^2h`$ term strictly vanishes (i.e. $`\nu _2=0`$ in Eq. (7)) by virtue of a topological symmetry  in the DT model. Thus the DT model in 2+1 dimensions may actually asymptotically belong (but with a very small value of $`\nu _2`$) to the EW universality class which according to Eq. (6) has an infinite value of the damping exponent $`\mu `$. It is therefore possible that the value $`\mu 2.5`$ in Fig. 5(b) may be indicative of a small correction in the 2+1 dimensional DT model arising from the EW term in the continuum equation. More work is needed to conclusively settle this issue. This would also explain why the scaling collapse for the damped oscillations in Fig. 5(b), particularly for the data in the kinetically rough ($`t>t_c`$) growth regime, is not as good as the corresponding 1+1 dimensional results shown in Fig. 3 — the asymptotic corrections to the LDV equation (Eq. (5)) arising from a small $`^2h`$ term (which is present in $`d^{}=2`$ DT simulation results of Fig. 5.

The noise reduced $`d^{}=2`$ WV model simulations, shown in Fig. 6, are even more surprising. The layer by layer growth oscillations in the initial transient time (upto $`t10`$ or so) are apparent in Fig. 6(a) for finite values of $`m`$ although the oscillations are already weaker than the corresponding DT results shown in Fig. 5. This is not expected based on the $`d^{}=1`$ results where Fig. 3 (DT) and Fig. 4 (WV) are essentially identical within our simulation sizes and times. Even more surprising is the scaling collapse in the $`d^{}=2`$ WV model shown in Fig. 6(b) which is clearly not a good scaling behavior — in fact, beyond the layer by layer regime (i.e. for $`t>t_c`$) there is essentially no scaling behavior in the $`d^{}=2`$ WV results. The scaling behavior of the WV model (Fig. 6) in the 2+1 dimensional WV model is clearly very different from (and worse than) the corresponding DT results. The best damping exponent value for $`\mu `$ we obtain from Fig. 6(b) is $`\mu 1.5`$, which is the same as the corresponding WV value in $`d^{}=1`$ as depicted in Fig. 4. We should emphasize that this estimate ($`\mu 1.5`$) for the damping exponent in the 2+1 dimensional WV model should be taken at best as a crude estimate for an effective exponent since there is no scaling behavior in the WV data shown in Fig. 6. Given the very similar $`d^{}=1`$ behavior in the DT and WV models (as shown in Figs. 3 and 4) it is very surprising that the $`d^{}=2`$ behavior in the two models (including the effective damping exponent values $`\mu 2.5`$ and $`1.5`$ respectively for the DT and the WV model) is so completely different.

What is the explanation for this striking difference in the $`d^{}=2`$ layer by layer epitaxial growth behavior in the DT (Fig. 5) and the WV (Fig. 6) model (particularly in view of their essentially identical behavior in $`d^{}=1`$ as seen in Figs. 3 and 4)? The explanation actually lies in the recently discovered fact  that, while the $`d^{}=1`$ dimensional WV model obeys  the continuum growth equations given in Eqs. (5) and (7) with $`\nu _4`$ and $`\lambda _{22}0`$ and $`\nu _2`$ very small but having a nonvanishing positive value, the $`d^{}=2`$ dimensional WV growth model is actually unstable in the sense that the WV morphology in the 2+1 dimensional growth forms a regular mounded structure with the mound edges having approximately constant slopes. Such an epitaxial mounding instability  in the 2+1 dimensional WV model becomes particularly manifest under the noise reduction technique as discussed in details in Ref. . This unstable mounded morphological growth  in the noise reduced $`d^{}=2`$ WV model leads to the peculiar behavior seen for late times in Fig. 6, where the epitaxial mounding instability prevents the usual layer by layer growth regime from behaving in the “usual” manner depicted in Figs. 2-5.

Finally, in Fig. 7 we show our 2+1 dimensional Family (F) model simulations under noise reduced ($`m>1`$, $`l=1`$) conditions. The F model, which by construction is designed to follow exactly the linear second order Edwards-Wilkinson equation (Eq. (7)) with $`\nu _4=\lambda _{22}=0`$ in Eq. (5) and $`\nu _20`$, has very simple growth rules (not shown in Fig. 1) : Each randomly deposited atom seeks to find the site of local height minina as the incorporation site. The F model is essentially the discretized version of the EW equation (Eq. (7)), and as such has the EW dynamical exponents : $`\beta =0`$ and $`\mu =\mathrm{}`$ in $`d^{}=2`$. The growth exponent $`\beta `$ (defines the presaturation kinetic roughening of $`W`$ as $`Wt^\beta `$ ) being zero (i.e. $`W\mathrm{ln}t`$) in the F model, growth is already very smooth because the roughening is only logarithmic in time. In the presence of noise reduction, therefore, the F model (Fig. 7(a)) shows persistent layer by layer growth oscillations in $`d^{}=2`$ with only logarithmic damping of the oscillation induced by kinetic surface roughening. Since the damping exponent $`\mu `$ in the noise reduced $`d^{}=2`$ F model is infinity (Eq. (6)), we cannot obtain any scaling collapse of the layer by layer growth simulation data of Fig. 7(a) which is obvious in the “scaling plot” shown in Fig. 7(b). A very large value of the exponent $`\mu `$ ($`100`$, for example) will of course produce trivial (and meaningless) data collapse, but we have checked that no finite reasonably small (upto $`\mu =10`$) value of the damping exponent $`\mu `$ produces scaling in Fg. 7(b). Thus, our $`d^{}=2`$ F results are consistent with the theoretical prediction  of $`\mu `$ being infinity in the $`d^{}=2`$ EW equation. We note that we have also carried out $`d^{}=1`$ dimensional noise reduced layer by layer growth simulations (not presented in this paper) in the F model, obtaining excellent scaling collapse with the theoretically predicted value of $`\mu =2`$. Our results for the $`d^{}=1`$ F model layer by layer growth are consistent with those reported in Ref. .

## IV Conclusion

We have presented numerical results for extensive computer simulations of various noise reduced limited mobility DT , WV , and F  growth models in 1+1 and 2+1 dimensions in order to study the damping of the layer by layer growth epitxy invariably induced by the shot noise inherent in the deposition beam fluctuations. We have used both multiple hit noise reduction technique and the long surface diffusion length method to obtain the layer by layer growth in our simulations, and have shown that these two different techniques for obtaining layer by layer growth are essentially equivalent. Our simulation results in general exhibit (with two exceptions noted below) very good scaling connecting the layer by layer growth regime with the kinetically rough growth regime. Our calculated damping exponents agree well with theoretical predictions where applicable. The two exceptions noted above are the 2+1 dimensional WV and the F model where scaling fails for different reasons. The 2+1 dimensional noise reduced WV model is known  to manifest unstable growth with spectacular epitaxial mounding, which inhibits layer by layer growth leading to the failure of scaling collapse. The 2+1 dimensional F model on the other hand exhibits very strong and persistent layer by layer growth oscillations (with little kinetic roughening) whose damping is expected on theoretical grounds to be extremely weak leading to an infinite value of the damping exponent, which is equivalent to saying that there is essentially no scaling since in the presence of noise reduction layer by layer growth regime lasts forever.

This work has been supported by the US-ONR.


format: text

text:
# STOKES’ DRIFT OF LINEAR DEFECTS

## I Introduction

Particles suspended in a viscous medium traversed by a longitudinal wave $`f(kx\mathrm{\Omega }t)`$ with velocity $`v=\mathrm{\Omega }/k`$ are dragged along in the $`x`$-direction according to a deterministic mechanism known as Stokes’ drift . As a matter of fact, the particles spend slightly more time in regions where the force acts parallel to the direction of propagation than in regions where it acts in the opposite direction. Suppose  that $`f(kx\mathrm{\Omega }t)`$ represents a symmetric square wave with wavelength $`\lambda =2\pi /k`$ and period $`T_\mathrm{\Omega }=2\pi /\mathrm{\Omega }`$, capable of entraining the particles with velocity $`\pm bv`$ (with $`0b(v)1`$ and the signs $`\pm `$ denoting the orientation of the force). Within one cycle, the particle velocity is positive for a longer time interval, $`T_\mathrm{\Omega }/2(1b)v`$, than it is negative, $`T_\mathrm{\Omega }/2(1+b)v`$; hence, the average drift velocity $`ub^2v`$.

Recent studies  have shown that thermal diffusion of the dragged particles can affect markedly the magnitude of the Stokes’ drift. A simple change of variables, $`xykx\mathrm{\Omega }t`$, makes it possible to map the Stokes’ mechanism into the more conventional problem of Brownian motion in a tilted washboard potential ; the dependence of the drift velocity $`u`$ on the viscous constant $`\alpha `$ and the temperature $`T`$ of the suspension medium follows immediately . No matter what the value of $`\alpha `$, the Stokes’ drift is characterized by an onset threshold that grows sharper and sharper as $`T`$ is lowered toward zero; at finite temperature a residual drift occurs in the forbidden parameter region, due to thermal activation (creep). Lately, the Stokes’s mechanism has been attracting the interest of an increasing number of experimentalists in view of its potential applications to single-electron transport in one-dimensional channels , light storage in quantum wells , optical tweezing of colloidal particles , and transport by capillary waves , to mention but a few.

In this article we consider the possibility that a directed linear defect confined to a planar substrate can undergo Stokes’ drift as an effect of a longitudinal wave travelling perpendicularly to its axis. By means of a convenient change of variable, like in Ref. , one expects to recast the question at hand into an instance of the more general problem of the diffusion of an elastic string on a tilted periodic substrate -. Such a mechanism has been used extensively in the literature to model the diffusion of dislocations in solids , the characteristic curve of Josephson junction arrays , the lubrication processes taking place between slipping macroscopic surfaces  and the transport of magnetic flux lines in type-II superconducting films . Here, as a major difference between point-like \[or (0+1)d\] and linear \[or (1+1)d\] defects, we show that: (i) In the limit of infinitely long strings, creep effects and critically sharp Stokes’ thresholds coexhist at any finite temperature; (ii) Such thresholds shift toward lower values with increasing the temperature, due to the configurational string fluctuations (or phonon dressing); (iii) The dispersion of the drifting string around its center of mass grows anomalously with time.

The paper is organized as follows. In Sec. II, we determine the dependence of the drift thresholds on the damping constant and the temperature of the substrate. The onset of the Stokes’ mechanism is sharp for infinitely long strings, whereas finite-length corrections may be computed numerically with no much effort. In Sec. III, the phenomenology of the Stokes’ drift for linear defects is analyzed in detail with particular attention to the excess diffusion of the string center of mass in the vicinity of the threshold and the anomalous dispersion of the string around its center of mass. In Sec. IV, we compare the effectiveness of the Stokes’ drift for (0+1) and (1+1) dimensional defects in solids.

## II Depinning thresholds

As anticipated in the Introduction, the dynamics of a directed elastic string dragged by a longitudinal wave propagating perpendicularly to the string in its glide plane, is intimately related to the diffusion of the very same string on a suitable periodic substrate subjected to an external static bias. Therefore, we address first the problem of the depinning of a damped string from a tilted washboard potential - with the purpose of better characterizing the threshold dynamics.

The perturbed sine-Gordon (SG) equation

$$\varphi _{tt}c_0^2\varphi _{xx}+\omega _0^2\mathrm{sin}\varphi =\alpha \varphi _t+F+\zeta (x,t),$$
(1)
represents the archetypal model of an elastic string in thermal equilibrium on a periodic substrate . Here, $`c_0^2`$ and $`\omega _0^2`$ characterize the string restoring force and the depth of the substrate valleys, respectively; correspondingly, $`c_0`$ is the maximum propagation (or sound) speed of a signal along the string and $`\omega _0`$ is the librational (or plasma) frequency of the string itself at the bottom of a SG valley. The coupling of the classical field $`\varphi (x,t)`$ to the heat bath at temperature $`T`$ is reproduced here by a viscous term $`\alpha \varphi _t`$ and a zero-mean Gaussian noise source $`\zeta (x,t)`$. The damping constant $`\alpha `$ and the noise intensity are related through the noise autocorrelation function

$$\zeta (x,t)\zeta (x^{},t^{})=2\alpha kT\delta (tt^{})\delta (xx^{}).$$
(2)
As the relevant SG substrate potential $`V[\varphi ]=\omega _0^2(1\mathrm{cos}\varphi )`$ is tilted by the bias term $`F\varphi `$, the resulting washboard potential retains its multistable structure only if $`|F|<F_3=\omega _0^2`$ (static depinning threshold ).

The elastic string can glide perpendicularly to the $`V[\varphi ]`$ valleys through three distinct mechanisms:

(a) Kink drift. – A string bridging two adjacent SG valleys bears at least one geometrical kink $`\varphi _+(x,t)`$, or antikink $`\varphi _{}(x,t)`$, with $`\varphi _\pm (+\mathrm{},t)\varphi _\pm (\mathrm{},t)=\pm 2\pi `$,

$$\varphi _\pm (x,t)=4\mathrm{arctan}\left[\mathrm{exp}\left(\pm \frac{xx_0ut}{d\sqrt{1u^2/c_0^2}}\right)\right].$$
(3)
In the presence of a positive bias $`F`$, the solitons $`\varphi _\pm `$ are driven apart to $`\mathrm{}`$, so that the string eventually advances by one step $`2\pi `$. SG solitons are known to behave like Brownian particles with finite size $`d=c_0/\omega _0`$, mass $`M_0=8/d`$ and average drift velocity $`\pm u_F`$, with  $`u_F=2\pi F/\alpha M_0`$. This mechanism of string transport works for $`F0`$, no matter how low the temperature.

(b) Thermal nucleation. – A string segment lying along one SG valley is unstable against nucleation of thermal kink-antikink pairs . Thermal fluctuations $`\zeta (x,t)`$ trigger the process by activating a critical nucleus. For $`FF_3`$ an unstable nucleus is well reproduced by the linear superposition of a kink and an antikink, each experiencing two contrasting forces: a restoring force due to the vicinity of its nucleating partner and a repulsive force due to either the tilt $`F`$ itself, or the preexisting (equilibrium) gas of kinks and antikinks . For $`FF_3`$ the size of the critical nucleus, obtained by balancing the two forces, is quite large, so that its energy is about twice the soliton rest energy $`E_0=M_0c_0^2`$. At low temperature $`kTE_0`$, the nucleation time is thus exponentially long, $`T_N\mathrm{exp}(2E_0/kT)`$. As a result, on establishing opposite, stationary kink-antikink flows, see also mechanisms (a), the string creeps in the $`F`$ direction with average center of mass velocity  $`u_{th}=2\pi n_0u_F`$, proportional to the density of thermal kinks (or antikinks)

$$n_0(T)=\frac{1}{d}\left(\frac{2}{\pi }\right)^{\frac{1}{2}}\left(\frac{E_0}{kT}\right)^{\frac{1}{2}}\mathrm{exp}\left(\frac{E_0}{kT}\right).$$
(4)
(c) Depinning from the substrate. – Even at low temperature, where thermal creep is negligible, the elastic string may unlock from the SG substrate. In the overdamped regime, $`\alpha \omega _0`$, this happens of course for values of $`F`$ above the static threshold $`F_3`$; in the underdamped regime $`\alpha \omega _0`$, the depinning threshold is determined dynamically : at $`T=0+`$, stationary string translation is energetically sustainable above the new threshold

$$F_2=3.36\mathrm{}\left(\frac{\alpha }{\omega _0}\right)F_3.$$
(5)
For $`F>F_2`$, the potential energy gained by the string during one jump forward is larger than the energy dissipated, in average, by the viscous force that opposes the process. \[In the noiseless case $`T0`$ , of no relevance to the present study, the string dynamics would be characterized by hysteretic loops with upper depinning threshold at $`F_3`$ and lower re-pinning threshold $`F_1=(4/\pi )(\alpha /\omega _0)F_3`$.\] On increasing $`\alpha `$ the depinning threshold grows monotonically according to the linear law (5), for $`\alpha \omega _0`$, and levels off at the constant value $`F_3`$ for $`\alpha {}_{}{}^{>}1.2\omega _0`$.

For tilt values above the depinning threshold $`F_2`$, an underdamped string slides downhill with average velocity $`u_\varphi `$; both the velocity of the string center of mass and the lateral motion of its kinks and antikinks are modulated periodically in time with period $`2\pi /u_\varphi `$. As the upper bound of $`u_\varphi `$ is $`F/\alpha `$ , the amplitude of the soliton sidewise oscillations is certainly smaller than $`(1/2)[2\pi /(F_2/\alpha )]c_0`$, where $`F_2`$ is given by Eq. (5) and $`c_0`$ is the limiting speed of any pulse propagating along the string, and therefore smaller than the soliton size $`d`$. Moreover, for $`kTE_0`$ the mean kink-antikink distance $`n_0^1(T)`$ is much larger than $`d`$, viz. $`n_0(T)d1`$, see Eq. (4); as a consequence, the spatial configuration of the falling string is only negligibly modulated by its passage over the substrate barriers. The same conclusion applies in the opposite regime $`\alpha \omega _0`$.

Such a periodic effect can be averaged out by adding explicitly a translational coordinate $`X`$ to the field $`\varphi `$, namely

$$\varphi (x,t)X(t)+\varphi (x,t).$$
(6)
On taking the spatial average

$$\mathrm{}_x=\underset{L\mathrm{}}{lim}\frac{1}{L}_{L/2}^{L/2}(\mathrm{})𝑑x$$
(7)
\[not to be mistaken for the stochastic average $`\mathrm{}`$ in Eq. (2)\] of both sides of Eq. (1) along a string segment of length $`L`$, one obtains

$$\ddot{X}=\alpha (\dot{X}+u_{th})U^{}(X,T)+F+\frac{1}{\sqrt{L}}\eta (t),$$
(8)
where we set

$$\varphi _t(x,t)_x=\dot{X}+u_{th},$$
(9)
$$\varphi _{tt}(x,t)_x=\ddot{X},$$
(10)
$$\zeta (x,t)_x=\frac{1}{\sqrt{L}}\eta (t),$$
(11)
$$V^{}[\varphi ]_x=U^{}(X,T).$$
(12)
In Eqs. (9) and (10) we have made use of the reflection symmetry $`_t\varphi _\pm =_t\varphi _{}`$ of the $`\varphi _\pm (x,t)`$ solutions. Note that for $`kTE_0`$ the $`\varphi _\pm `$ center of mass $`x_0`$ fluctuates  with velocity variance $`u^2=kT/E_0c_0^2`$. In Eq. (9) $`u_{th}`$ denotes the creep velocity of the string due to the lateral drift of kinks and antikinks, see mechanism (b) above.

The spatial average of the random force in Eqs. (1)-(2) is statistically equivalent to the Gaussian noise (11), where the overall intensity is inverse proportional to the string length $`L`$, $`\eta (t)=0`$ and

$$\eta (t)\eta (t^{})=2\alpha kT\delta (tt^{}).$$
(13)
The effective potential $`U(X,T)`$ can be computed only under certain low-temperature approximations. Of course, it retains all the spatial symmetries of $`V[\varphi ]`$, namely $`U(X,T)=U(X,T)`$ and $`U(X+2\pi ,T)=U(X,T)`$. In view of the rigid translation (6), the identities $`V^{}[\varphi ]_x\omega _0^2\mathrm{sin}(X+\varphi )_x=\omega _0^2[\mathrm{sin}X\mathrm{cos}\varphi _x+\mathrm{cos}X\mathrm{sin}\varphi _x]`$ and $`\mathrm{sin}\varphi _x=0`$, lead to

$$U(X,T)=V(X)\mathrm{cos}\varphi _x$$
(14)
with $`V(X)=\omega _0^2(1\mathrm{cos}X)`$. \[Note that the averages (7) $`\mathrm{cos}\varphi _x`$ and $`\mathrm{sin}\varphi _x`$ are to be taken along the kink profile.\] The $`T`$dependent configurational factor $`\mathrm{cos}\varphi _x`$ deviates from unity because of (i) the local $`\varphi `$ steps (solitons), (ii) the extended $`\varphi `$ excitations (phonons) sustained by thermal fluctuations. The first correction is proportional to the (anti)kink density $`n_0(T)`$ and therefore exponentially small; neglecting the soliton contribution corresponds to replacing $`\mathrm{cos}\varphi _x`$ with $`\mathrm{cos}\psi _x`$, where $`\psi (x,t)`$ represents the equilibrium phonon gas of a string sitting at the bottom of one SG valley. The field $`\psi (x,t)`$ can be regarded  as a linear superposition of phonon modes $`\psi _k(x,t)`$, with dispersion $`\omega _k^2=\omega _0^2+c_0^2k^2`$ and average amplitude set by the equipartition theorem $`\psi _k^2_x=kT/\omega _k^2`$; hence

$$\psi ^2_x=\frac{kT}{2\pi }_{\mathrm{}}^+\mathrm{}\frac{dk}{\omega _0^2+c_0^2k^2}=4\frac{kT}{E_0}.$$
(15)
At low temperature the distribution of the phonon field is well approximated by a Gaussian function $`P[\varphi ]=\mathrm{exp}(\psi ^2/2\psi ^2_x)/\sqrt{2\pi \psi ^2_x}`$, so that

$$\mathrm{cos}\varphi _x=\mathrm{}\left\{e^{i\psi }P[\psi ]𝑑\psi \right\}=e^{\psi ^2_x/2}.$$
(16)
On combining Eqs. (14) and (16) we conclude that

$$U(X,T)=\omega ^2(T)(1\mathrm{cos}X)$$
(17)
with

$$\omega (T)=\omega _0e^{kT/E_0};$$
(18)
hence, the $`T`$ dependence of the string depinning threshold (5).

Let us discuss now the diffusion of the string center of mass described by the reduced Eq. (8) or, after the Galileian transformation $`Y=X+u_{th}t`$, by

$$\ddot{Y}=\alpha \dot{Y}\omega ^2(T)\mathrm{sin}(Yu_{th}t)+F+\frac{1}{\sqrt{L}}\eta (t).$$
(19)
The processes $`X(t)`$ and $`Y(t)`$ represent a driven Brownian motion on a tilted washboard potential . In the limit of an infinitely long string the effective temperature of the process (8),(19) $`T^{}=T/\sqrt{L}`$ tends to zero. This implies, for instance, that the depinning threshold of an underdamped string is well approximated by Eq. (5) after $`\omega _0`$ was replaced by $`\omega (T)`$. As shown in Risken’s textbook , the depinning mechanism is abrupt; for $`F=F_2(T)+`$ the average drift velocity of the string $`u_\varphi \dot{X}`$ jumps all of a sudden from its creep value $`u_{th}`$ up to close to its asymptotic value $`F_2(T)/\alpha `$. In Eqs. (8) and (19) the dependence on the string temperature $`T`$ is restricted to the effective potential amplitude (18) \[as a result the threshold value $`F_2(T)`$ gets lowered\] and to the creep velocity $`u_{th}(T)`$. Note that even the dispersion of the creep velocity tends to vanish in the limit $`L\mathrm{}`$, as the creep mechanism through pair nucleation is essentially a shot noise process .

Finite-length corrections to the picture outlined above can be computed explicitly by applying the numerical recipes of Ref.  for a finite effective temperature $`T^{}`$. Depinning curves for vanishingly small $`T^{}`$ values are displayed in Ref. ; anticipation of the depinning thresholds with the string temperature is hinted at, though not explicitly discussed, in Refs. .

## III Stokes’ drift

We go back now to the initial problem of a directed elastic string (say, a lattice dislocation) confined to a planar substrate and subjected to a longitudinal wave $`f(\varphi vt)`$ (a stress wave) propagating perpendicularly to the string direction. Let the string be coupled to a viscous heat bath at temperature $`T`$; the corresponding field equation reads

$$\varphi _{tt}c_0^2\varphi _{xx}=\alpha \varphi _t+f(\varphi vt)+\zeta (x,t),$$
(20)
where $`\zeta (x,t)`$ is defined in Eq. (2) and $`v=\mathrm{\Omega }/k`$ denotes the propagation speed of a travelling wave with wavelength $`2\pi /k`$. By means of the Galileian transformation $`\varphi \psi \varphi vt`$, Eq. (20) can be rewritten as

$$\psi _{tt}c_0^2\psi _{xx}=\alpha \psi _t+f(\psi )\alpha v+\zeta (x,t).$$
(21)
For a sinusoidal wave with $`k=1`$, $`f(\psi )=f_0^2\mathrm{sin}\psi `$, this equation coincides with Eq. (1) under the correspondence rule

$$\omega _0f_0,F\alpha v,$$
(22)
and $`\alpha ,T`$ and $`c_0`$ unchanged.

### A Onset thresholds

The center of mass of the string $`\psi (x,t)`$ diffuses according to Eqs. (8) and (19), after transformation (22). As long as $`\psi (x,t)`$ is locked to the bottom of a (creeping) SG valley, see Eq. (19), the initial string $`\varphi (x,t)`$ drifts with velocity $`u_\varphi v`$; upon depinning, the center of mass of $`\psi (x,t)`$ travels to the left with negative velocity $`u_\varphi F/\alpha =v`$; correspondingly, $`\varphi (x,t)`$ gets traversed by the propagating wave without undergoing any appreciable drag, were not for the residual creep velocity $`u_{th}`$.

The relevant Stokes’ threshold $`v_S`$ admits two simple limiting approximations; in the overdamped regime $`\alpha \omega _0`$

$$v_S(\alpha \mathrm{})=v_3=\frac{F_3}{\alpha }=\frac{f^2(T)}{\alpha },$$
(23)
in the underdamped regime $`\alpha \omega _0`$

$$v_S(\alpha 0)=v_2=\frac{F_2}{\alpha }=3.36\mathrm{}f(T),$$
(24)
where $`f(T)=f_0\mathrm{exp}(kT/E_0)`$. Stokes’ drift requires that $`v<v_S`$, i.e. travelling waves must be selected with sufficiently low frequency or short wavelength. At low damping the threshold velocity $`v_Sv_2`$ is insensitive to $`\alpha `$ and the drift effect is more pronounced. In the limit of infinitely long strings, $`L\mathrm{}`$, the threshold mechanism grows sharper and sharper, as to be expected for the diffusion of a Brownian particle in a tilted washboard potential at zero temperature $`T^{}`$. At variance with the (0+1) dimensional case , the onset of the Stokes’ drift for an elastic string in thermal equilibrium is thus a critical transition, even at finite (but relatively low) temperature.

### B Excess diffusion

As proposed in Ref. , a dynamic threshold can be characterized in terms of the diffusion coefficient

$$D=\underset{t\mathrm{}}{lim}\frac{X^2(t)X(t)^2}{2t}$$
(25)
of the process at hand. Here, the function $`D(v)`$ tends to an exponentially small value $`D(0)\mathrm{exp}(kT^{}/f(T)^2)`$ for $`v0`$, and to the Einstein constant $`D(\mathrm{})=kT^{}/\alpha `$ for $`v\mathrm{}`$ . Both limits vanish for an infinitely long string, $`T^{}=T/\sqrt{L}0`$. More notably, in the vicinity of the Stokes’ threshold $`v_S`$, $`D(vv_S)`$ is expected to shoot up drammatically; such an effect, know as excess diffusion, may be enhanced either by lowering the temperature or by weakening the viscous coupling, independently. In the overdamped regime (better said, in the Smoluchowski limit $`\alpha \mathrm{}`$) and for $`T^{}0`$, the diffusion constant $`D(v_3)`$ decays much slower than $`D(\mathrm{})`$, that is

$$\frac{D(v_3)}{D(\mathrm{})}G_0\left[\frac{f^2(T)}{kT^{}}\right]^{\frac{2}{3}}L^{\frac{1}{3}},$$
(26)
with $`G_00.90`$. In the underdamped regime, the ratio $`D(v_2)/D(\mathrm{})`$ turns out to grow even faster, although an appropriate scaling law has not been determined, yet . In both instances the excess diffusion peak is predicted to be centered very narrowly around the relevant threshold velocity $`v_S`$.

### C Anomalous dispersion

Away from the Stokes’ threshold $`v_S`$, the string either drifts along with the travelling wave, if $`v<v_S`$, or creeps forward in the same direction, but with marginal velocity $`u_{th}v\mathrm{exp}(E_0/kT)`$, if $`v>v_S`$. The spatial configuration of $`\varphi (x,t)`$, or equivalently of $`\psi (x,t)`$, is markedly different in the two opposite dynamical regimes. The creeping string $`\varphi (x,t)`$ corresponds to the case of the string $`\psi (x,t)`$ falling down a tilted periodic substrate (21) ($`v>v_S`$): As shown in Ref. , under stationary periodic boundary conditions the dispersion

$$S(t)=\psi ^2(x,t)_x\psi (x,t)_x^2,$$
(27)
of an unlocked string segment of finite length $`L`$ is insensitive to the substrate, $`S(\mathrm{})=kTL/4c_0^2`$. Accordingly, the spatial dispersion of a creeping string with $`L\mathrm{}`$ is expected to be the same as for a directed elastic string diffusing in a plane, free of periodic modulations, either in space or in time, i.e.

$$S(t)=\frac{kT}{c_0}\left(\frac{t}{\pi \alpha }\right)^{\frac{1}{2}}.$$
(28)
The counterpart of a drifting string $`\varphi (x,t)`$ is a string $`\psi (x,t)`$ locked to the tilted substrate (21); therefore, the spatial configurations of both strings result in a linear superposition of diffusing kinks and antikinks of size $`d=c_0/f_0`$, separated by an average distance $`n_0^1(T)`$. Here the soliton width depends on the wave amplitude $`f_0`$ and on the string elastic constant $`c_0`$, while the soliton height concides with the wavelength $`2\pi `$ (for $`k=1`$). The question as how the dispersion $`S(t)`$ of an infinitely long string pinned to its SG substrate evolves with time, due to thermal pair nucleation, has been addressed in Ref. : Under the further assumption that kink-antikink collisions are destructive, namely for

$$v2\left(\frac{2\alpha }{f_0}\right)^2v_3,$$
(29)
the creep dispersion grows anomalous with fractional exponent larger than 1/2, that is

$$S(t)0.32(u_{th}t)^{2/3}.$$
(30)
Condition (29), clearly satisfied in the overdamped regime, may apply to lower damping values, too, provided that $`v`$ is taken small enough.

In conclusion, while for $`v>v_S`$ the string $`\varphi (x,t)`$ is almost transparent to the travelling wave $`f(\varphi vt)`$, for $`v<v_S`$ the string spreads out anomalously with exponent 2/3 and its center of mass drifts with velocity $`u_\varphi `$ close to $`v`$.

## IV Conclusions

Longitudinal waves in solids may drag along point-like and linear defects, alike, a mechanism known as Stokes’ drift. The onset of such a transport process, though, exhibits distinct properties depending on the dimensionality of the movable defects involved:
(i) The depinning threshold is very sharp in the case of a long defect line, no matter what the substrate temperature. Increasing the temperature may push indeed the Stokes’ threshold $`v_S`$ toward lower values, but causes no broadening of the activation region, at variance with the depinning of point-like defects. In view of this remark, we conclude that Stokes’ drift of linear defects is characterized by a clear-cut threshold signature.
(ii) Due to their spatial extension, linear defects subjected to Stokes’ drift undergo a substantial deformation in the form of an enhanced transversal dispersion. Moreover, a dragged defect line acquires a multi-kink profile, where the geometry of each individual (anti)kink component depends on the incoming wave parameters. The ensuing internal dynamics of the transported linear defects adds to the dissipative properties of the propagation medium overall .