PHYSICS REPORTS (Review Section of Physics Letters) 216, No. I (1992)! —61. North-Holland
PH Y S I CS R E PORTS
Propagation and scattering of nonlinear waves in disordered systems Sergey A. Gredeskula and Yuri S. Kivsharh* ~Departmentof Physics, Ben-Gurion University of the Negev. P.O. Box 653, 84105 Beer-Sheva, Israel hlnstitute for Low Temperature Physics and Engineering, 310164 Kharkov, Ukraine Received December 1991; editor: A.A. Maradudin
Contents: 1. Introduction 2. Linear problems 2.1. Lyapunov exponent and localization length 2.2. Plane-wave transmission through a disordered segment 2.3. Transmission of wavepackets 2.4. Generation of dark solitons by a random input pulse 3. Stationary nonlinear transmission 3.1. Preliminary remarks 3.2. Fixed output problem 3.3. Fixed input problem. Multistability and disorder
3 4 5 7 10 12 15 16 17 29
4. Modulational instability and wave propagation 4.1. Modulational instability in homogeneous models 4.2. Influence of inhomogeneities on modulational instability 5. Soliton propagation through disordered media 5.1. Dynamical solitons 5.2. Envelope solitons 5.3. Topological solitons (kinks) 6. Concluding remarks References
33 33 34 37 38 41 46 5$ 5()
Abstract: We review recent numerical and analytical results related to the study of the coexistence of nonlinear wave propagation and disorder effects. mostly for soliton bearing systems. First of all, we briefly summarize scattering problems in the linear approximation which show many features related to the Anderson localization. Considering wavepackets instead of standard problems of the plane wave propagation, we show that the change of the asymptotic behaviour of the averaged transmission coefficient as a function of the disordered segment length may be demonstrated even at the level of the linear approach. For nonlinear problems we begin by discussing the stationary case which may be analyzed in detail, The main effect of nonlinearity, bistability, may be drastically modified by disorder in this case. However, for nonstationary problems disorder may give rise to more intensive development of nonlinear properties of physical systems via modulational instability, which is enhanced by inhomogeneities. We point out that localization effects created by disorder may vanish in the presence of strong nonlinearity, and we present examples showing that nonlinearity leads to an actual improvement of the transmission when it contributes to create soliton pulses. We analyze soliton propagation through disordered and inhomogeneous media and demonstrate that nonlinear transmission strongly depends on the soliton type: results are different for dynamical, topological, and envelope solitons.
*
Present address: Institut fur Theoretische Physik I, Heinrich-Heine-Universität Düsseldorf 1, D-4000 Düsseldorf, Germany.
0370-l573/92/$15.OO
©
1992 Elsevier Science Publishers By. All rights reserved
PROPAGATION AND SCATTERING OF NONLINEAR WAVES IN DISORDERED SYSTEMS
Sergey A. GREDESKUL Department of Physics, Ben-Gurion University of the Negev, P.O. Box 653, 84105 Beer Sheva, Israel
and Yuri S. KIVSHAR Institute for Low Temperature Physics and Engineering, 310164 Kharkov, Ukraine
NORTH-HOLLAND
1. Introduction Wave propagation in nonlinear disordered media has become an extensively studied subject in recent years [Bass et al. 1988; Bishop et al. 1989; Abdullaev et a!. 1991; Sanchez and Vázquez 1991]. In linear systems, disorder is generally connected with Anderson localization, a concept that applies to a number of different waves like phonons, acoustic and electromagnetic waves, etc. (see, e.g., Lifshitz et al. [1988],Ping Sheng [1990]and references therein). Localization means that the transmission coefficient T of a plane wave decays exponentially with the system length L and, in particular, that a positive finite number exists, the so-called localization length A(k) where k is the wave number, such that for L ~ A(k) very little transmission is allowed. Narrow linear wave packets behave in an analogous fashion [Kivshar et al. 1990, 1992]. The study of the influence of nonlinearity on the transport properties of disordered systems was started by Payton et al. [1967],who investigated numerically the role of nonlinearity on heat flow in a discrete system. They emphasized that heat flow is enhanced by introducing nonlinearity into the system (see also the recent interesting papers by Bourbonnais and Maynard [1990a,b] on the same subject). However, at the time they studied the problem, the concept of soliton had not yet been established and, as a matter of fact, the relation between disorder and soliton propagation had not been observed. In the framework of stationary nonlinear equations, the action of nonlinearity seems to be opposite to that of disorder. Nonlinearity may change the dependence of the transmission coefficient on the length, that still tends to zero as the size of the system increases, but following a power law [Devillard and Souillard 1986; Douçot and Rammal 1987a, b; Knapp et al. 1989, 1991] instead of the exponential behaviour exhibited by the linear case. Nonlinearity may also produce multistability in the wave transmission through a disordered system [Knapp et a!. 1989, 1991]. In nonstationary problems, nonlinearity may lead to modulational instability which is enhanced by disorder [Caputo et al. 1989; Kivshar 1991b]. Modulational instability is important because it is the factor that causes appearance of solitons instead of plane waves. The most remarkable manifestation of the influence of nonlinearity on the properties of disordered systems is related to the fact that many nonlinear systems support undistorted propagation of localized waves, usually called solitons. As is well known, there are solitons of three general types: dynamical solitons, envelope solitons, and topological kinks. The amplitude of a dynamical soliton (or a dynamical kink) is proportional to its velocity, and the scattering of such an object in a disordered medium is similar to that of a linear wave. This means that the scattering of such solitons may be treated as the scattering of a superposition of linear waves with different wave numbers that usually yields power-like dependence of the transmission coefficient on the system length; nonlinearity only slightly modifies this dependence (see Li et al. [1988b],Ishiwata et al. [1990]).Envelope solitons are two-parameter ones and they exhibit more complicated behaviour in disordered media. As was demonstrated by Kivshar et al. [1990], in the case of the envelope soliton, localization effects may vanish drastically above a certain threshold in nonlinearity intensity, which implies that the influence of disorder on the propagation of such intense pulses can become negligible. A similar conclusion follows from the numerical study of the energy transport in anharmonic lattices with disorder: a regime of nonlinear diffusion is characterized 3
4
.S’.A. Gredeskul. Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
by an exponent which depends strongly on anharmonicity but it is insensitive to disorder [Bourbonnais and Maynard 1990a, b]. Another type of solitary waves is topological kinks, which may exist only in nonlinear systems with two (or more) equivalent ground states, i.e.. when nonlinearity is strong. In a number of cases kink dynamics in disordered systems may be described, in the framework of a collective-variable approach. by an equation for an effective (classical or relativistic) particle, its coordinate being the kink position [Pascual and Vázquez 1985; Bass et al. 1987, 1988; Pascual et al. 1989; Biller and Petruccione 1990; Petruccione and Biller 1990; Gredeskul et al. 1992]. This situation is similar to that for perturbationinduced dynamics of regular nonlinear systems supporting the kink propagation [Fogel et al. 1976, 1977; Kivshar and Malomed 1989]. However, the kink scattering may demonstrate some resonant effects which are related to resonant excitation of impurity modes [Kivshar et al. 1991b], and this mechanism may induce additional energy loss of the kink in a disordered medium. The present review aims at demonstrating a number of features in common to both nonlinear dynamics and disorder effects, mostly for soliton bearing systems. Chapter 2 briefly summarizes scattering of linear waves and linear wavepackets which we will use later for comparison with the scattering of nonlinear waves. In ch. 3 we present results for stationary nonlinear problems. Chapter 4 gives examples demonstrating that the typical approach generally used in the theory of localization in disordered linear systems is not suitable for some nonlinear models due to the modulational instability phenomenon. We also point out that modulational instability may be enhanced by disorder, and that it gives rise to generation of solitons. This is the reason why we study soliton propagation extensively. Chapter 5 reviews the soliton propagation through disordered media for different types of solitons; this propagation has a number of distinctive features for each type of nonlinear localized waves. Chapter 6 concludes the paper and presents some prospects in the study of the effects produced by a strong nonlinearity in the presence of disorder,
2. Linear problems Transmission problems are the most important and interesting in the theory of one-dimensional disordered systems. The main reason is that transmission characteristics are closely connected with various kinetic coefficients or other interesting physical quantities. Indeed, let us mention a few important physical examples demonstrating such a connection. Let T= t~2 and R = r~2be the transmittivity and reflectivity of a finite one-dimensional disordered system. t and r being the transmission and reflection amplitudes, respectively. These data are enough to compute some transport characteristics. Indeed, the conductance G of such a system may be represented in the form (the so-called Landauer formula [Landauer 1970]) e2 T G = 2~F11— T’ e being the electron charge. As is well known (see, e.g., Brekhovskikh [1973]), in a layered medium the field F(R, R 0) of a point source placed at the point R0 is given by F(R, R~~) = ~
f
11(kpsin O)f d~~—r(z0)r (z1) H~1
1
2(0,
z) sin 0.
S. A. Gredeskul. Yu. S. Kivshar, Propagation and scattering of nonlinear waves in di.sordered systems
5
where, in particular, r,(z0) are the reflection amplitudes from two semispaces separated from each other by the plane z = z0. The mean heat flux J through a disordered segment is expressed by the formula [Keller et al. 1978], J=2cldk (T)~(kc), where ~T) stands for the averaged transmittivity. The intensity 1(z0) of a wave incident from the right on the disordered layer inside the layer equals 2I~1 r(z 2 ‘(2.1) 1(z0) = I(0)~1+ r(z0)~ 0)~ where r(z) is a reflection amplitude from the semispace z < z 0. At last, the rate Q of a sound pulse —
,
dissipating in a pure superconductor has the form (see Bratus’ and Shumeiko [1985]),
Q
=
~N(0)
~f ~f VF
dE (R)(r~
-
6
)[n~(r~) n(e)], -
~
where R is the reflectivity of an electron from the deformation potential created by the sound pulse in the reference frame moving with the pulse. Moreover, transmission characteristics are sensitive to the dynamical properties of physical systems and, therefore, they contain information about the nature of states within a disordered segment. The present chapter reviews different transmission problems in the framework of linear disordered models. In the first section we briefly discuss the Anderson localization phenomenon in onedimensional systems, and define several important characteristics of the scattering as the Lyapunov exponent and the localization length. The question about quasistationary states in open disordered systems is discussed too. The problem which is under consideration in the second section is plane-wave transmission through a qisordered segment. After a qualitative description of the problem, we obtain the imbedding equation for the reflection amplitude and then, using the technique of the Fokker— Planck equation, we find the averaged transmittivity of the system. The main feature of the transmittivity in linear disordered systems is its exponential decrease with the growth of the thickness of the disordered layer. The third section contains results for the wave packet transmission. In this case the asymptotics of the transmittivity as a function of the length of the disordered segment may show not an exponential but a power-law behaviour, like the asymptotic behaviour of certain relations between the parameters of the system and those of the wave packet. In the last (fourth) section, creation of dark solitons by a random input pulse is considered for the case of nonlinear wave propagation in optical fibers. We have decided to include this section into this chapter, because this nonlinear problem may be reduced to the linear transmission problem for the Dirac equation with a random potential. 2.1. Lyapunov exponent and localization length To demonstrate the genesis of such a phenomenon as the Anderson localization in one-dimensional systems and to show the influence of localization effects on the transmission properties, we will consider the simplest (and well-known) example, namely, the stationary Schrodinger equation, 2tfrldx2+U(x)tfr=Ei/i, E~k2. (2.2) —d
6
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering
01
nonlinear waves in disordered systems
Let us look for a real solution of eq. (2.2) and parametrize it by the formulas, ~=e~sin~,
d~/dx=ke~cos~.
(2.3a,b)
Then, the functions ~(x) and ~(x) satisfy the following system of equations: d~Idx= k
—
[U(x)/k] sin2~,
d~/dx= [U(x)/2k] sin(24).
(2.4a, b)
The function ~(x) is the phase of the solution, and it is responsible for its oscillations, but the function ~(x)describes the amplitude. We will assume that the potential U(x) in eq. (2.2) is a spatially homogeneous random function with correlations vanishing at infinity. Then the following statement is valid. For any boundary value çb(O) = ~ [whatever the value ~(0) may bel a nonrandom positive limit exists (with probability 1),
~L~’i ~(x)/x=
lim K~(x))Ix=l/21(k2)(>0)
(2.5)
-
(The angular brackets ( ~t always denote averaging over all realizations of the random potential.) This property means that any solution of the Schrödinger equation (2.2) with a fixed logarithmic derivative is exponentially increasing outside the initial point (all details concerning the linear disordered ~ystemsmay be found, e.g., in the book by Lifshitz et al. [1988]).The quantity 1(k2) is called localization length. More precisely, the r.h.s. of eq. (2.5) is the so-called Lyapunov exponent of eq. (2.2). This quantity describes very well localization properties of solutions at high energies, i.e., for E = k2 ~ U, where U is a characteristic value of the potential U(x). In the intermediate region E U the quantity 1(k2) from eq. (2.5) more or less coincides with the localization length defined by the Anderson function. However, on the fluctuation spectrum the function 1(k2) has nothing in common with the true localization length because the latter quantity has to be connected with the size of the localization region, whereas the former one defines only the tails of the solution. Nevertheless, we will call the quantity 1(k2) defined in eq. (2.5) the localization length. lithe random potential U(x) in eq. (2.2) is a Gaussian white noise, i.e., ...
—~
(U(x))
=
0,
(U(x)U(0)~t= 2D6(x),
then the localization length in the high-energy region, E l~(k2)—~2k2!D,k2>>D2’3.
(2.6) ~‘
D2’3, is proportional to the energy k2, (2.7)
This result is valid also in the more general case when the characteristic parameters k~and r~of the correlation function, which are defined as KU(x)U(0)~= k~f(x/r~) [where in the function f(x) we have introduced both the amplitude and distance on which it changes], satisfy the inequality k~r~ ‘~81. In this case the energy values are also bounded from above, so that the region of validity of eq. (2.7) may be defined as D2~3~ k2 ~ r~2.
(2.8)
S.A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
7
If we deal with the one-dimensional Helmholtz equation, d2uldx2
+
k2[1
—
s(x)]u
=
0,
(2.9)
then it is easy to notice that after the change u ~s, k2r(x) U , we will obtain the same Schrödinger equation (2.2). Moreover, the inequality (2.8) is always satisfied when k is sufficiently small, because in the case of the Helmholtz equation D = k4d, where 2d = f dx (s(x) e(0)). As a result, the localization length for the Helmholtz equation is inversely proportional to k2, —~
—*
—
l~(k2)= 21k2d.
(2.10)
This result and its connection with the previous one (eq. 2.7) may be easily understood taking into account that in the vicinity of the stable (not fluctuating) boundary of the spectrum the localization length will diverge. In the case of the Schrodinger equation the stable boundary is E~= ~, but for the Helmholtz equation it equals E~= 0. To conclude this section we would like to note that the positiveness of the Lyapunov exponent for all values of energy means that the spectrum of eq. (2.2) with spatially homogeneous potential and with correlations disappearing at infinity is purely point-like, i.e., all states are localized. We will comment on this result in more detail below (see subsection 2.4.2). 2.2. Plane-wave transmission through a disordered segment 2.2.1. Qualitative picture Let us consider now the Schrodinger equation (2.2) with a random potential which differs from zero only within the interval 0< x < L. If the plane wave with wavenumber k is incident from the right of the disordered segment, then the wave outside the segment may be presented in the form,
I e~~*~) + r eik~_I~ x> L x<0 ,
~IJ(X)G(1~,te_ikx
(2.11)
r and t being the reflection and transmission amplitudes. So, at the boundaries of the segment the solution of eq. (2.2) satisfies the following conditions:
~~-1
~f
x=o =
—ik,
d.ix=L
~
—ik
=
~—~--~
(2.12)
.
We can express the solution inside the segment in terms of two linearly independent solutions
i/ia
and ~i 5
of the type (2.3), where the functions (~, 4~)and (~,~) satisfy the system of equations (2.4) and the boundary conditions ~‘
~s’
4~s~=~ = 0 ,
=
4~c~x=0
As a result, for the transmittivity T, defined as 2,
(2.13)
T= t~ we may obtain the expression [Pastur and Fel’dman 1974]: T= 4{2
+
exp[2~~(L)]+ exp[2~
1. 5(L)]}
(2.14)
8
S.A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
From eq. (2.14) and the statement (2.5) it follows that for a long segment, the transmitti~vityon a typical potential realization decreases exponentially, with the decrement coinciding with the inverse localization length, T’—exp(—L/l),
(2.15)
L~’1.
It is necessary to note that the qualitative estimate (2.15) is valid only with logarithmic accuracy because the quantity exp(2~)fluctuates exponentially. So, in a rigorous sense eq. (2.15) has to be written in the form
L~~lnT=1
(2.16)
,
with probability 1. Besides, the transmittivity (2.13) is not a self-averaging quantity but its decrement (2.16) is self-averaging for large L. From the estimate (2.15) it follows immediately that on a typical realization, the difference between unity and the modulus of the reflection amplitude is also exponentially small for large L, i.e., 1— r~-’exp(—LI1)
(2.17)
.
Notice that’ such an exponential [in the sense of eq. (2.17)] behaviour of the transmittivity on a realization does not mean that the averaged transmittivity also decreases exponentially, due to the absence of the self-averaging property. 2.2.2. Embedding approach and averaged transmittivity The problem of average transmittivity asymptotics for large L has been considered in many papers, among which we would like to mention those by Gerzenstein and Vasilyev [l959a, b] and Papanicolaou [1971] (see also the books by Klyatskin [1980]and Lifshitz et al. [1988]). Let us introduce a new dependent variable in the Schrodinger equation (2.2), namely, the logarithmic derivative of the wavefunction, z
=
~
d~Idx.
(2.18)
This quantity satisfies the nonlinear Ricatti equation, dzldx= —(z2 + k2)+ U(x),
(2.19)
and the boundary conditions, z(0)
=
—ik ,
z(L) = —ik[1
—
r(L)1/[1
+
r(L)]
-
(2.20)
Then, let us consider the function r(x)
[ik+ z(x)] I[ik
—
z(x)1 ,
which satisfies the Ricatti equation
(2.21)
S.A. Gredeskul, Yu,S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
dridx
=
2ikr
—
(iI2k)U(x)(1 + r)2,
9
(2.22)
as follows directly from eq. (2.19). If we demand that r(x) satisfy the boundary condition r(0) = 0, which is a consequence of the first condition given in eq. (2.20), then the meaning of this function at point L coincides with the reflection amplitude r(L). If we introduce the modulus p and the phase 0 of the reflection amplitude, i.e., r = p eiO, then the functions p and 0 satisfy the system of equations dpidx = —(U12k)(1 dOidx = 2k
—
—
(U12k)[2
p2) sin 0,
(2.23)
(p
(2.24)
+
+
lip) cos 0].
In the standard case the potential U(x) is considered as a Gaussian white noise, and the transmission problem in the high-energy region (2.8) can be reduced to a Fokker—Planck equation for the probability density of the reflection amplitude modulus p (see, e.g., the books mentioned at the beginning of the subsection). The simplest way to obtain this Fokker—Planck equation was firstly proposed by Antzygina et al. [1980](recently, in the paper by Doucot and Rammal [1987b]the same method was proposed once again). We will not discuss technical details of this approach, but we will present only the final result for the average transmittivity, (T~ ~ir512(lIL)312 exp( —
—
L141).
(2.25)
It is seen from this formula that the average transmittivity decreases also exponentially with the growth of the length L of the disordered system. But the connection between this exponential behaviour and the exponential decrease of the transmittivity on a realization is not straightforward, because the main contribution to the averaged quantity (2.25) is due to resonant transparent realizations with T~1, and the exponential dependence in eq. (2.25) reflects only the measure of such realizations. 2.2.3. Open systems and quasistationary states In the case when the potential U(x) in the Schrödinger equation (2.2) differs from zero only in the segment [0, U, the continuum spectrum of the problem coincides with the positive semi-axis 0< E < Each positive value of the energy E > 0 is two-fold degenerate. One of the states corresponding to it has at x = +°~only an outgoing wave, ~.
~(z, E)
[1—r(E)r+(E)] eikx + rt(1 t~(E)t 2(E)~
—
e~x,
X <0,
(2.26)
U
Here r~(E)[r_(E)] is the reflection amplitude of a plane wave incident from the left (right) with energy E from the semisegment [x0, L] ([0, x0]); x0 is an arbitrary point within the interval [0, L]; r1(E) is the reflection amplitude of a plane wave incident from the left with energy E from the semisegment [0, x0]; t12(E) are the transmission amplitudes of the semisegments, and the asterisk stands for the complex-conjugate value. It is seen from eq. (2.26) that on a typical realization when the moduli of the reflection amplitudes r±(E)are exponentially close to unity, the amplitude 1 r_(E)r~(E)of the initial wave becomes exponentially small at the points E~for which the product r_(E)r÷(E)is real and positive. If the —
l()
S. A. Gredeskul, Yu. S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
segment length L is sufficiently large in comparison with the localization length 1(E,), i.e. L then we may obtain for the phase of the reflection coefficient the result E~= (nITIL)2[1 + o(1)].
>
(2.27)
The second term in the r.h.s. of this equation essentially depends on the specific realization and it is small in comparison with the universal first one, but it is not small in comparison with the distance (2~r2n/L2)between the neighbouring points E, and E~.In the 1i~’-complexplane in the vicinities of the points E, of the real axis there exist points ~,, such that Im ~, <0, ~,— El exp(— LII), and they satisfy the equation -~
1
—
r(~, 7)r5( ) = 0.
(2.28)
Then in expression (2.26) for the wavefunction only outgoing waves are present. So, the function t/i(x. ~‘,) corresponds [Freylikher and Gredeskul 1991] to the well-known quasistationary states in the quantum mechanical sense [Landau and Lifshitz 1977]. Some simplifications may be obtained in the case of a semi-axis with respect to which the solution of the problem satisfies the boundary condition ~
(d~Idx)~= —ik(l
—
r )I( I
+
r ),
and the continuous spectrum is non-degenerate. In particular, when r = 1 (total reflection) it may be shown that for an arbitrary positive E the envelope of the state ~1i is localized within the segment [0. LI (but usually it is not normalized because of oscillating tails). Such “random” quasistationary states play an essential role in the waveguide propagation of waves in a randomly stratified medium [Freylikher and Gredeskul 1990, 19921. As we will see below, in another physical problem these states lead to creation of the so-called dark solitons (see section 2.4 below). 2.3. Transmission of wavepackets Let us now consider transmission of wavepackets through a disordered medium according to the paper by Kivshar et al. [1992].In a general form the incident wave may be presented as ~(x) =
f
~(k)
ikx
dk.
(2.29)
In this case the average transmittivity of the wavepacket is also given by an integral over the spectrum, 2~T(k)) SdkiW(k)L KT ~ = Idki~1’(k)i
(2.30)
where T(k) is the transmittivity of the plane wave with wavenumber k. We start from the case of the Schrödinger equation with a Gaussian white-noise potential (2.6), where the function (T(k)) is described by eq. (2.25), and we assume that the spectral density 1P(k)1 of the packet is taken in the form
S. A. Gredeskul, Yu. S. Ktvshar, Propagation and scattering of nonlinear waves in disordered systems
11
11P(k)12 = {cosh[IT(k k0)12’q]}2 , i.e., the packet spectrum has its centre at k
(2.31)
—
0 and its width is of order of The calculations with the help of eq. (2.30) show that the behaviour of the average transmittivity 3, where k is ofa the packet essentially depends on the dimensionless parameter ~ defined as ~ = (kik0) characteristic value of the wavenumber, built from parameters of the random potential and the system length as well as the width of the spectral density (2.5) of the packet, k = (D’qLI4ir)113 In the limit case ~ 1 the spectral density is a narrow function and the average transmittivity (2.30) of the packet differs only slightly from the “monochromatic” average transmittivity (2.25), ~.
.
.~
~
~1.
Thus, the L-dependence of the transmittivity in this case is also exponential, (T~~) —..exp[—UI4l(k~)]. But in the opposite case, i.e. for ~ > 1, this dependence becomes slower. In particular, when the inequalities
3
\1~L>>r
1>>L’, 1—kik0>>(~iV~L)” are valid, then the average packet transmittivity takes the form (T~~) (ir2i2V’~)(~Ik) exp[
—
(3irkI2~)(1
—
~_h/3)]
and, therefore, we obtain —ln(T~~3iki2i~U”3.
(2.32)
We do not claim that the analogous procedure is also valid for the transmittivity on a realization. But let us suppose that the corresponding results make sense with logarithmic accuracy. Replacing in all formulas of this subsection the quantity ~T(k)) by the quantity T(k) ~ exp(— Lii), and the quantity (T~ 0)by the quantity T~0,we find
In ~
~—(3~ki2”~~)(1
—
2”~I3~”~), ~ >1.
Thus, the packet transmittivity on a realization decreases with the growth of the disordered 113. In thesegment case of length U more slowly than the “monochromatic” transmittivity, because —ln ~ L the Helmholtz equation (2.9) the decrease of the average wavepacket transmittivity turns out to be much slower. Indeed, the integrand in the numerator of eq. (2.30) is proportional to the quantity exp(—k2Ldi4) [see eqs. (2.30) and (2.10)]. Thus, when L—~c’~, only small values kci~U”2contribute to the average transmittivity and, therefore, we obtain (TWP)HO~L”2,
i.e., the law of the decrease of (T~~) with U is not exponential but power-like.
(2.33)
S. A. Gredeskul. Yu.S. Kivshar, Propagation and scattering of nonlinear is’aves in disordered testems
12
2.4. Generation of dark solitons by a random input pulse 2.4.1. Dark solitons In this section we will study a nonlinear problem which reduces to a linear transmission problem in a disordered system. From the physical point of view the problem under consideration is the propagation of short optical pulses in nonlinear single-mode optical fibers. This propagation is described by the well-known nonlinear Schrodinger (NLS) equation [Hasegawa 1988; Hasegawa and Tappert l973a, b] which in an appropriate system of normalized coordinates takes the form i aulax
=
—~riuI~t
—
21u1u.
(2.34)
Here u is the (complex) amplitude envelope of the pulse, x is the distance along the fiber, and the time variable t is a retarded time measured in the reference frame moving along the fiber with the group velocity. The value o~= ±1in eq. (2.34) corresponds to the negative (+) or positive (—) group-velocity dispersion (GVD). In the framework of the approach based on the inverse scattering transform (see. e.g., Zakharov et al. [1980]) the solution of eq. (2.34) for u(x. 1) with the initial condition u(0,
t)
=
u(i)
(2.35)
reduces to the problem of finding scattering data for an associated linear equation. We will discuss below the case of the positive GVD, i.e., the case = 1. It is known [Zakharov and Shabat 1973] that the NLS equation with positive GVD is exactly integrable. and the spectrum of the associated problem is described by the Dirac-type equation. if
—iu~di~iIdt+ [v(t)x~ =
(~),
v(t)
=
—
—
w(t)a~]~= A~.
Re[u5(t)].
w(t)
=
(2.36)
lm[u(t)1.
o~,o~,o~being the Pauli matrices. In the case I u0(t)~1,.=
u
=
const.
(2.37)
the associated linear problem contains discrete levels. Each discrete level of the scattering problem corresponds to a dark soliton. The dark soliton is described by the expression 2. u(x.
t)
= U11
(A—i~+expZ exp(2iu~x+ i~1),
Z
=
2~u1(r
—
t~—
2Aux).
v=
—
A
and it corresponds to the eigenvalue A = A e Here 00 is the phase of the “initial” condition (2.35) (we assume that lim 11~u(0, r) = u1 e’~’).and t~,~ are arbitrary constants. However, in real optical experiments the condition (2.37) is never satisfied. Moreover, the limit values of the initial condition are always equal to zero. This means that in a rigorous sense the dark soliton does not exist at all. Nevertheless, recently dark solitons were observed by three experimental
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
13
groups [Krökel et al. 1988; Weiner et al. 1988; Grudinin et al. 1988]. In the first two of these papers the input pulse was an extended constant wave (cw) background pulse with a rather sharp and short-time perturbation. The authors of the third paper treated their experimental results assuming that a random input pulse may generate dark solitons. Let us discuss, first of all, the results of the American groups [Krökel et al. 1988; Weiner et a!. 1988]. Indeed, the associated linear problem in these two cases has no discrete spectrum at all but it has quasistationary ones (see Gredeskul Ct al. [1990], Kivshar [1992]). It means that solutions of the NLS equation cannot contain exact solitons which are stable at infinite distances, but soliton-like excitations with large (but finite) lifetime. Then the possible mechanism to generate such excitations in the case of a random input pulse [Grudinin et al. 1988] may be the following. If all the states of eq. (2.36) with the random potential (v, w) on the whole axis are localized, then the estimate (2.17) for the modulus of the reflection coefficient from the disordered segment described by the finite random input pulse, 0(1) + iw(t) with the duration T 0 is also valid in this case. (Note only that in eq. (2.17) instead of the localization length the localization time r has to be used, and in addition T0 has to replace U.) Therefore, eq. (2.28) has to possess solutions with real parts like eq. (2.27) and imaginary parts exponentially small (in the parameter T01r). These solutions may correspond to soliton-like excitations. 2.4.2. Localization properties of the Dirac equation First of all, we will prove the positiveness of the Lyapunov exponent y( A) for eq. (2.36) with a Markovian random potential on the whole axis [Gredeskul and Pastur 1991]. Using the current conservation J = (.Ji, cr5.Ji) = const., we may parameterize the state by the formula iji
with real-valued
~,
~/i
and 4. It may be shown that the Lyapunov exponent defined as
‘y(A)= lim (~(t))It
(2.38)
is given by y(E)
=
(v sin 4’
—
w cos
where the symbol (...) ~ stands for averaging with the help of the stationary distribution p(4, v, w) of the random functions 0(t), w(t) and ~(t). The phase ~(t) satisfies the equation ~ d95/dt= ucos~+ wsin~—A. If the random potentials 0(t) and w(t) are Markovian then the triad (4, v, w) forms a three-dimensional Markovian process. Its probability density p(~,u, w) satisfies the Fokker—Planck equation ap/at+(a/a4)[2(ucos4+wsin~—A)]+~~p=0, where the operator
.~
(u
(v, w)) acts only on the arguments u, w. So, for the Lyapunov exponent
14
S. A. Gredeskul, Yu.S. Kivshar. Propagation and scattering of nonlinear waves in disordered systems
(2.38) after simple calculations (see, e.g., Bratus’ et al. [1988])we obtain the result .~
p~(u)1~p~,(u)
~( A) =4 ~
p11(u)
where we have introduced the Fourier transform of the stationary distribution p(~,u).
Pm(~’)
=
Jd~ e’”~p(~,u).
Following the proof scheme described in the book by Lifshitz et al. [1988] (see pp. 136—139) it is possible to show that the Lyapunov exponent y( A) is positive with probability 1 for all values of the spectral parameter A. The positiveness of the Lyapunov exponent is only a necessary condition for the pure point character of the spectrum. It was shown by Kotani [19861that for complete localization it is sufficient if the density of states is a continuous function. For the simplest case [when u(t) and v(t) are Gaussian white-noise potentials and, therefore, they are very sharp functionsi this statement may he easily checked (see Bratus’ et al. [1988] where an even more complicated case was analyzed). 2.4.3. Quasistationary states and soliton-like excitations Now we can calculate those real values A,, of the spectral parameter A which are close to the corresponding complex values A,, describing the so-called quasistationary states [Gredeskul et al. 1990]. For simplicity, we will consider only input pulses, so that u11(t) = v(t), and w(t) = 0. In this case the embedding equation for the reflection amplitude r has a form similar to eq. (2.22) in the case of the Schrodinger equation, 2). dridt = 2iAr iv(t)( I + r It follows from this equation that, for a sufficiently large pulse duration T and sufficiently weak amplitudes of the potential u, we may write A,, = (irnIT 11)[1 + o(1)]. The difference between the last formula and the result (2.27) is connected with the fact that the spectral parameter A in the Dirac equation (2.36) 2. is analogous to the wavenumber k in the Schrödinger equation (2.2) rather than to the energy E = k As for the localization time, it can be easily calculated for the case of the Gaussian white-noise potential v(t) when —
~v(t)~= 0,
(v(t)v(t’)~= 2~(t
—
t’).
It is found that in the quasiclassical region ~i A the localization time is equal to (2.39) and to leading order it does not depend on the spectral parameter. The result (2.39) is also valid for random potentials v(t) with nonzero (but short) correlation time r~.More precisely, if the correlation function of the random potential,
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
(o(t)o(t’))
=
w~f((t—
15
t’)Ir~)
has a characteristic scale w~,then in the case ~ ~ 1, the result (2.39) is valid in the spectrum region 4 T,~. The spectral density a( A) I 2 of the generated excitations coincides with the transmittivity of the spectral problem (2.36) with the finite random potential 0(t), 0(t) 0 if t E [0, T0]. Therefore, as ~
‘~
in the case of the Schrödinger equation, we can represent the spectral density in terms of the reflection and transmission amplitudes r~and t12 (see subsection 2.2.3), 2 = t 2111 r_r+12 . (2.40) a(A)L 1t21 Furthermore, using the estimate (2.18) we can write the moduli of the reflection amplitudes p±as = e~, where the values ~ are exponentially small with respect to the parameter T 01r. Then, the spectral density (2.40) may be written in a simple form, such as 2 = 4~~6_i(ö~ + (2.41) —
~_)2.
Ia(A)L
But both exponentially small values ~ fluctuate exponentially. Therefore, the spectral density (2.41) is also, as a rule, exponentially small. The exceptional values of A~are only those ones for which the values &~and ~ coincide with exponential accuracy. Then we may obtain for these A~the estimate a( A~ ) = T( A~) -~ 1. It means that such values of the spectral parameter correspond to the resonant transparency of the realization [Lifshitz and Kirpichenkov 1979]. Thus, only those solitons are effectively excited by the random input pulse 0(t) which correspond to resonant transparent states with T( A~) 1, but not to all quasistationary states. -2
—
3. Stationary nonlinear transnussion The simplest generalization of transmission problems to the nonlinear case is to consider transmission through a finite nonlinear disordered segment. Such a formulation of the problem gives rise to the possibility of obtaining well-defined transmission characteristics with the same inter-relations as in the case of linear problems. That is the first simplification. The second one is connected with the approximation of taking into account only one temporal harmonic and neglecting the higher ones. In special cases, i.e. for the NLS equation, such a procedure is not an approximation approach, but only a choice of the class of solutions. At this point we would like to make two important remarks. First of all, one-frequency solutions cannot provide a complete analysis of the nonlinear transmission problems, but in some cases they may demonstrate only the main features of the phenomena. Secondly, one of the main questions which have to be answered when we deal with nonlinear problems is the stability of the obtained solutions, even if they are exact. We will discuss the stability in the next chapter demonstrating that some solutions analyzed in the present chapter are not stable for time-dependent problems. Before begining qualitative investigations of stationary scattering, we will briefly discuss general properties of transmission through a nonlinear disordered segment. As we saw in the previous chapter, in the linear case Anderson localization leads to exponential decrease of transmittivity. But the main reason of this effect is the inference phenomenon that is supressed in the nonlinear case due to absence of the superposition principle. Hence, it is natural to expect slower than exponential decrease of the transmittivity in nonlinear problems.
16
S. A. Gredeskul. Yu.5. Kivshar. Propagation and scattering of nonlinear waves in disordered systems
Another new property related to the nonlinear transmission is connected with the difference between two formulations of the scattering problems. The first is the fixed output problem, when the outgoing current or the outgoing intensity is fixed. The second formulation is the fixed input problem, in which it is assumed that the intensity of the input wave is fixed. In the linear limit both of these formulations are the same. But when nonlinearity is large enough, the solution of the fixed input problem can be nonunique due to the bistability phenomenon (see, e.g., Delyon et al. [1986],and also Flytzanis [1986] for a review of some applications, and references therein). The fixed output problem is simpler than the fixed input one. its properties for different kinds of randomness will be the subject of our investigations in the main part of this chapter, section 3.2. In subsection 3.2.1 we discuss analytical and numerical results obtained in the pioneer paper by Devillard and Souillard [19861, their main result being the power law of decreasing transmittivity when the length of the segment increases, T ~ L The next subsection is connected with another approach developed by Knapp et al. [1989,1991] and their results confirm the previous ones. In the last subsection the embedding approach (see, e.g., Babkin and Klyatskin [1980], Klyatskin [1988]) to the problem is presented. This approach was successfully used by Doucot and Rammal [1987a,b] for nonlinear transmission problems. Their main results are also in good agreement with the previous conclusions. The fixed input problem is considered in section 3.3. Our attention is paid to the bistability phenomenon and the influence of disorder on this effect. -‘.
3. 1. Preliminary remarks The subject under consideration in this chapter is the stationary NLS equation. = k2çb. —d2~/iidx2+ U(x)~/i aI~I2~ —
(3.1)
We assume that both nonlinearity a) and disorder U(x) are concentrated within the finite-length interval [0, L] so that outside the interval our system is described by the linear Schrödinger equation. Therefore, we can consider the transmission problem in a standard formulation (2.11) when there are an incident wave from the right, a reflected wave for x> L, and a transmitted wave for x <0, (—
fA~(e_ xL) ~i(x) {A 11t e~, — —
+
r e’”~~’1) , if x> L if x <0.
3 2 (..
The only (but essential) difference between eq. (2.11) and eq. (3.2) is that we have introduced the wave amplitude A11 because, due to nonlinearity, all transmission characteristics will depend on this value. As a result, we have well-defined scattering characteristics, i.e., the transmission and reflection amplitudes r and t, respectively. In the case of the stationary equation (3.1) there is only one variable x and we will treat it as a temporal variable for convenience. The NLS equation (3.1) leads to current conservation, J = —(1 i2i)(~j*dç~i/dx
—
~(i
d~/J*Idx) = const.
(3.3)
[This definition differs from the usual one only in the sign, because in eq. (3.2) it is assumed that the current goes to the left.] If the random potential is absent, i.e. U(x) = 0, the energy is also conserved,
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems = const. E = ~(Id~IiidxI2+ k2I~pI2+ ~a~~p~4)
17
(3.4)
In the linear case, the conservation laws (3.3), (3.4) and the form of the solution yield the evident result ti = 1 (i.e., the free solution with only an outgoing wave on the left semi-axis). But in the nonlinear case such a procedure leads to a nontrivial result. Indeed, substituting eq. (3.2) into eqs. (3.3), (3.4), we obtain for the quantity F = E/kJ + 1 the following relations: alA i21t12 2 alA 1211 + r14 4k2 =F=1T+ 4k21t12
2+
(3.5)
This means that for the nonlinear segment the reflection amplitude is not an arbitrary value bounded by the condition ri 1 only. As follows from eq. (3.5), in the fixed output problem, when the current is fixed, i.e. when J = kiA 2iti2 = const., or, equivalently, the output intensity w 01 0 is fixed, 21t12, (3.6) w0 = A01 the reflection amplitude must satisfy the equation,
i1+rl 4
8k 3 ri 2 aJ(1
—
‘
ri2) +
(1
—
ri2)2 =1.
In the second case, when the intensity w of the initial wave is fixed, w
1A 2, (3.7) 01 the reflection amplitude has to satisfy another equation, which follows from the previous one after the change J—* w. Note that the input (w) and output (w 0) intensities defined in eqs. (3.6) and (3.7), respectively, satisfy the equation =
2) w = w01(1 ri 3.2. Fixed output problem .
(3.8)
—
3.2.1. Stepwise random potential To the best of our knowledge, the first investigation of wave transmission through a nonlinear disordered segment was carried out by Devillard and Souillard [1986].In that study, wave propagation inside the segment is described by the NLS equation (3.1) with a stepwise random potential U(x). More precisely, this potential is assumed to be equal to a uniformly distributed random quantity U~,so that Umin < U(x) < Umax
18
S.A. Gredeskul. Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered sysiem.s
U
Umax
~ Urnj,v
Fig. 1. A realization of a stepwise random potential.
dzIdx
=
—(z2 + k2)
+
U(x)
+
aJIIm z,
with the boundary conditions (2.20). The real (—u) and imaginary (—v1) parts of the impedance, z = —(u + iiv), satisfy on each interval [xn, xn+i] the system of equations. duidx=u2—11v2+(k2— U~)+aJv~3~~’i3v , dvldx
=
—2uv
—3~/Ci3u.
(3.9)
(3.10)
It is easy to see from these equations that the relation (3/au) du/dx + (3/3u) dv/dx
=
0
(3.11)
is also valid, and hence the motion defined by eqs. (3.9), (3.10) conserves the measure d~= du dv. Then it follows that the dynamical system described by the equations (3.9), (3.10) is a Hamiltonian one, with Hamiltonian
2
2
1 2
~/{~(u,v)=[u +(k —U~)]v+1iv+~aJv-
(3.12)
This fact will be widely used in the next subsection. Integrating the system (3.9), (3.10) over the interval (xn+t,xn) we can obtain the connection between z~ 1= z(x~) and z,, = z(xn÷i) [note that the function z(x) is continuousi, z,, = ~J3(Zn ~ y~).The necessary number of iterations gives rise to an expression for the total impedance ZN z(U) in terms of the initial condition Z0 = z(0) = —ik and the realization parameters 1’
S. A. Gredeskul, Yu. S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
(N; y1,.
yN~ U1,.. . ,
. .
19
UN). According to eq. (2.21), the total reflection amplitude equals
r(U) = [ik + z(U)]/[ik
—
z(U)]
(3.13)
.
As is well known in the theory of disordered systems, all simple and observable results are related to the so-called self-averaging as well as mean values. But the reflection amplitude is not a self-averaged value, so that the quantity r~defined by the formula
which is similar to eq. (3.13), is strongly fluctuating. Therefore, it is more convenient to investigate another quantity with more predictable behaviour. 2 and its One of for the large main lengths problems of this subsection is to investigate T = 1 to zero ri when behaviour of the nonlinear disordered segment. Ifthe the transmittivity transmittivity tends the length tends to infinity, the expression 1 ri has the same behaviour. So, the subject of our interest now will be the arithmetic mean value of this quantity, —
—
(1- jr~~
(3.14)
Now, the first result of Devillard and Souillard [1986]may be formulated as follows. For any r >0 the value of T described by eq. (3.14) tends to zero with probability 1 when L tends to infinity. So, the transmittivity of the nonlinear disordered segment averaged according to the definition (3.14) tends to zero when the segment length is large enough. The second result is that for J small enough and L not so large the transmission decreases exponentially. This result is more or less evident because in the fixed output problem with the small current J and not so large L the nonlinear term in the stationary NLS equation can be neglected. The most interesting problem is to analyze the decreasing rate of transmittivity in the nonlinear regime. The third result of Devillard and Souillard [1986]is related to this problem, and we are going to discuss it in detail. The result may be formulated as follows. The transmittivity cannot decrease faster than L when U ~. To prove this statement let us estimate, first of all, the rate of the energy increase using the formula (3.12). As we will see in this subsection and in the next one, the main part of such an analysis is based on the fact that it is possible to find a quantity which is conserved in the absence of disorder, but is not conserved in the presence of disorder, so that its behaviour would help us to understand the transmittivity dependence when U tends to infinity in the disordered case. Using eq. (3.12) and the continuity of the wavefunction for the value E~ E(x~~ 1) we may obtain the following result: 2 E~~1 E~= ~(U~ U~+1)ic1s(n)i It follows from this equation that the estimate
-2
—~
i~,
—
—
—
E~i~
is valid, where LW
=
Umax
—
Umin.
But from eq. (3.12) it follows that k’(n)i2 ~2~1E7~,
so that
20
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
E~ + E~~ z~ U\[E~7~. If we replace this inequality by an equality and the finite difference by the derivative with respect to the corresponding variable, then after solving the resulting equation we will + E 2. As a result, the estimate obtain En = (1Ia)(~j~U)2(n 0) —
En ~(1/4a)(z~U)2(n
E 2 (3.15) 0) is valid, and this result may be rigorously proved by induction. In eq. (3.15) only the constant E 0 depends on the output current J. On the other hand, from the definition of the energy (3.12) we may obtain another estimate, 2(k2 Umax). (3.16) flA0~ Since we deal with the fixed output problem when the current J = kTIA 2 is a constant, from eqs. 0I is bounded from below by (3.15), (3.16) it follows that the transmittivity of the whole segment (n = N) the value, +
—
1 (3.17) (N+E 2 0) This is the end of the proof because N ~ U when U tends to infinity. The estimate (3.17) is mathematically correct but there is a flaw in its physical sense. Indeed, neglecting in eq. (3.12) all terms except the nonlinear one, in contrast to the case considered above, we may obtain, instead of eq. (3.17), the following result: T4aJUma~
k(z~U)2
T
(aJ/k ~U)I1
+ rNi/(n +
E 0).
(3.18)
From the mathematical point of view this inequality is not informative because its r.h.s. may be equal to zero. But it seems very natural to assume that the contribution of the realizations with small values of 1 + rNi is not so large. Thus, inequality (3.18) allows us to hope that the inverse proportionality between the transmittivity T and the segment length U is correct for large U. To conclude this subsection it is necessary to mention that the numerical calculations made by the authors of the paper cited above have shown that in the case of small current J and not too large U the exponential decrease of the transmittivity indeed takes place (fig. 2) while for sufficiently large values of the current the decrease is according to a power law (fig. 3), i.e., T ~ U -
3.2.2. Markovian potential Here we will consider the case when the potential U(x) in the stationary NLS equation (3.1) is a Markovian one. This problem was considered by Knapp et al. 119891, the Helmholtz 2~iui2u,so that who in eq.started (3.1) from we have to replace equation (2.9) with the additional nonlinear term k
2
2—
U(x)—~—ke(x), a—*kw.
(3.19)
We will deal with the thus modified version of eq. (3.2) and use the boundary conditions following from the form of the solution [see eq. (3.2)] outside the nonlinear disordered segment, i.e., ~/i
S.A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
21
10_a
+ 3
~
5x10
Umm
4++H.+~~4ALA
‘
200
A
400
600
800
N
Fig. 2. The regime of exponential decay of the transmission for small current (or small nonlinearity) and short enough length; YL = L -, logI ~t’(— L)/ !P(0)I. The + corresponds to the case al = iO-’~,the L~to aJ = 10_20. In both cases nmn = (k2 — Um=)”2 = 1.2, and nma, = (k2 — Urn 2 = 4.4 1,)” lDevillard and Souillard 1986].
(di/i/dx + ik~Ii)I~... 0 0,
=
(d~iIdx— iki~J)IX_L= —2ikA0.
(3.20)
Let us introduce now two functions q and 0 instead of the complex function ~/i(x), =
A0q(y)~t~ e’°~
102
(3.21)
io~
iO~
N
Fig. 3. The regime of slow polynomial behaviour for large current or long systems: C~stands for the Cesaro sum N’ E (1— Ir,/rol). Here aJ = 1, n,,,1, = 1.2, and nrn,,, 12 lDevillard and Souillard 1986].
5. A. Gredeskul, Yu. S. Kivshar. Propagation and scattering of nonlinear waves in disordered systems
22
where y = kx is just a renormalized coordinate. The boundary condition for the function 0(y) is chosen in the form t~e’~°~ = t, and then the function q(y) satisfies the boundary conditions at x = 0. q(0)
=
1 ,
q’(O)
(dq/dy)(0)
=
0.
(3.22)
As for the second boundary condition from Eq. (3.20), it now has the form [q’ i(q —
+
1/q)]~~ = —2i/~t~,
(3.23)
where I = kL is the dimensionless length of the segment. Substituting the solution (3.21) into the current-conservation relation and the NLS equation modified with the help of eq. (3.19), we may obtain the equations 0’=1/q2,
q”—1/q3+q(1+r(y)+f3q2)=0,
(3.24a,b)
where the new dimensionless parameter /3 describing nonlinearity is /3 w3J/k. The transmittivity T of the segment can be expressed via q(I) and q’(I) using the boundary condition (3.23) at the point y = I, T4{2
+
[(q’2 + 1/q2)
+
q2]~ (3.25)
1}~
The r.h.s. of this equation is similar to that of eq. (2.14). It is easy to show by direct calculations that, indeed, the following equations hold: 2
+
1 /q2)i~= exp[2~~(L)],
q2~ 1 = exp[2~~(L)].
(q’
In the absence of disorder, when U(x) = 0. eq. (3.24b) leads to the conservation of the dimensionless energy ~‘(q’, q), 2 + 1/q2 + q2 + ~/3q~) = 1 + ~/3. (3.26) ~(q’, q)~~(q’ This quantity is connected with the energy E defined above by E = kJ~.Using the conservation law (3.26) and the expression (3.25) for the transmittivity, we may obtain explicit expressions for the transmittivity and the input intensity w = iA 2 via the solution q(y) of eq. (3.24b) with the boundary 11 conditions (3.22), T
{l
+
~f3[1q4(I)j}~, -
w
=
(/3/~){1+
~f3[1
-
q4(I)]}.
It is not difficult to see that eq. (3.24b) is equivalent to a Hamiltonian system (p p’—0~/C/3q,
(3.27a, b) q’)
q’—0H/ap,
with the Hamiltonian ~C 0(p,q) ~f’(p,q) [see eq. (3.26)1. But in the presence of a random potential, eq. (3.24b) can be also reduced to Hamiltonian form with “time”-dependent Hamiltonian function
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
~(p, q)
=
~‘(p, q) + ~q2~(y),
~ss
e(y/k).
23
(3.28)
In this case the energy (3.26) is not conserved and expression (3.27a) for the transmittivity is no longer valid. Now, we have to substitute into eq. (3.25) the result following from (3.26), T= 2[~(p, q) + 1— ~/3q4]’.
(3.29)
The remainder of this subsection is devoted to the calculation of the asymptotic value of the energy ~ as a function of the segment length. If such a problem is solved, the energy for large U will dominate, properly describing the transmission decrease. Let us consider, first of all, the unperturbed system (i.e., without random potential) and let us rewrite the energy definition (3.26) in the form ~(p,q)=~p2+V(q), V(q)—=~(1/q2+q2-l-~f3q4). The unperturbed motion is periodic and, if we introduce action—angle variables I and q5 according to the standard definitions, q
1 .1 1= ~ p dq,
fdq 4) = —w(I) j —,
&(I)
=
-~-j,
(3.30)
then they satisfy the equations of motion dI!dy
=
0,
d4)/dy = —w(I).
(3.31)
According to eq. (3.31) the action I is conserved (as well as the energy ~‘), and the angle grows linearly. The former dynamical variables are connected with the new ones by the relations, q
=
Q(I, 4)),
p
=
{2[~’(I)
—
V(Q(I,
4)))]i/2}1/2
where the function Q(I, 4)) is defined by the second equation of (3.30) with q = Q. Note that the 4)-dependence may be found explicitly in terms of Jacobian elliptic functions. The existence of the random potential changes the equations of motion (3.31) and their new form is dlldy=r(y)ah(I,4))/a4),
d4)Idy=—o.(y)—~(y)ah(I,4))/3I,
where in accordance with eq. (3.28) we have defined h(I, 4)) = ~Q2(I, 4)), so that now the “time”dependent Hamiltonian function in terms of the action—angle variables has the form ~‘C(I,4)) = + ~‘h(I,4)). To analyze the energy (or action) evolution with increasing system length, we assume that the wavelength k’ and the correlation radius r~of the random potential e(x) are of the same order. Thus, the correlation radius of the function ~(y) is of order unity. To take into account these scales we will introduce a new “temporal” variable t = 82y, and a new dynamical variable (t) f(tI~2),where ö 1 for all values of I, ç1~ 4) T [where dT/dy = —w(I)] and T. Then, the limit case ~—+0will correspond to the scale connections, namely k1 r~ ~ U, because t will be fixed (or it will tend to f~
f:
~‘,
—
‘~
24
S. A. Gredeskul. Yu . S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
infinity) and the wavelength and the correlation radius will be negligibly small. In order to obtain a nonzero correlation function of the random potential in the limit 6 —*0, it is necessary to replace also the function ~(y) by the scaled function 6~(t). After such a procedure we obtain a system of ordinary differential equations of order 1 for the scaled functions I~,r~, and the scaled random potential 6~ plays the role of a coefficient in those equations. So, if the potential is Markovian, then the set of variables (~,I~,~ forms a four-dimensional Markov process and we can write the Kolmogorov equation for the conditional average uB of any (bounded smooth) function of these variables. Finally, in the limit case 6 —-*0 we obtain the diffusion-type equation ~
~i~)
8u/8t=Iu.
(3.32)
The operator I in eq. (3.32) has the form
~2~’
~
I D,1
ah
Iah
=j ds R(s)\~j-(~T
+
~2)
(3.33)
\ (~~+ ~2+ sw(I))/,
(3.34)
and R(s) is the correlation function of the random potential ~(y). The analysis shows (see details in the original paper by Knapp et al. [19891)that for large I the term —D11 corresponding to diffusion in the i-space dominates in the expression5,the (3.33) for thecoefficient diffusion diffusion operator. In particular, for the exponential correlation function R(s) = C e~ is D 213/a413, and, therefore, the diffusion equation (3.32) has the form, 11 = A1 3u/8t
=
(8/8I)[(A/a413)12131 0u/01.
(3.35)
To obtain time asymptotics of the action we note that if 1(0) = I and t(I, M) is the first passage time of the interval [I, M], 1< M [e.g., t(I, M) is the smallest root of the trajectory equation 1(t) = M, where 1(0) = I], then the mean value (t(I, M)) tM(I) satisfies the equation (ä/0I)[(A/a413)1213] 8tM(I) /31 = —1
,
(3.36)
and the boundary condition (3.37)
tM(M)=O.
The simplest (but not rigorous) way to proceed is to treat the function u in eq. (3.35) as a conditional average of 3(t t(1, M)), i.e., as a probability density of the first time passage. Then, multiplying eq. (3.35) by t and integrating it by parts over t from 0 to infinity, we obtain eq. (3.36). The solution of eq. (3.36) with the boundary condition (3.37) has the form —
tM(I)
= (3a ‘~3/4A)(M413
—
j413)
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
25
which corresponds to a time growth of the action I proportional to t3~, i.e., M ~ [tM(I)]314.But the analysis of the first of eqs. (3.30) shows that for large I the energy ~ is proportional to j4/3~Thus, we finally obtain the result (3.38)
~ t,
that again leads to the law T ~ L ~‘for the transmittivity if the energy 1~’is the leading term of the denominator in eq. (3.29). Numerical simulations by Knapp et al. [1989]confirm that picture (see figs. 4 and 5). In fig. 4 we can see that the mean ratio ~‘Iyhas a tendency to be stabilized and the stabilization happens earlier when the nonlinearity parameter becomes larger. On the other hand, in the case when the nonlinear parameter is sufficiently small, we obtain the result of the linear theory of disordered systems, i.e., the mean value of the energy grows exponentially (Fig. 5). 3.2.3. Weak Gaussian white noise Embedding approach. Another approach which is based on the invariant embedding method was used by Doucot and Rammal [1987a,b]. They started from the well-known embedding equation for the reflection amplitude r [Babkin and Klyatskin 1980; Klyatskin 1988] which in the case of the NLS equation (3.1) takes the form 0r(L, w)/t3L = 2ikr
—
(i/2k)[U(x)
—
awll + ri2][(1
+
r)2
+
(r
—
r*)w 3r(L, w)/t3w].
(3.39)
In the nonlinear case the reflection amplitude depends on the initial wave intensity w as well as on the system length L, and now the embedding equation involves partial derivatives. In the linear case this equation coincides with eq. (2.22). We can use now the method of characteristics (see, e.g., Courant and Hilbert [1962]) and replace eq. (3.39) by the system of ordinary differential equations,
—
0.25
~=O.25 arno.1
— — —
._,‘._.~....
~ 0.20 I.’. I..
~ 0.15
,‘
~ o.io
-
0.05 0.00~ 0
I
2000
I
4000
6000
I
8000
y Fig. 4. Stabilization to linear growth for different values of a. The curves represent the averaged energy divided by y for 1000 realizations with e = 0.7 and mean correlation length 1 = 0.5 lKnapp et al. 1989].
26
S.A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered sstems Linear (a=0)
— — -
—
2
a=10
—
io~
a=10
-
—
a=1
be V
a ~1O6.
—-
F
io~~ ~
——---
1Oo~ 0
I
I
I
2000
4000
6000
8000
Fig. 5. Comparison of energy growth with the linear case for small y. The curves mean the same as in fIg. 4.
dw/dL
=
—(w/k)[U(L)
dr/dL =2ikr— (i/2k)(1
—
awil + r121 im r.
(3.40)
+
r)2[U(L)
(3.41)
—
awll + r~2]
[cf. eq. (2.25)]. Both values, the reflection amplitude r and the intensity w, are the embedding parameters. But the input intensity w(L) is connected with the output one by the relation (3.8), i.e., w(L)
2), (3.42) 0/(1 r~ and, therefore, in the fixed output problem we have to solve the differential equation for the reflection amplitude only, =
w
—
=
2ikr
~
-
(1
+
r)2(U(L)
-
aw)) 1+rV)
(3.43)
and to use the natural boundary condition r(0) = 0. For the modulus p and the phase 0 of the reflection amplitude (r = p e’~)we obtain from eq. (3.8) the following equations: ~=_l_P~sin0(U(L)_awl+2PcosO+P2)
~dL
k
(3.44)
i—p
(3.45)
In this subsection we will consider a weak Gaussian white-noise potential like the one in eq. (2.8). So, in the main approximation it is natural to neglect disorder completely, and consider the dynamics of a
S. A. Gredeskul, Yu. S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
27
nonlinear ordered system. We know that in the absence of the potential U(U) the conservation law H0=2F=const.
(3.46)
is valid. Here the function F is defined by eq. (3.5) and the Hamiltonian H0 coincides with the Hamiltonian defined in eq. (3.12) at U = 0. In terms of p and 0 the integral of motion F is
2)2
(3.47)
4k1—p2 F= l—p~ 2awo(1+2~cosO+~
The orbits corresponding to the unperturbed motion always pass through the point p = 0 due to the initial condition r(0) = 0. They can be examined with the help of the conservation law (3.47) or of the parametrical equations (3.44), (3.45) at U(U) = 0. But from the preceding subsection we know that the random potential U(U) leads to very small values of 1 p(L) for sufficiently large U. As a result, up to large “time” U the point (p, 0) is found on an unperturbed orbit and it does not go through the zero point. Let us parametrize this orbit by employing the meaning of the transmittivity T 0 corresponding to the phase 0 = Then, due to the conservation law (3.46) all other transmittivities T(O, T0) can be found, too. Using eqs. (3.44), (3.46) and a perturbation expansion, one can find the following result of the second-order approximation: 2T (z~lnTIi~U)= —Const.’tfk 0/aw0. (3.48) —
ii~.
Here the angular brackets stand for averaging over all realizations of the random potential, and Const. in the r.h.s. of eq. (3.48) is proportional to the value D characterizing the intensity of the potential. The term i~ln T is additional to the logarithm of the transmittivity for the period of the motion along the unperturbed orbit, and it is calculated in second order of perturbation theory, t~Ubeing the length of the period. If we now express the transmittivity T0 = T(IT, T0) in the r.h.s. of eq. (3.48) via T(O, T0) and replace the finite difference in the l.h.s. of this equation by the corresponding derivatives, then eq. (3.48) becomes a differential equation describing the decrease of the transmittivity T(O, T0) according to a power law. In particular, when 0 = i~, then T0 = 2.T(~,T0) andcontrary, it follows0 from eq. (3.48) that Tfor If, on the = 0, then T 2(0, such Owe will have the relation Tmax(U) = T,~(U)ixU 0~T 0) and the result is Tmjn(U) = T0(U)x L’. Thus, the transmittivity orbital motion 7 with a power ‘y,averaged such thatover 1< yan<2. decreases according to the power law, Tav(U) m U Hamiltonian approach. In the quoted paper by Douçot and Rammal [1987a,b] another approach was also used, based on eq. (3.47) which we will write now in the form F
2/T + (awo/4k2)il
+
r12/T2.
(3.49)
This value is conserved in the absence of a random potential. But, if U(x) ~ 0, then F changes as U decreases. This change can be easily analyzed using the “time”-dependent Hamiltonian (3.12). Indeed, let us replace there u and v by the new dimensionless variables b = kv and a = k~u.Besides that, it is necessary to replace the values U~by U(L) and ~‘~(u,v) by ~l’(a,b, L)/k. Then eq. (3.12) has the form ~C(a,b, U)/k=2F(a, b)—[U(U)/k2]b,
(3.50)
28
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
2)b2. (3.51a, b) 0/2k Equations (3.51a, b) when written in appropriate notation turn out to be of a form similar to the r.h.s. of eq. (3.49). Because a and b are canonically conjugated variables, they satisfy the Hamiltonian equations 2F(a, b)
=
k1 da/dU
a2b
=
+
G(b),
0(2F)/ab
G(b) = b + 1/b
U(U)/k2,
—
+
(aw
k’ db/dL
=
—0(2F)/0a.
Therefore, F in the “time”-dependent case satisfies the equation d(2F)/dU
=
—2abU(U)/k2.
One can see from eq. (3.49) that for very large U, when the transmittivity becomes extremely small, or for a large nonlinearity parameter, i.e. for aw 2, the second term in the r.h.s. of the equation ~k dominates. So, in this case we can expect that the 0estimate for the transmittivity Tx F1 /2 will be valid. Using a perturbation expansion for F up to second order we obtain a differential equation for ~F) which is qualitatively similar to that for T(L) obtained in the previous subsection. Its solution with the boundary condition F(L = 0) = 2 + a w 2, has the form 0/4k /~F)= C
2)[1 1(aw0/k
+
C
2, 2(D/aw0)L]
where C 1 and C2 are some dimensionless constants, and D is the coefficient in the correlation function (2.6). From this result and formula (3.26) the following estimate is immediately obtained: 1.
(3.52)
Tx L Nevertheless, if the nonlinear parameter a w 2 is extremely small, then for not very small lengths L 0/k In this case the nonlinearity is negligibly small, and the first term in the r.h.s. of eq. (3.49) dominates. the transmittivity decreases exponentially, so that F grows exponentially, and the estimate lnFxL/l
(3.53)
is valid, I being the localization length. It is an interesting problem to find the crossover between linear and nonlinear regimes because for very large U we have to obtain again the power law (3.52). To solve this problem, let us note that the nonlinear part of the function F defined by (3.51a) becomes essential at b b~where b~is of the order of the inverse nonlinear parameter, b~ k2/aw 0. Therefore, a strong influence of nonlinearity on the behaviour of F will arise only in the case when 2/crw the inequality F(L) ~ b~ holds. Combining this result with eq. (3.53), we obtain the estimate U~ 1 ln(k 0) forantheexponential length U~ 2, we have where theofcrossover happens. nonlinearity, decrease transmittivity for UThus, forL aw0 k It is interesting to note that in the case of a small output power, disorder helps nonlinearity. Indeed, the fast growth of the reflectivity ri2 with increasing U leads, according to eq. (3.28), to an increase in w(L) that corresponds to an effective growth of nonlinearity. —
-~
—~
‘~
S. A. Gredeskul, Yu. S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
29
3.3. Fixed input problem. Multistability and disorder Let us discuss now the fixed input problem. We will start from the embedding equation (3.40), (3.43). A possible procedure is to solve eq. (3.43) with the fixed output power w0 and the initial condition r(0) = 0 and then to substitute that solution into eq. (3.40) to find a solution with the initial condition w0. After such a procedure we have to equate the result to the fixed input power w and find appropriate solutions for the output power w0. The result of such a procedure is expressed, in fact, by formula (3.42) with r = r(L, w0). However, the equation 2] (3.54) w = w01[1 ir(L, w0)i may have more than one solutions w 0(w) implying multistability which is a specific feature of the nonlinearity effect. As is well known (see, e.g., Delyon et al. [1986] and references therein), multistability is a nonlinear phenomenon which exists even without disorder. Following Doucot and Rammal [1987a,b], let us consider the value w(L) which is related to the initial value w(0) = w0. Such a characteristic of the scattering is a periodic function and its period L0 is a function of w0, U0 = L0(w0). If we denote as wm(wo) the maximum value of the function on the characteristic, then it is clear that the input power w satisfies the inequality —
w0 < w < wmax(wo).
(3.55)
On the other hand, this inequality means that only values w0 such that Wmmn(W) < w0 < w, where wmin is the inverse function to wmax, can lead to the input power w. There is a very simple way to show which solution of eq. (3.54) is a function of the system length L. Let us consider the characteristic (fig. 6a) with the starting point w0 satisfying the inequality (3.55), and let us cross it with the straight line w = const. Each point of intersection of this line and the characteristics (w, U) corresponds to a point with coordinates (w0, L) in fig. 6b. If we carry out this procedure for all possible values of w0, we obtain the resulting dependence w0(L) (see fig. 6b). The points corresponding, e.g., to w0 = w and w0 = Wmjn(W) belong to different characteristics and so they have different periods. As a result, the curve w0 versus L has a form similar to that presented in fig. 6b, and for L large enough there are more than one values of the output power w0 corresponding to the same input power w. In the paper by Doucot and Rammal [1987a,b] it has been shown that the degeneracy of the solutions of eq. (3.54) grows linearly with the length of the disordered segment, and that a weak random noise does not essentially change this degeneracy. In the framework of the approach developed by Knapp et al. [1989]we may search for multistability in a simpler way. Indeed, it is necessary to solve eq. (3.46) with e(y) =0 for q( y) and substitute this solution into eqs. (3.27a, b) which gives the “transmittivity versus input power” dependence in a parametric form. Using numerical procedures, Knapp et al. [1989]observed multistability and they noted that the threshold value /3~of the dimensionless decreases as the 3cr may be nonlinearity obtained for parameter large values/3 of U by means of inverse of thetheory system [Knapp length. The perturbation et al.critical 1989],value /3~/ (16I3)~3(kL)1°.However, a more important and interesting effect discussed by Knapp et al. [1989]is related to the threshold intensity for the bistability phenomenon. This threshold intensity is much smaller when randomness is present in the system. For example, when maxi ~(y) i 0.1, kL = 1000 and the dimensionless nonlinearity parameter is aw 2 01k =-
30
S. A. Gredeskul, Yu.S. Kivshar. Propagation and scattering of nonlinear waves in disordered sy~e,r~
Wmaz
wmmn(W)
(b)
L Fig. 6. (a) A characteristic starting from the point
w,,
versus L. and (b) the solutions of eq. (3.54) versus L.
0.008, the bistability phenomenon occurs for about 25% of all realizations of the random potential studied by Knapp et al. [1986]. However, in the problem without the random potential with r(y) = 0 bistability is not observed up to aw0/k2 —0.3 (the other parameters being the same). In the recent paper by Knapp et al. [19911 a new quantity characterizing multistability and localization in a nonlinear disordered segment was introduced. That quantity is the average density of solutions w 0 of eq. (3.54) as a function of transmittivity T = W0/W which the authors denote by D(T; U, w). Let N(T~,T2 L, w) be the expected number of solutions w01 of eq. (3.54) lying between T1 and T2. Then, the average density of solutions is given by D(T;L,w)=lirnN(T, T+z~T;L,w)/AT. Therefore, the total number of the solutions is given by the simple relation. N(0,1; L, w)=JdTD(T;L, w).
(3.56)
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
31
This is a measure of how likely multistability is to appear in the nonlinear case. The averaged transmittivity (T) can be found with the help of the density (3.56) as (T)=
TdTD(T;L,w)
J0tdTD(T;L,w)
J~
(3.57)
.
The averaged multiplicity density D gives a lot of information about the localization phenomenon. When there is localization, the density function will have a peak near the point T = 0 and most of the mass will be near zero. In the opposite case, when the peak is not near zero and the mass is spread out, it is no longer reasonable to say that there is localization. Note that in the nonlinear case (w > 0) the behaviour of (T) for large U cannot possibly resemble the linear case (2.25). One might think that this is true because, once the fields become small deep inside the scattering region, the linear theory takes over and, hence, we have localization. However, as U increases, more solutions are created and they are localized even at large L. All these facts are reflected in the behaviour of (T) which is both a statistical average and an average over all realizations. Its behaviour at large U is very important but until now it has remained an unsolved problem. Let us now briefly discuss numerical simulation results carried out by Knapp et al. [1991]in the diffusion (white noise) limit (see subsection 3.2.2). Figure 7 shows the function D( T; U, w) for L = 30 and w = 0 (i.e., in fact, the linear case) and, for comparison, the case at w = 0.16. The solid curve for w = 0 has a peak and most of its mass is near zero. This is an indication of localization. However, the mass for the dotted curve (at w = 0.16) is shifted from zero to the right, indicating delocalization in the nonlinear case. Notice also that the curve at w = 0.16 is always higher than that at w = 0. This is because in the nonlinear case there may be multiple solutions. The average multiplicity of solutions is depicted in fig. 8, which shows the total multiplicity N(0, 1; U, w) as a function of the input intensity w for U = 7.5, 15.0, 22.5 and 30.0. These curves appear to be increasing exponentially fast and can roughly be interpreted as predicting the likelihood of multistability [e.g., if N(0, 1; U, w) = 5.0, a typical realization would have five possible solutions, or 20.0
I
I
w=0.0 w=0.16 15.0
i~10.0
0.00.~
.‘
i..,..
~2~4
.....~
L.~,..
0.8
1.0
Fig. 7. Average density of solutions, D(T; L, w), as a function of transmittivity, T, for w = 0 and w = 0.16 for L = 30.0. The density was compiled using 5000 realizations [Knappet al. 1991].
32
S.A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems 15.0
I
—~
I
o *~-~-*
L=7.5 ~ L=15.0 0 L=22.5 L=30.0
0.0 0.00
0.04
0.08
0.12
0.16
w Fig. 8. Average multiplicity of solutions, N(0, 1; L, w). as a function of w. The averages are for 5000 realizations [Knapp et al. 1991}.
output states, for given L and WI. Figure 9 shows the weighted averaged transmittivity (T) as a function of w for U = 7.5, 15.0, 22.5 and 30.0. Note that for the larger lengths where localization was being felt in the linear case, these curves are increasing as functions of the input intensity w. This is further evidence for delocalization for nonlinear wave propagation. Finally, fig. 10 provides a direct comparison of the decay of the weighted averaged transmittivity for the linear and nonlinear cases (here the scaling parameter 6 from subsection 3.2.2 differs from zero). The solid curve represents an exponential decay of the linear localization. For small L, (T) is nearly identical for the linear and nonlinear cases. In this regime nonlinear effects are not strong and there is no multistability. However, above a certain value of U, the decay slows down for the nonlinear case and it eventually appears to stop. This indicates that there are a significant number of states which are delocalized. 1.0
I •
I
0.8
0.6 A
-
—.
a
~
A
4
4
0
~
0.4
D———-~
a-
p
0.2
0.0 0.00
~
~p ~-*
•
0.04
I
0.08
L=7 .5 L=15.0 L=22.5 L=30.0
,.
0.12
0.16
w Fig. 9. Weighted average transmittivity. (T), as a function of a. The averages are for 5000 realizations lKnapp et al. 1991].
S. A. Gredeskul, Yu. S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems 1.0
0.8
•
I
N,,
.
---
0.00
10.00
33
20.00
30.00
w=0.00 w=0.04 w=0.08 w=0.12 w=0.16
40.00
L
Fig. 10. Weighted average transmittivity, (T), as a function of L, for different values of w. The averages are for 1000 realizations with e = 0.25 [Knappet at. 1991].
The results of Knapp et al. [1991]presented above give strong evidence that there exist delocalized transmission states for nonlinear wave propagation through a random medium. However, since only a one-harmonic approximation was used, not all these states are physical due to the effect of modulational instability. We will discuss the problem in the next chapter.
4. Modulational instability and wave propagation 4.1. Modulational instability in homogeneous models Equation (3.1) was introduced by Devillard and Souillard [1986]and Doucot and Rammal [1987a,b] as a stationary version of the time-dependent NLS equation
i4)1+4)~~+ai4)I24)r(x)4).
(4.1)
Let us now look for solutions of eq. (4.1) in the form 4)(x, t) = exp(—ik2t)4i(x). Then, substituting this expression into eq. (4.1) yields eq. (3.1) for the function ~i(x).However, unlike the linear case, solutions of nonlinear problems have to be investigated for stability properties. As is well known, the homogeneous NLS equation (e = 0) with a >0 exhibits modulational instability of a nonlinear plane wave [Benjamin1967; Benjamin and Feir 1967]. To describe this phenomenon, let us consider the exact solution of eq. (4.1) with e(x) = 0, 4)(x,t)=Aexp(iO),
O—kx—.wt,
(4.2)
for which the following dispersion relation holds: ü.=k2—aiAi2.
(4.3)
S. A. Gredeskul, Yu. S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
34
To analyze the stability of the solution (4.2), (4.3), we have to derive equations for small deviations from its exact form. Substituting 0 = 0 + p(x, t) and A = A + a(x, t), we may obtain a system of the linear equations for the functions p(x, t) and a(x, t), ~
+
2aA2a
—
Ap,
—
2kAp~= 0,
a, + Ap 11
+
2ka~= 0.
(4.4a, b)
If we look for solutions of eq. (4.4b) in the form (a, p) = (a0, p0) exp(iqx iI’lt), 12 and q being the frequency and wave number of small periodic perturbations, then we find the following dispersion relation: —
2 = q2(q2 2alAi2). (4.5) (12 2kq) For a > 0 eq. (4.5) has a negative solution for Q2 which corresponds to an instability. For example, at k = 0 the region of instability, 0< q2 <2aIAl2, has the maximum growth rate = aiA~2. In a discrete chain the modulational instability may be affected by discreteness [Kivshar and Peyrard 1992]. To demonstrate the main difference between the continuous and discrete cases, let us consider the discrete version of the NLS equation [cf. eq. (4.1)], —
—
~
(4.6)
which has the exact solution 4~(t)=Aexp(ikan—iwt),
(4.7)
for which the following dispersion relation is valid: w =4sin2(~ka) aiAl2.
(4.8)
—
Looking for evolution of linear waves excited around the solution (4.9), we may compute the dispersion relation for the wave number Q and frequency Ii of such linear waves, [12 —2 sin(ka) sin(Qa)]2
=
4 cos(ka) sin2( ~Qa) [4cos(ka) sin2( Qa) ~
—
2a Al2].
(4.9)
As may be seen from eq. (4.9), at k = 0 the wave (4.7) will be unstable for any Q provided the wave amplitude exceeds the critical value A~r= 2/a. Thus, in some cases even a regular periodicity of the lattice may enhance modulational instability; then it is quite clear that this effect may be stronger in a disordered system. 4.2. Influence of inhomogeneities on modulational instability In the case of inhomogeneous systems the problem of modulational instability is more complicated. To understand the main features of modulational instability in inhomogeneous media, let us consider the simplest case of a single point impurity (see fig. 11). The impurity generates a reflected wave which interacts with the incident wave. Therefore, in the region before the impurity position we have two
S.A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
0 Fig. 11. Generation of a counterpropagating wave due to the wave scattering by an impurity (x scattering and reflected waves enhances modulational instability.
35
x =
0). In the region x <0 interaction between the
coupled nonlinear waves. Modulational instability of two coupled waves has been extensively investigated in nonlinear optics and plasma physics (see, e.g., McKinstrie and Luther [1990]for a review and references therein). For example, for two counterpropagating waves with amplitudes A and B of the same order the maximum growth rate is (see, e.g., Barday and Remoissenet [1991]), (4.10) ~a{iAi2+ iBi2 + [(IAi2+ iBi2)2 + 6iAi2iBi2]112} , i.e., it is always larger than the value cr~= alA!2 for a homogeneous case, and in the case A! = iBi =
eq. (4.10) yields the result crcr = 3o~.In a general case, the effect of wave scattering includes three simple waves (see fig. 11), but the main properties come from the interaction between the incident and reflected waves. So, the arguments presented above show that modulational instability is always enhanced by inhomogeneities and this phenomenon has a simple physical explanation. As an example of this, we would like to present some numerical results by Peyrard and Bishop [1990].They considered a system of rotators given by the Hamiltonian (in fact, the discrete sine—Gordon model) (4.11) where the coupling constant between rotators S(n) was changed in 2.5 times in 10 central cells (the total chain contained 200 cells). Figure 12 depicts the behaviour of a plane wave in such an inhomogeneous case, the rest of the parameters being the same as in the homogeneous case. It is clear that a small disorder generates modulational instability much faster (cf. scales in fig. 12). The modulational instability of nonlinear plane waves suggests that the results of Devillard and Souillard [1986]and Doucot and Rammal [1987a,b] (see also the previous section) for eq. (3.1) have a questionable applicability for the time-dependent case, at least for a >0 in eq. (4.1). To analyze the scattering in the framework of the time-dependent NLS equation, Caputo et al. [1989]have performed numerically scattering experiments in the full time-dependent, discrete NLS equation, iA~
+
A~~ 2(A~+ + 1 A~_1 2A~+ aiA~i +1 A~_1)= V~A~ , —
(4.12)
the potential V,, being equal to zero for n <0. As is well known, the discrete Ablowitz—Ladik model (4.12) is exactly integrable at V, = 0, so that solitons may propagate in the discrete chain without distortions [Ablowitz and Ladik 1976]. As follows from the numerical simulations, there is not much difference between the cases a = 0 (linear) and a = —10 (nonlinear medium without modulational instability), and the scattering of the waves looks very similar (fig. 13). The field reaches a stationary state, and its values are comparable for the two cases. On the contrary, for a = +10 (fig. 14), the field
36
S. A. Gredeskul, Yu. S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
5
I I
I
I
I
I
I
I
I
I
I
500.0 a=1.0 ip= 50
t=
~
I
0
C
a1 O•
1
_~
!!!!1tlh1hui~1h!hhu!!h!1u1u!u~!!ut1
—
0
10 I
0
~p’~ 5~•
I
I I
I
I
50
I
I
50
I
I
00
I
100
50
I
150
I
200
I
200
cells Fig. 12. Base rotation ~, along an inhomogeneous chain at t = 1500 for an initial wave amplitude of 1.0. For comparison, the homogeneous case at the same parameters but at t = 2000 is also presented [Peyrard and Bishop. 1989].
does not decay as fast as for cases mentioned above and, additionally, it does not tend towards a stationary state; there are some pulses which are created and propagate through the medium, the latter being the solitons of the NLS equation (see fig. 15). Hence, from the comparison between the two cases a = —10 and a = +10 it is clear that these coherent structures (in fact, solitons) are responsible for the increase in energy transport for the nonlinear case. Therefore, solitons play an important role in nonlinear wave scattering and that is why it is necessary to study the soliton scattering in detail.
-ø.løøe+Ø1
[og(]UI) .
-ø.500e+ø1 ø.øøøe+øø Fig. 13. The scattering experiment for a nonlinear wave,
f3 =
.....~.
ø.l5øe+04 n
___
—10. Shown is log
ø.300e+04
10~Ujversus n for a given time ICaputo et al. 1989].
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
-ø.løøe+ø1
-ø .
37
I
3øøe~øl
log(
[UI
—ø.500e+ø1 ________ I ø.øøøe+øø ø.150e+04 n Fig. 14. Same as fig. 13 for
0.255e+0ø
$
=
ø.300e+04
10 [Caputoet al. 1989].
I
-
ø.150e+04 n
ø.300e+04
0. 1 27e+0O
[UI
ø.øOøe+O0 ø.300e+ø1 Fig. 15. The scattering of a nonlinear wave at f3
= +
10. Shown is [U,]versus n for a given time [Caputo et al. 1989].
5. Soliton propagation through disordered media We will analyze the propagation of solitons through disordered media separately for different types of solitons. Our main goal is to point out the qualitative difference between three types of solitons scattering in a disordered medium. In all problems related to the soliton propagation we will assume that the medium is nonlinear so that the soliton is created before the scattering, but disorder is considered to be concentrated in a finite region.
38
5. A. Gredeskul, Yu. S. Kivshar. Propagation and scattering of nonlinear waves in disordered systems
5.1. Dynamical solitons 5.1.1. Models One of the systems where dynamical solitons may propagate is a chain of particles with nonlinear nearest-neighbour interactions (see, e.g., Flytzanis et al. [1985]). The model can be described by the Hamiltonian (5.1)
where y~(t)is the displacement of the nth particle from its equilibrium position, = dy~/dtand m~are the velocity and mass of the nth particle, respectively. In the simplest case the interaction potential V(r) is assumed to be polynomial, ~P
V(r)
~Gr2+ ~Ar3 + ~Br4.
=
(5.2)
The equation of motion of the nth particle is given by the expression, ~
—V’(y~ y~
=
—
(5.3)
1) + V’(y~~1 — y,).
In the continuum limit, eq. (5.3) can be transformed into a generalized Boussinesq (Bq) equation which takes into account isotopic disorder (see Kivshar et al. [1991b]), 2u~ —p(u~)~ q(u~) f(x)u,, c 1 hurxx2 =0, (5.4) —
—
—
where a is the lattice spacing, u
2
2
y~and
3
4
—
c ~Ga/m, p~Aa/m, q~Ba/m, h=Ga/12m (5.5) are constants; the function f(x) describes inhomogeneities [in the homogeneous case, when m0 = m, we have f(x) = 1]. In the absence of any perturbation and when there is only a cubic nonlinear term (q = 0) (the standard Bq equation), eq. (5.4) has a kink-like solution of the form u(x, t)
Vt)/2L], 2 c2), L \/h/(V2 Am [(2 sgn h)/p]\/h(V V being the kink velocity; the kink energy is given by =
Am tanh[(x
=
—
Ek
(5.6) (5.7)
—
—
c2),
= j~g(mA~,/aU)(4V2 + c2).
(5.8)
On the other hand, when only the quartic nonlinear term is included (p the form u(x, t)
=
Am tani[(x
—
Vt)/L],
Am
=
±2\/i~7~,
=
0), the soliton solution has
(5.9)
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
39
where L is given by the same expression as above. The total energy of the soliton (5.9) is Ek
=
~(mA~/aL)(2V2
c2).
—
After some transformations, it is possible to derive from eq. (5.4) the Korteweg—de Vries (KdV) equation (see, e.g., Kivshar et al. [1991b]),which is integrable by inverse scattering transform, and it is useful to apply perturbation theory for the soliton transmission. Another example of the system (5.1) is the Toda lattice, when V(r)
=
(a/b)[exp(—br)
—
1] + ar,
(5.10)
for which the discrete equation (5.3), (5.10) becomes exactly integrable [Toda 1967]. 5.1.2. Soliton transmission coefficient
Let us assume that a disordered segment embedded between two homogeneous semi-infinite chains has randomly distributed impurity masses m = -ym 0 with the probability p in the disordered part. A soliton is incident on the disordered segment, and we are interested in the transmission coefficient, which is the transmitted energy E1 over the soliton energy E1~.This problem was studied numerically by Li et al. [1988b] for the potential (5.2). Some results for a quartic nonlinear disordered chain are presented in fig. 16. The behaviour of the soliton may be analytically studied in two limit cases [Li et al. 1988b]: the independent scattering approach, and the linear limit. The independent scattering approach is valid in the dilute case, when we assume that the concentration p of the mass impurities is low, i.e., the mean distance between two nearest-neighbour impurities is larger than the soliton size. In this case we can write T=JJ T~(E~),
(5.11)
where N is the number of impurities inside the disordered segment, E1 is the incident kink energy for
£
0.2
£ 0
£ ~
Io~oISnSrqY
£
0.1
£
0.05
100 Nj Fig. 16. Transmission coefficient of the dynamical soliton with the incident energy E,, the quartic nonlinearity [Li et aI. 1988b].
1000 =
0.66 as a function of the disordered segment length N, for
40
S.A. Gredeskul, Yu.S. Kivshar. Propagation and scattering of nonlinear waves in disordered systems
the ith impurity, and T.(E,) is the transmittivity for the ith impurity. The transmitted soliton for the ith impurity is then the incident one for the (i + 1)th impurity. So, we may write the relation E,~1= E~T~(E~) .
(5.12)
Multiplying these equations for i = 1,2, , N and denoting T= EN±i/El, we again obtain eq. (5.11). The recurrence relation can be written in the form .
=
E1~1 E,
=
—
—E1[1
—
.
T,(E,)]
.
—E1R1(E)
(5.13)
,
where R1(E,) is the reflectivity for a single impurity, which has been previously calculated numerically [Li et al. 1988a] and analytically [Li et al. 1988a; Kivshar et al. 1991b], and it is given by R(E)
=
a0E~(y
—
1)2
(5.14)
,
where a0 ~0.66, a = 2 for the quartic potential (p = 0), and a0 = 8.0, a = 2/3 for the cubic potential (q = 0). Since on the average there are (~x)pimpurities in the interval between x and x + 6x, from eq. (5.13) a differential equation can be written down, taking the continuum limit dE/di~= —pE(x)R[E(x)].
(5.15)
Note that eq. (5.15) is an approximate one; moreover after direct averaging eq. (5.13) leads to another equation, d(E)/dx= —p((E)
—
(ET(E)~),
which is not the same as eq. (5.15); the latter may be obtained only under the assumption that the following decoupling is valid: ~ET(E) ~ = ~E) T( ( E)). Such a decoupling means that for macroscopic scales the energy E is not a fluctuating quantity. After integrating eq. (5.15) with the help of eq. (5.14) we obtain T(U)
E(U)/E1~= (1
+
c0L~~,
C0 =
paa0E~(y
—
1)2.
(5.16)
Therefore, for lengths L ~ 1/c0 the independent scattering limit gives for the quartic potential (p = 0) T(L) U /2 and for the cubic potential (q = 0) T(L) U -3/2 i.e., this dependence is implicitly dependent on the type of nonlinearity through the reflection coefficient. The result T(U) L i/2 is in agreement with numerical simulations for the quartic potential (see fig. 16). As was mentioned by Li et al. [1988b],the same behaviour follows from the averaged transmission coefficient (2.35) with the localization length 2 cot2(~k), (5.17) A(k)= 2(lyyp) (y—l)p(l—p) -
-
i.e., formula (2.35) with eq. (5.17) yields [Li et al. 1988b] [cf. eq. (2.38)],
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering ofnonlinear waves in disordered systems
41
T—
(l—p+yp) 1 L~1 518 (y l)[p(l p)]”2 VL’ So, we have seen that the independent scattering picture gives T U h/2 for a quartic potential (p = 0) and T(L) -~L -3/2 for a cubic potential (q = 0). However, for large enough distances, when nonlinear terms can be neglected, the total transmission must decay as L 1/2 i.e., much slower than exponential for T in the case of a plane-wave. For the quartic potential these dependences are the same. Numerical results for the disordered lattice with the Toda potential (5.10) have been obtained by Ishiwata et al. [1990](see also the paper by Tsuboi and Toyama [1991],where equivalent solitons were analyzed). The authors were interested in the same problem considering mass disorder in a nonlinear chain of particles interacting via the Toda exponential potential (5.10). The main result of the paper by Ishiwata et al. [1990] is the empirical dependence of the soliton amplitude A~at the nth particle, -
-
-
—
-
=
A 01(1 + a\fA~n)” ,
(5.19)
where p = 1—1.2 and a 0.009, which was established numerically. As we can see, the dependence (5.19) is also power-like and very close to the one previously discussed. Therefore, the transmission coefficient of a dynamical soliton decays much slower than the exponential law for a plane wave, in fact, according to a power law. On the other hand, similar dependences arise in the linear theory of disordered systems for wavepackets; in particular, for the quartic interparticle interaction potential they are the same. So, we may conclude that dynamical solitons in disordered systems propagate in the same fashion as nonlinear wavepackets for which the function T(U) is modified due to nonlinearity, but, in general, this dependence is power-like, as is the case for the averaged transmittivity of wavepackets in the linear theory. To close this section, we would like to point out that similar results, i.e., power-like decay of the soliton amplitude (which is connected with the soliton energy), follow also from the analysis of Wadati [1983] and Wadati and Akutsu [1984] (see also discussions of the results and other examples in the review of Bass et al. [1988]) for an exactly solvable stochastic KdV equation. However, the physical sense of the averaged soliton shape calculated by Wadati [1983]and Wadati and Akutsu [1984]is unclear when the results are compared with either the independent scattering approach results or direct numerical simulations, where a unique realization is actually observed. 5.2. Envelope solitons 5.2.1. Physical models We will start our discussions of envelope soliton propagation from two physical examples in which these solitons naturally arise. First of all, let us consider the nonlinear Klein—Gordon model for the wave variable U(x, t), 3=o, (5.20)
u,~—u+~u—.i~u
w 0 being the frequency in the linear limit. Equation (5.20) arises in a number of different physical problems (see, e.g., Sanchez and Vázquez [1991]and references therm) and it may be obtained also as an expansion of the widely applied sine—Gordon equation when U ~ 1. In the framework of eq. (5.20),
42
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
we will consider the effects of changes of the natural frequency, in the form of random point impurities with intensities and random positions x,,, E~
+ r(x),
r(x)~ ~
s~6(x xv). —
With this change, eq. (5.20) takes the form (5.21) U,,— U1~+o4~U—f3U3=—s(x)U, and it may be thought of as a perturbed nonlinear Klein—Gordon equation. To study the properties of nonlinear wavepackets as described by eq. (5.21), we will derive an equation for the envelope u(x, t), assuming U(x,t)=u(x,t)e
101)2
*
+u (x,t)e
loll) (
(5.22)
,
the asterisk standing for the complex conjugation. In this way, we are separating the system dynamics into relatively fast oscillations with linear frequency w 0, and a slowly varying nonlinear envelope u(x, t). Substitution of eq. (5.22) into eq. (5.21) leads, in the lowest order of the oscillations, to the equation 2u = —~(x)u. (5.23) + 2iw0u, u~1 3/3!u! If we assume a slowly varying u, the inequality u,, w 0u, is valid and, therefore, eq. (5.23) may be transformed into the well known NLS equation for the dimensionless envelope, —
—
-~
+
+ 2I~!2c~ =
r=—t/2w0,
r(x)q~,
(5.24)
4~=\/~/3u,
(5.25)
where r(x) is the same as in eq. (5.21). For another example in which a NLS equation can be derived, we consider vibrations of a one-dimensional lattice composed of atoms with atomic masses m0 in which each atom interacts only with its nearest neighbours. Let u~(t)be the displacement of the nth atom from its equilibrium position and let k2 and k4 be nearest-neighbour harmonic and quartic anharmonic force constants, respectively. The potential energy W of the chain is taken to be 2+ ~k W= ~k2 ~ (u0÷1 u~) —
4.
(5.26)
4~ (u~~1u~) —
If the chain has isotropic impurities with the masses m~,then the equations of motion are given by m~ü~ = k
3J (5.27) 2(u~~ + u~ 2u~)+ k4[(u~+ u~)3 + (u,, u,2 ) Without the isotopic impurities, the nonlinear chain (5.27) was analyzed in a number of papers (see, e.g., Takeno et al. [1988],Burlakov et al. [1990a,b], Kiselev and Rupasov [19901,Takeno and Hon [19901),and it is the model for which the intrinsic localized modes were first discussed by Sievers and Takeno [1988]. .
i
—
—
—
S.A. Gredeskul. Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
43
We will assume that the disorder is weak, i.e., the relative value (m~ m)Im is small enough, and we will apply the so-called quasi-continuum approximation. To do this, first of all we note that the linear spectrum of the homogeneous chain (m~= m) has a cut-off frequency t0m = 4k2/m and, besides, —
at this frequency a particle in the chain oscillates out of phase with its nearest neighbours. Therefore, we may look for solutions with w 2> w in the form ~,,
u~(t)= ~(—1)°V’~(t) &olmt + c.c.
(5.28)
,
where c.c. stands for the complex conjugate quantity. Assuming that the function 0I’,~(t)2is slowly varying may in time, i.e. (~P~)~ ~as ~a function 111, and that in the vicinity of the cut-off w~,,,it be considered of the continuous argument na frequency, x, we maywhere derivew from eq. (5.27) the following NLS equation: —
-~
2~Wm11ç + (a2k
2~P= _w~,F(x)lP. (5.29) 7/m)1I’~~ + (12k4/m)[~1’I The function F takes into account the isotopic disorder, and, in particular, in the case of a single impurity of mass M this function takes the form F(x) = a[(M m) /m]6(x). Equation (5.29) describes in a rigorous way localized oscillations in the chain when their characteristic length is large enough in comparison with the lattice spacing a. After an appropriate scaling we may get from eq. (5.29) the NLS equation (5.24) exactly. Equation (5.24) is the exactly integrable NLS equation [Zakharov and Shabat 1972], and, aside from the application considered here, it is important in its own right, for it has to do with a number of problems related to solid state physics (see, e.g., Kivshar and Malomed [1989]).The function e(x) may then represent a structural disorder of the associated system. The most remarkable feature of the homogeneous nonlinear system described by eq. (5.24) at s(x) = 0 is that it allows the distortionless propagation of localized excitations in the form of envelope solitons, —
exp[+~iVx i(~V2 a2)T] —
u 5(x,r)=a
cosh[a(x
—
—
Vr)]
(5.30)
,
where a and V are the soliton amplitude and velocity, respectively. These are the objects we are interested in, and whose propagation along the disordered medium we intend to describe in the following subsection. 5.2.2. Independent scattering approach
We next move to the problem of NLS soliton scattering by a random system of point impurities with equal intensities e. The soliton is incident on the disordered layer form the left, scatters, and decomposes into reflected (r) and transmitted (t) parts. We assume that, after passing through each impurity, the wavepacket may be described again as a soliton plus some small-amplitude waves (radiation). To begin with our approach to the calculation, we recall that the NLS system is characterized by two integrals of motion, namely the energy E and the “number of quasiparticles” N defined by E= j1 dx(!u~!2 +e(x)tu! 2 —lul 4 ),
N=
j1 dx!uI 2
.
(5.31a,b)
44
S. A. Gredeskul, Vu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
During the scattering the values of these quantities are changed due to the emission of radiation, but they are conserved for the total system. So, we describe the scattering process by means of two quantities: the total energy transmission coefficient T~ = E,/E1, that is, the transmitted energy E~over the incident one E1, and the “number of particles” transmission coefficient T~ = N1!N1, where N1 and N1 are defined analogously. It is to be understood that the constraints E = E, + Er = const. and N. = N, + Nr = const. must hold. When the concentration p of impurities is low, it is quite reasonable to suppose that the average distance between two neighbouring impurities is larger than the soliton size. In this limit we may treat the scattering by many impurities independently, defining a transmission coefficient for the whole layer as T~H. T1, T1 being the transmission coefficient of the jth impurity. In this so-called independent scattering approach, back scattering interference effects are neglected, and hence the transmitted soliton for the jth impurity is then the incident one for the + 1)th scatterer, and we can write (cf. the equations of the previous section) (j
E1~,=E1T~(E1,N1), E1~, E1
=
—
=
N1~,=N~.T~(E1,N,),
—E1R~(E1,N1),
(5.32a,b) (5.33)
N1+, N1 = —N1R~(E1,tv~,), (5.34) 1~ = 1 T~’~ stand for the energy reflection coefficient and the “number of particles” where thecoefficient. R~” reflection It is very important to note that, unlike the linear case, in the nonlinear problem we have to distinguish two independent characteristics of the scattering, e.g., two reflection coefficients. Since on average there are (i~x)pimpurities in the interval ~x, we can derive from eqs. (5.33) and (5.34) the following differential equations: =
—
—
dE/dx = ~pE(x)R~[E(x), dN/dx
=
N(x)J,
(5.35)
—pN(x)R~[E(x), N(x)],
(5.36)
assuming that during the scattering the soliton parameters have small changes. In order to deal with the system (5.35), (5.36) it is now necessary to compute the exact soliton reflection coefficients R~(E,N) and R~(E,N) defined for a single inhomogeneity. 5.2.3. Evolution of soliton parameters Soliton reflection coefficients may be calculated via soliton perturbation theory (see e.g., Kivshar and Malomed [1989]). The reflection coefficients are defined by the reflected wavepackets in the form of linear waves emitted by the soliton during the scattering. The details may be found in the papers by Kivshar et al. [1987a, bJ and Kosevich [1990], and here we will give only final results, R(N)=
6~~JdyF(y,a), 2[(~/4a)(y2+a2 F(y, a)
=
cosh
R~ ~~Jdyy2F(y,a),
-
1)1’
a ~N/V.
(5.37a,b) (5.38)
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
45
As a consequence, eqs. (5.35) and (5.36) give rise to the following system of integro-differential equations [Kivshar et al. 1990]:
=
4JdyF(y,a),
(5.39)
~
(5.40)
where the distance is measured in units of x0 = 64/irpe2, i.e., z = xIx0. In the remainder of this subsection, we describe the variation of N and V along the z axis. By carefully studying the system (5.39), (5.40), it turns out that a is the very relevant parameter of the problem. This quantity is physically significant, because it is related to the nonlinearity of the incoming soliton; the soliton will be “more nonlinear”, so to speak, for large values of a. This can be understood by realizing that for each given value of V(0), the greater a is, the larger the number of quasi-particles contained in the soliton becomes, and, in addition, its spatial extent becomes narrower; this is obvious from formula (5.30). On the contrary, if a is small, it can be easily seen that the wave looks very similar to a linear wavepacket: Once we have clarified the physical meaning of a, we begin the analysis of the two coupled integro-differential equations (5.39) and (5.40). We first consider the linear limit, a 4 1; in this simplification, the system (5.39), (5.40) can be solved analytically in a quite approximate fashion, yielding (note that if a 4 1 the derivative of V becomes negligible, for the integrals are never much greater than unity) V(x) = V(0) = const., and therefore, the two transmission coefficients become the same, (5.41)
40,
T~~(x) N(x)IN(0) = E(x)/E(0) = e_X/ A 2(0) /pe2 = 1/pR, , R~ e2IV2(0), 0 V
(5.42)
where R, is the reflection coefficient of a single impurity. This result demonstrates the exponential decay of the transmission coefficient. It is a consequence of the quasilinear features of solitons with small a, that closely resemble those of linear wavepackets, for which A(k 0) = A0, and k0 = V/2 has the sense of a carrier wavenumber of the packet [see eq. (5.30)]. In the opposite limit, a — 1 or above, the system of integro-differential equations (5.39) and (5.40) has been solved numerically, and the possible evolution of the system can be described as follows: the dynamics depends essentially on the value of the parameter a(0) = N(0) /V(0), because the system given by (5.39) and (5.40) has a fixed point for a, namely, the solution of the transcendental equation 2 —2 + G(a~)= 0, a~
G(a) =
_____________
,~ ~
~
.‘
.
(5.43)
~y,a,
This value was computed both by numerically solving this equation, and by direct integration of the system (5.39) and (5.40), getting a perfect agreement between the two results, and concluding that = 1.285 05(4) is an unstable fixed point.
46
S.A. Gredeskul, Yu.S. Kivshar. Propagation and scattering of nonlinear waves in disordered systems
A straightforward calculation of the derivative of a(z) is the only requirement to establish that only two kinds of behaviour are allowed for the system (5.39), (5.40). (i) a(0) < a0. The system evolves to a final state in which N tends asymptotically to zero while V goes to a constant, nonzero value, and hence a(c~s)= 0, as imposed by the above mentioned considerations. This behaviour corresponds to a decay of the transmission coefficient, and, as we saw earlier, for small enough a this decay is exponential and is quite correctly described by formula (5.41). If 1 a < a~the decay of N(z), though still present, is more complicated; it is slow at first, and after an interval, it undergoes a crossover to a fast exponential decrease, faster than the one than takes place when a 4 1. (ii) a(0)> a0. This kind of initial condition leads to solutions such that both N and V suffer a very rapid decay at first, up to a point at which the two functions become almost constant. Obviously, the same happens to a, that tends to some value around a 10 when z is large enough. However, the general property of these curves is the same; both transmission coefficients tend to certain asymptotic values, and, henceforth, the rearranged wavepacket is transmitted without any reflection, i.e., localization is inhibited. For different sets of N(0) and V(0) the final dependences when they give rise to equal values of a(0) are almost the same (see pictures in the original paper by Kivshar et al. [1990,1991]). 5.3. Topological solitons (kinks) 5.3.1. Effective equation of motion The last type of solitons is that of topological solitons or kinks. Such a soliton exists only in a nonlinear system with two (or more) equivalent ground states. The simplest example is the sine— Gordon (SG) model, u,~ u~+ sin u —
=
ef(x) sin u
(5.43)
,
which in the homogeneous case, i.e. when r = 0, has a kink solution of the form, uk(x,t)—4tan{exp[(xX)/V1V]}
,
(5.44)
X = Vt being the kink coordinate. As is known (see, e.g., Kivshar and Malomed [1989]), for a number of cases the kink emission is exponentially small, so that the kink dynamics may be studied in the so-called collective-coordinate approach. In the framework of this approach the evolution of the kink’s coordinate X is described by a simple equation of classical mechanics (see, e.g., Currie et al. [1977],McLaughlin and Scott [19781).To derive such an equation, we note that the unperturbed SG system has an infinite number of quantities that are conserved in the evolution, among which the momentum, P~—JdXUtUx.
(5.45)
For the kink2)”2. of eq. (5.44), eq. (5.45) takes the the formmomentum of the well-known In the presence of perturbations, is no longerrelativistic conserved;expression using eq. P = 8V1(1 V —
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
47
(5.43) it is possible to show that it varies according to dP
I
J
e
=
dx f(x)(cos u)1,
provided the boundary condition u —~0(mod 2ir) holds at x —~ ±co. The adiabatic approach is now defined by the assumption that, for e small enough, the kink shape will not be affected and only its coordinate X and its velocity V will become slowly varying functions of time. this hypothesis it 2 4 1, the kink center obeysWithin the following evolution can be then shown that, in the nonrelativistic limit V law: d2X
f(x) sinh(x
f
—
X)
—
d
8~d~~r=4ej dx cosh3(x—X) =—~-~U(X),
U(X)
—2e
J
dx
-=
f(x) cosh (x X)
,
(5.46)
(5.47)
—
where we have used the approximate expression P 8 dXldt, valid for small velocities. Thus, in the framework of the adiabatic approach, the motion of the SG kink can be thought of as that of a particle with mass 8 in the external potential U(X) defined by eq. (5.47). The same can be shown for relativistic kinks [Bergmann et al. 1983]. The following two cases arise naturally from eqs. (5.46) and (5.47). If f(x) varies rapidly over distances of the order of the kink length, then e has to be small for our approximation to hold. For example, in the case f(x) = ~(x), we have [McLaughlinand Scott 1978] U(X) = —2e sech2X. On the other hand, if f(x) varies slowly, i.e., its characteristic length (say L) is much larger than that of the kink, it is not necessary for e to be small, because all the parameters of perturbation theory are of the order of L 1, and we are left with U(X)
—2e
J
dx f((x + X)/L)sech2x~s4gf(XIL).
5.3.2. Kink transmission; a general stati.stical approach In the simplest situation, the function f(x) may be taken to be a sum of delta functions, i.e., ef(x)
=
~
s,,5(x
—
an),
(5.48)
where the numbers e,, and a~are chosen at random, and it is assumed that the distances b~= any, are identically distributed random variables with probability density p(b)
=
—
a~
b~1exp(—b/b 0).
(5.49)
Then, the equation for the kink coordinate X takes the form, 2X/dt2 8d
=
—dU(X)/dX,
(5.50)
48
S. A. Gredeskul, Vu. S. Kivshar. Propagation and scattering of nonlinear waves in disordered systems
U(X) = ~ u~(X),
u~(X) u(X
—
an),
u(X)
—2r sech2X.
(5.51)
Here, as above, we have approximated P=8V=8dX/dt, which is valid for V2 41. Thus, in the collective coordinate framework, the motion of the SG kink can be interpreted as the motion of a nonrelativistic particle with unit mass in an effective, random potential defined in eq. (5.51). In the paper by Gredeskul et al. [1992] eqs. (5.50), (5.51) have been analyzed for the case when disorder appears as randomly distributed point impurities with equal intensities, i.e., for e,~= Since the method of the Fokker—Planck equation — which is valid for time-dependent random perturbations studied, e.g., in the papers by Bass et a!. [1988] and Pascual and Vazquez [1985] (see also Biller and Petruccione [1990]) cannot be directly applied to the problem [note that the potential (5.51) is not a Markov one], a statistical procedure has been developed to compute the mean characteristics of the kink propagation. Here we will briefly summarize the main results and ideas of that study. To simplify the computation of the kink statistics, let us consider the potential U(X) as a finite one, with a width Li, i.e., ~.
—
U(X)lzOVx,
x!>Li/2.
(5.52)
Now, we assume that the distances y,, bn Li = a~ ~, a~ Li > 0 are identically distributed random quantities with the probability density [cf. eq. (5.49)] p(y) = y~’exp(—y/y 0). Here y0 has the meaning of the mean length of the empty interval between the effective scatterers (5.52). From a rigorous point of view, such a statistics of the random potential differs from the exact one which follows from eq. (5.49). But in the low-concentration limit, when y0 ~‘ Li, this difference is negligible [unless of course one is not specially interested in the effects induced by overlapping of the potentials U(X)]. Then, provided the above mentioned assumption holds, the relation y0 b0 is also valid. We are interested in the calculation of the mean values (and their dispersions) of the kink parameters, namely the mean coordinate and velocity of the kink at each time instant t. To this end, it is actually much more convenient to describe the process in terms of the temporal, rather than the spatial, distribution of the scatterers; i.e., strictly speaking, to study the problem by means of the time intervals r,, elapsed when the kink is between the nth and (n + 1)th scatterers (n = 1, 2, . . .). Each one of these intervals is separated by equal times with duration —
—
15/2
f
—
—
15/2
dX
~(X)
f
dx
-~/2 ~V~-2U(X)’
where X(X) is the local kink velocity. The time t0 has the meaning of the time it takes the kink to be transmitted through one scatterer. With respect to the time intervals ‘Tn~ they are equal to y~/v0,and subsequently they have the probability density p(T) = r0 exp(—T/r0), where r0 y0/v0. Therefore, instead of the sequence of random spatial intervals of the form Li, y,, Li, y2, Li, y3,. , Li, y,,, Li, ~ , we will consider the corresponding sequence t0, r1, t0, T2, t0, T3, , t0, T~,t0, ~ The next point to take into account is that there are configurations of scatterers of two different types. Configurations of the first kind (hereafter referred to as configurations I) are defined in such a way that the interval [0, t] includes n full intervals T1,. , T~,where the kink moves in the region where the potential is zero, i.e., freely, The probability of configurations I can be calculated as .
.
.
.
.
.
.
.
.
.
.
.
S.A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
p~(t)= (6(r’) ~
—
T’)),
t — nt~
—
~‘
49
(5.53a, b)
~,
where 0(T) is the Heaviside step function, and the angular brackets in eq. (5.53a) stand for the average over the n + 1 random variables T1,. . . , T~~ Equation (5 .53a) has a simple physical interpretation. The first 0 function means that n scatterers are contained in the interval [0, X], where X = X(t), and the second one accounts for the fact there are no other scatterers in this interval. After some straightforward calculations it is possible to obtain for P~(t)the exact result, P~(t)= (1/n!)[(t
—
nt0)/T0]~ exp[—(t
—
nt0)I’r0]
n = 0, 1,
,
. . . ,
N,
(5.54)
where N [tIt0], [.“] standing for the integer part. The configurations of the second type (we will term them configurations II from now on) consist of n full intervals, like configurations I, but the last t’ is not enough to cross the time duration of one scatterer. The probability density p,, t’) for such configurations, with fixed n and t’, can be shown to bep~(t~t’)= (~(T’ — t’)O(t’)), with the same as in eq. (5.53b). The appearance of a delta function in this expression is connnected with the fact that, at moment t, the kink has moved during a time t’ under the influence of the last [the (n + 1)th] scatterer, whereas the theta function ensures that this time interval t’ is positive. Again, after the calculation of the corresponding integral, we obtain the result (tI
~‘
p~(t~t’) = (lIn!)[(t
—
t’
—
nt~)/;]~ exp[—(t
—
t’
—
nt~)/r~] .
(5.55)
Concerning the normalization of both probability densities, (5.54) and (5.55), they have to satisfy to
J
P~(t)+ ~
1~
J
dt’ p~(t~t’) + dt’ pN(tlt’) =1,
(5.56)
where t~ t~{t/t~},{x) standing for the fractional part of x. It is easy to check directly that condition (5.56) is fulfilled. Now, we come to the next step. It is clear that, following the same line of reasoning, any dynamical variable F(t) can be considered either as a function of n, i.e., F~(t),for configurations I, or as a function of n and t’, F~ (tI t’), for configurations II. For instance, the following relationships are valid for the kink coordinate X(t): X~(t)= v0t n~0(for configurations I) and X~(t~t’) = v~t n~0 (for configurations II), where we have introduced the notation vet’ x’, ~ v0t~— Li, and t’ is connected with xl by —
—
~‘
—
~‘
—
—4/2+x’ t’=
I I
J
—15/2
dx X(X) .
so tJ—4t t’(x’) is the part of the scatterer length (in the temporal scale) that is contained in the whole motion interval [0, t], or [0, x(t)]. Analogously, for the kink velocity, the expressions v~(t)= v~(for configurations I) and v~(tIt’)= v(t’) u[x’(t’)] (for configurations II) hold. The above mentioned results allow us now to present the mean value of any time-dependent function, say F(t), as
50
S. A. Gredeskul, Vu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems I~
(F(t))
=
~ F~ (t)P0 (I) + n=()
~
f
dt’
~
(t~t’)p,, (t~t’) +
f
dt’
FN(tj t’ )p\(t~t’),
(5.57)
where N N(t) = [t/t~]and t0 is defined above. The procedure described in this subsection to obtain the mean kink parameters is quite general, and it may be applied to a number of different soliton-bearing systems. Finally, we would like to point out that to present a global picture of kink propagation in disordered media, it is necessary to compare the results obtained by means of collective coordinates to direct numerical simulations of the model. This numerical work is currently in progress. 5.3.3. Resonant kink—impurity interactions As has been demonstrated above, the kink—impurity interaction may be described by a simple picture where the inhomogeneity gives rise to an effective potential for an effective particle (kink). However, the model of a classical particle is valid only in the case when the impurity cannot support an impurity mode, a local oscillating state at the impurity site (see, e.g., Braun and Kivshar [1991] and references therein). This impurity mode may be excited due to the kink scattering and it may change the result of the kink transmission. Recently, the importance of impurity modes in the kink—impurity interactions has been pointed out in the papers by Fraggis et al. [1989], Kivshar et a!. [1991a], and Fei et al. [1991a,b]. In particular, Kivshar et al. [1991a]have found that a kink may be totally reflected by an attractive impurity due to resonance energy exchange between the kink translational mode and the impurity mode. This effect is quite similar to the resonance phenomena in kink—antikink collisions in some nonlinear Klein—Gordon equations (see Campbell et al. [1983, 1986], Peyrard and Campbell [1983])and it drastically differs from the previous theoretical approaches which took into account only radiative losses (Malomed [1985], Kivshar and Malomed [1989]). To demonstrate the origin of the resonant interactions, we will start from the SG model (5.43) which includes a local point-like impurity, i.e. f(x) = 5(x). When the impurity is absent, the model (5.43) supports the topological kink (5.44) propagating without distortion. In the presence of the s-impurity the potential (5.47) takes the form 2X, (5.58) U(X) = —2e/cosh i.e. for e >0 the impurity gives rise to attraction. According to the result of Malomed [1985],the kink may be trapped by such an attractive impurity only due to radiative losses, and there exists a threshold velocity (exponentially small in r) V,h 0(e), such that if the kink initial velocity is larger than ~hr(E), it will pass through the impurity and escape to infinity, otherwise it will be trapped by the impurity. In this consideration the reflection of the kink is impossible. Kivshar et al. [1991a]have studied the kink—impurity interactions for r >0 by numerical simulations. They found that there are three different regions of initial kink incoming velocity, namely, region of passing, of capture, and of reflection (see fig. 17). They also found that a critical velocity V0 0.2678 (for r = 0.7) exists, such that if the incoming velocity of the kink is larger than 1/c, the kink will pass the impurity mnelastically and escape without change of direction, losing part of its kinetic energy through radiation and excitation of an impurity mode. In this case, there is a2 linear relationship between the = a(V~ V~),a ~0.887 being squares of the kink initial velocity V and its final velocity V~:V, —
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems 8
51
0.14
~rTT
~.°°
50
100
150
200
&25
T
0.26
0.27
0.28
0.29
INITIAL KINK VELOCITY
Fig. 17. The kink coordinate X(t) versus time for initial velocities V situated in three different regions: the region of passing (solid line,
Fig. 18. Final kink velocity as a function of the initial velocity (e = 0.7). Zero final velocity means that the kink is captured by the
V
1 = 0.268), of capture (dotted line, V = 0.257), and of reflection (V = 0.255) lKivshar et al. 1991a].
impurity lKivshar et al. 1991a].
constant. If the incoming velocity of the kink is smaller than V~,the kink cannot escape to infinity from the impurity after the first interaction, but will stop at a certain distance and return (due to the attractive force of the impurity) to interact with the impurity again. For most of the velocities, the kink will lose energy again in the second interaction and finally get trapped by the impurity (see fig. 17). However, for some special incoming velocities, the kink may escape in a direction opposite to the incident one after the second interaction, i.e., the kink may be totally reflected by the impurity (see fig. 17). The reflection is possible only if the kink initial velocity is situated in some resonance windows. By numerical simulation, Kivshar et al. [1991a]have found a number of such windows and the results for = 0.7 are presented in fig. 18. Using the idea of the resonant energy exchange between the kink translational mode and the impurity mode, it is possible to predict the centers of the resonance windows [Kivshar et al. 1991a], -
=
V~— 11.0153/(nT + 0.3)2,
n
=
2,3,.
(5.59)
. . ,
where T is the period of the impurity mode oscillation, and V~is the critical velocity. This formula has been shown to provide very good predictions (see the table in the original paper by Kivshar et a!. [1991a]). Following the paper by Kivshar et al. [1991a],we will give a brief analysis of resonance structures in the kink—impurity interactions. The main initial point is observation that the nonlinear system (5.43) supports a localized mode. By linearizing eq. (5.43), the shape of the impurity mode can be found analytically to be 2, (5.60) uim(x, t) = a(t) e~jxj/2 , a(t) = a0 cos(L2t + 0~), 12 = ~s where 11 is the frequency of the impurity mode, and 00 is an initial phase. As a matter of fact, the impurity mode (5.60) can be considered as a small-amplitude breather trapped by the impurity, with energy Eim = Q2a~Ie. Now we may analyze the kink—impurity interactions through the collective-coordinate method taking into account two dynamical variables, namely the kink coordinate X(t) [see eq. (5.44)] and the amplitude of the impurity mode oscillation a(t) [see eq. (5.60)]. Substituting the ansatz —
=
Uk + Uim =
4 tan~exp[x
—
X(t)] + a(t) e~~2
(5.61)
S. A. Gredeskul, Vu. S. Kivshar. Propagation and scattering of nonlinear waves in disordered systems
52
into the Lagrangian of the system, and assuming that a and e are small enough so that the higher-order terms can be neglected, it is possible to derive the following (reduced) effective Lagrangian: L = 4X2 + (lIe)(à2 — u22a2) — U(X) — aF(X), where U(X) is given by eq. (5.58), and F(X) two dynamical variables are 8X+U’(X)+aF’(X)=O,
=
(5.62)
—2r tanh X/cosh X. The equations of motion for the
(5.63)
ä+122a+ ~rF(X)=0.
The system (5.63) describes a particle (kink) with coordinate X(t) and mass 8 placed in an attractive potential U(X) (r>0), and “weakly” coupled with a harmonic oscillator a(t) (the impurity mode). Here we say “weakly” because the coupling term aF(X) is of order 0(r) and it falls off exponentially. The system (5.63) is a generalization of the well-known equation 8X= —U’(X) describing the kink—impurity interactions in the adiabatic approximation. The dynamical system (5.63) can describe all features of the kink—impurity interactions. Firstly, it may be used to calculate the threshold velocity of kink capture, which is given by the relation [Fei et a!. 1991a, b] ~rr sinh[QZ(Vhr)I2V~hr] V0hr
= \/‘~
CO5h(Q1T/2V~hr)
-
, Z(V)
=
cos
I
2
[(2V
2 —
s)/(2V + r)]
.
(5.64)
Comparing the analytical results with the direct numerical simulations of eq. (5.43), Fei et al. [1991a] found that the perturbation theory is valid only for very small e (r ~ 0.05), while formula (5.64) gives good estimates of Vhr(E) for r over the region (0.2, 0.7). As was pointed out by Kivshar et al. [1991a] and Fei et a!. [1991a],eqs. (5.63) can be used as a qualitative model to explain the mechanism of resonant energy exchange between the particle and the oscillator. The resonant reflection of the particle by the potential well corresponds to the reflection of the kink by an attractive impurity. Therefore, the collective-coordinate approach can give a qualitative explanation of the resonance effects in the kink—impurity interactions. It is important to note that similar resonance phenomena have been observed in the kink—impurity interactions in the 4~’model, -
+
[1- r8(x)](-~+ ~3)
=
0.
However, the resonant structures in the ~ kink—impurity interactions are more complicated than in the SG model because the 4 kink has an internal mode which also can be considered as an effective oscillator. Fei et al. [1991a]have developed a collective-coordinate approach taking into account three dynamical variables, and they have found that due to the joint effect of the impurity and the kink internal-mode oscillation, some resonance windows may disappear. The resonant interactions described above have been analyzed numerically and analytically for a single impurity, and we can say that the physical mechanism of this effect has been well understood. It is clear that in the case of a few impurities the well-defined energy-exchange process is likely to be destroyed, especially for the case of a random lattice of impurities. However, the possibility of exciting impurity modes during the kink propagation will lead to the additional and, as we have seen for a single
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering ofnonlinear waves in disordered systems
53
impurity, more important energy loss of the kink. So, a possible mechanism of the kink decay in disordered media is the excitation of localized mode vibrations due to impurities but not radiation of small-amplitude waves. This mechanism was mentioned in the earlier paper by Tsurui [1973] who discussed soliton propagation in the nonlinear lattice with isotopic disorder. However, to the best of our knowledge, there are no more results in this direction, so that it is still an open problem to be explored in the near future. 5.3.4. Interference effects The collective-coordinate approach considering the kink as a particle does not take into account the wave origin of solitons. The wave properties of solitons manifest themselves in radiation emitted due to scattering by impurities. Interference phenomena might be the most remarkable effects for studying the wave properties of solitons. Interference may arise when the characteristic distance between impurities becomes commensurable with the characteristic length of the linear waves emitted by the soliton. As a consequence, the simplest way to observe interference is to consider a system with two point impurities [Kivshar et al. 1987b, 1991b]. Solitons emit a rather wide spectrum of linear waves, and the problem is far more difficult than the scattering of a monochromatic plane wave. As a matter of fact, the above mentioned condition that the two characteristic lengths must be commensurate has to be valid for an averaged spectral structure of the soliton emission. When a soliton has an internal frequency (like, e.g., an envelope soliton), resonant scattering is naturally expected [Kivshar et al. 1987b], while for topological kinks this phenomenon does not seem to be possible. However, resonant scattering has been predicted [Kivshar et al. 1991b] for a slowly moving SG kink when the spectral density of the emission generated by it has a narrow maximum; interference effects should appear as oscillations of the soliton reflection coefficient as function of the distance between impurities. The simplest model to describe interference effects is given by the equation (5.65) which takes into account two point-like impurities. Because the homogeneous SG equation is integrable, some theoretical analysis of the problem (5.65) is possible through perturbation theory for solitons based on 1ST [Kaup 1976; Karpman and Maslov 1977, 1978; Kivshar and Malomed 1989]. This analysis has been recently carried out by Kivshar et al. [1991b],who were interested in the influence of the parameter D on the scattering properties of SG kinks, or more precisely, in their reflection coefficient. This coefficient, R, may be defined as R E~/Ek,where Ek = 81(1 v2)”2 is the energy of a SG kink with velocity v far away from the inhomogeneous region, and ~ is the energy the kink emits backwards as radiation, due to its interaction with an impurity. The emitted energy ~ and, consequently, the reflection coefficient R can be calculated analytically by means of the aforementioned 1ST perturbation theory, the main restriction of it being the necessary assumption that the kink velocity is not changed during the scattering (the so-called Born approximation). By this means, the emitted spectral density can be shown [Kivshar et al. 1991b] to be given by —
e
2{(D12v)[kv 2(k)
=
[k
2 2
e 1(k)
=
—
~o(k)]},
(5.66)
4e,(k) cos
—~
8v
(1
—
v )
—
2[irVl cosh
—
v2w(k)12v]
,
(5.67)
54
5. A. Gredeskul, Vu. S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
where w(k) (1 + k2)”2, and n = 1,2 stand for the case when one and two impurities are present in the system, respectively. Having these expressions in mind, it is straightforward to obtain the corresponding reflection coefficients, which turn out to be R~=W1_v2Jdkr~(_k),
(5.68)
where n = 1, 2 stands for the case with one or two impurities, respectively. It can be seen from eqs. (5.66) and (5.67) that for small v, v2 4 1, the spectral density r,(k) has a single maximum at k = 0, with quite a narrow peak of width of order 2vlir -= v. As a consequence, such maximum will provide the main contribution to the emitted energy and should give rise to a resonant dependence. On the other hand, if the velocity is large, there are two maxima at ±km,km 2vRr(1 v2)”2 (1 v2)”2, and the function r,(k) is not exponentially small in the region IkI < km. Hence, after averaging over all wavenumbers there is no leading contribution, the oscillatory dependence disappears, and resonant scattering is not to be expected. It must be noted that, of course, the results for v2 4 1 make sense only from a theoretical viewpoint, because due to total reflection the kink velocity cannot be less than a certain threshold, Uthr below wicich the kink is reflected by the impurity. This is the fundamental reason for the necessity of numerical simulations. The numerical results for two values of the initial kink velocity in the SG mode! (5.65) are shown in figs. 19 and 20. The shown ratio R 212R, is a suitable quantity to use in order to search for interference effects. Indeed, when the impurities are far from each other, i.e., when D ~, the reflection coefficient has to coincide with 2R,, and the above mentioned quotient must go to unity; any difference from unity may be treated as coming from interference. In both figures, the full lines correspond to the analytically computed dependence given by eqs. (5.66)—(5.68). The main differences arise at v = 0.4 (fig. 20); in particular, the asymptotic behaviour differs between the analytical and the numerical curves. This disagreement, that in principle should not occur, can be explained in a ~iaturalway [Kivshar et a!. 1991b]. Recall that the analytical results for two impurities were obtain~dunder the assumption that the kink does not change its velocity during the —
—
—~
Fig. 19. Dependence of the reflection coefficient on the distance between impurities for a kink with v = 0.4 when the impunty is repulsive (e = 0.1, empty circles) or attractive (e = —0.1, full circles). The so)id line is an analytical result [Kivsharet al. 1991b1.
Fig. 20. The same as in fig. 19, but for v = 0.9. Only reflection coefficients for an attractive impurity are shown [Kivshar et al. i991b1.
S. A. Gredeskul, Yu. S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
55
scattering (the so-called Born approximation). However, as a matter of fact, after the first scattering the kink loses some part of its kinetic energy, so that it interacts with the second impurity at a smaller velocity, say v — z~v.Hence, the ratio R212R, does not go to 1 when the distance between impurities goes to infinity and, as was shown by Kivshar et al. [1991b],these arguments allow us to explain the difference between the theory and numerical data. 5.3.5. Two-component models In the previous subsections we considered one-component models in which topological kinks may exist. However, there are a lot of systems where the wave field supporting a kink excitation is coupled with other degrees of freedom. A very simple (and well-known) example of such models is a hydrogen-bonded (HB) chain. The main idea to model the HB chain is that protons move in double wells due to hydrogen bonds with a heavy-ion lattice (oxygen lattice) which is deformable [Antonchenko et a!. 1983; Zolotaryuk 1986]. The local distortions of the oxygen lattice can lower the activation barrier for the protons and thus promote their motion. In order to describe this effect, models must include two coupled degrees of freedom at each lattice site, so that the proton sublattice supports topological solitions (kinks) while, in the absence of coupling, the equations for the oxygen sublattice are linear ones. The different models proposed (see, e.g., Antonchenko et a!. [1983],Zolotaryuk et a!. [1984], Zolotaryuk [1986], Peyrard et al. [1987], Hochstrasser et a!. [1988], Pnevmatikos [1988], Halding and Londahl [1988]) seem to provide an effective description of the proton mobility in hydrogen-bonded chains, and they may also play an important role in interpreting certain biological processes. To describe the proton mobility in hydrogen-bonded systems, it is also necessary to estimate the influence of inhomogeneities which always exist in real systems and may sufficiently change the steady-state motion of solitons. Propagating in an inhomogeneous medium, a soliton emits radiation in the form of linear waves. This radiation is usually small for one-component soliton-bearing systems because the on-site potential supporting the topological solitons (kinks) yields a gap in the spectrum of linear excitations. However, in the case of a two-component HB system, emission of a kink moving with a variable velocity may be not small. It is easy for it to generate linear waves in the oxygen subsystem, because the latter has no gap. Therefore, the kink emission in HB (as well other two-component) models has to play an important role in the proton mobility in real systems. To describe this new type of effects analytically, we consider the model which consists of a diatomic chain of proton and oxygen atoms (see, e.g., fig. 1 of the paper by Peyrard et al. [1987]). The Hamiltonian consists of three parts. The proton part is given by 2+ ~ —p,,)2 + U(p,,)], (5.69) H~= ~ [~m(dp~Idt) .
with the on-site potential U(p) = U 2]2. The proton coordinate p,, is the displacement of the hydrogen from the middle of the0[1bond(pip0) in the static case. The oxygen Hamiltonian is written as ~
—
2 + ~kq(q,,÷, q,,)2], (5.70) 0 = ~ [~M(dq,,/dt) where q,, describes the displacement from the equilibrium position. The interaction Hamiltonian is H
H,~,=
—
~
(q,,+,
—
q,,)[1
—
(p,,1p 2]. 0)
(5.71)
56
S.A. Gredeskul, Vu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
The physical content of this interaction is the lowering of the double potential barrier due to the oxygen displacement. In the continuum limit, which is valid for the long-wave approximation, when the variable na z is considered as continuous, the Hamiltonian (5.69)—(5.71) yields the equations of motion u where 0
2) —
u~ u(1 —
—
=
—auv~,
v
u
(5.72a, b)
2 v~= a/3uu~, 0
—
s
2/m, s~= kqa2/M s~=s~/c~ , c~= k~a a = (aq 2, w~= (U~/mp~), 0W0w0iU0c0), f3 = (m/M)(p01q0) and a is the lattice spacing. In eqs. (5.72a, b) we have used the dimensional variables, T u =pip
0
,
v
=
q/q0
,
=
w
0t, (oi0/c0)x. The parameter a stands for an effective coupling between the hydrogen and oxygen sublattices. In the absence of impurities, the homogeneous equations (5.72a, b) support the steady-state propagation of the kink in the hydrogen subsystem together with its “shadow” in the oxygen subsystem, z
=
u(x, t) = ±tanh[~/l(V)], 4 = x Vt, 2 V2) cosh2[~/l(V)], dv(fl/d~= af312(s where the parameter 1(V) given by
(5.73) (5.74)
—
—
12(V) = 2(1 V2)(V2 ~ (V2—s2+ a /3)
(5.75)
—
has the meaning of the soliton width. A simple analysis demonstrates that the kink (5.73)—(5.75) exists provided 0
or s2
(5.76a,b)
and it can be proven that in the region (5.76a) the kink is a stable solution [Zolotaryuk et al. 1984]. The standard mode! for the hydrogen-bonded chain briefly presented above assumes homogeneous parameters along the chain and, as a result, the motion of soliton excitations is possible in the steady-state regime. Let us investigate the influence of a local inhomogeneity of the oxygen subsystem on the soliton motion. To obtain equations of motion in the inhomogeneous case, we assume that in the vicinity of the impurity the parameters of the oxygen subsystem are changing, so that M—* M
+
a ~M 5(z)
=
MEl + EM~(x)l,
k~ kq —~
+
a ~kq 6(z)
=
kq[l
+
Eq6(x)I,
a(w
EM
01c0)(AM/M) ,
Eq
a(wo/c0)(z~kqikq)
When the influence of the impurity is taken into account, the equations of motion become u~, u~ u(1 —
—
—
u 2 )= —auv~,
(5.77)
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems = af3uu~ sM3(x)vl, + s2eq[6(x)vx]x s2v~~ —
—
57
(5.78)
.
The more general case which describes the change of all the parameters of the model, i.e. the hydrogen as well as oxygen subsystems, has been considered by Kivshar [1991a]. In the HB chain two subsystems are coupled, and the oxygen subsystem has no gap. It means that radiation-induced effects in the two-component hydrogen-bonded chains play a much more important role than in one-component chains. For example, the threshold velocity of the radiation-induced capture of the soliton by a local attractive inhomogeneity is not small in the two-component models [Kivshar 1991a]. To demonstrate some features of radiation-induced effects in the model, we will present the results for the spectral density and the total energy emitted by the soliton scattering by the local impurity. In the zeroth approximation we have to consider emission generated by the solution (5.73), (5.74) with = x Vt, assuming, for simplicity, the case of a slowly moving kink. The energy of waves emitted in the oxygen subsystem may be presented as —
~
Eem=
~(k)
J
dx(v~+s2v~)= dk~(k),
J
(5.79)
J f
dx F(x, t) e~~~t) 2
(5.80)
dt
=
The result (5.79), (5.80) provides a basis to calculate radiation-induced losses of the soliton due to emission in the oxygen subsystem. For V2 ~2 a 2/3 the simple calculations give the spectral density in the form [Kivshar 1991a] ~
~(k)
kl ~ (EMV+es)2 = 122/3(0)
—
l~ l(V=0).
(5.81)
The total emitted energy may be obtained integrating eq. (5.81) over all k (V >0), Eem
=
i~g1T(a2f3V/s5lO)(rMV+
(5.82)
eqs)2.
As can be seen from eq. (5.82), the emission is determined by the soliton velocity and it is not so small as for one-component models. Let us consider now a disordered HB chain with local impurities. In the independent scattering approach, when interference effects are ignored, we can write the relation E. = E 1 T.(E,), where T,(E,) is the soliton transmission coefficient of the ith impurity. Then, as in the case of dynamical (nontopological) solitons, +,
=
E~4., E1 —
=
—E1[1
—
T,(E1)] =
where R.(E1) is the reflection coefficient. Since on average there are (í.~x)pimpurities in the interval between x and x + i~x,p being the impurity concentration, the following differential equation can be derived:
58
S.A. Gredeskul, Vu.S. Kivshar. Propagation and scattering of nonlinear waves in disordered systems
dE/dx
=
—pE(x)R[E(x)}.
(5.83)
The soliton reflection coefficient R(E) may be defined as a ratio of the energy, emitted by the soliton during its scattering in the direction opposite to its propagation direction. In the case when the soliton emits radiation in another (oxygen) subsystem, we may use for the reflected energy the total energy emitted by the soliton in this subsystem. Then the quantity R(x)E has the meaning of the emitted energy Eem defined above as the soliton radiation-induced losses during its scattering by the inhomogeneity. If, instead of E, we put into eq. (5.83) the kinetic energy of the soliton defined as Ek = 2m~fXO, where meff is the soliton effective mass, then eq. (5.83) gives rise to the law of radiation-induced damping of the soliton velocity in the case when impurities are installed only into the oxygen subsystem, ~
dV/dx = —pB(EMV
2
+ Eqs)
,
B = ~ ir(cr213/s51 0).
(5.84)
The solution of eq. (5.84) with the initial condition V(0) = V0 has the form r MO V—esC’x q /1
+ i...X)
4
2
C=~7Tpema/3
EMO V+rs
5
s l(,mCff
.
(5.85)
Equation (5.85) allows us to estimate the propagation length L of the soliton in the disordered system, L
-~-EMVOIEq5C~
(5.86)
where C is defined above. It is interesting that, according to eq. (5.86), L tends to infinity for EM i.e., for a pure isotopic disorder.
—~
0,
6. Concluding remarks We have briefly reviewed the problems related to propagation of nonlinear waves through disordered media for stationary and nonstationary regimes. In the case when the nonlinear medium shows modulational instability, the nonlinear constant-amplitude wave decays into a set of solitons. So, nonlinearity leads to an improvement of the transmission only when it contributes to the creation of soliton pulses. Solitons are more robust than linear waves or wavepackets and that is why their creation acts against localization effects giving a nonzero contribution to the transmission coefficient outside of the localization region. Moreover, above a certain threshold in amplitude, the envelope soliton described by the Schrödinger equation has very low radiative losses, and it propagates almost undistorted. Using the simplest independent-scattering approach, we have pointed out the difference among three types of solitons scattering in disordered media. We have shown that topological solitons (kinks) are more robust objects than other solitons, and in typical situations they may be described in the framework of a suitable collective-coordinate approach. As can be concluded from the paper, one of the main questions in the study of the influence of nonlinearity on disorder is modulational instability, which may generate solitons to improve transmission. We believe that the study of modulational instability in inhomogeneous media will improve our knowledge about nonlinear effects in disordered systems. It seems to us that another possible line interesting for subsequent studies in this field is direct numerical simulations of different soliton-bearing
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
59
systems, including disorder, and their comparison, because the propagation of solitons through a disordered medium essentially depends on their type. One of the important problems in soliton propagation is to understand the role and manifestation of a competition between the soliton length and the correlation radius of disorder. In any case, the study of soliton propagation is very important because solitons give the most relevant contribution to the transmission characteristics of a disordered nonlinear system and, in particular, they may increase substantially the energy transmission in nonlinear systems.
Acknowledgements It is a great pleasure for us to thank all the people with whom we have discussed and worked on this subject for the last years. We would like to name especially Leonid Pastur (Institute for Low Temperature Physics and Engineering, Kharkov), Zhang Fei, Angel Sanchez and Luis Vázquez (Universidad Complutense, Madrid), Alan Bishop and David Campbell (Center for Nonlinear Studies, Los Alamos), Michel Peyrard (University of Dijon), and Stephanos Pnevmatikos (Research Center of Crete) (it is very regrettable that he tragically died while we were writing this review). We are deeply indebted to Alexej Maradudin for his current and stimulating interest in the work. We thank Angel Sanchez for critical reading of the review, and Oksana Chubykalo for her help in the preparation of the manuscript. We also thank all our colleagues for permission to use figures from their original papers, and Jean Guy Caputo for original figures from his paper. This work was partially supported by the Isreaeli Ministry of Science and Development.
References Abdullaev, F.Kh., AR. Bishop and St. Pnevmatikos, eds, 1991, Springer Proceedings in Physics: Nonlinearity with Disorder (Springer, Berlin), in press. Ablowitz, M.J. and J.F. Ladik, 1976, J. Math. Phys. 17, 1011. Antonchenko, V.Ya., A.S. Davydov and A.V. Zolotaryuk, 1983, Phys. Status Solidi b 115, 631. Antzygina, TN., L.A. Pastur and V.A. Slyusarev, 1980, Fiz. Nizk. Temp. 6,3. Babkin, G.I. and V.!. Klyatskin, 1980, Soy. Phys — JETP 52, 416. Barday D. and M. Remoissenet, 1991, Phys. Rev. B 43, 7297. Bass, F.G., Yu.S. Kivshar and V.V. Konotop, 1987, Zh. Eksp. Teor. Fiz. 92, 430 lSov. Phys. —JETP 65, 2451. Bass, F.G., Yu.S. Kivshar, V.V. Konotop and Yu.A. Sinitsyn, 1988, Phys. Rep. 157, 63. Benjamin, TB., 1967, Proc. R. Soc. London Ser. A 299, 59. Benjamin, T.B. and J.E. Feir, 1967, J. Fluid Mech. 27, 417. Bergmann, D.J., E. Ben-Jacob, Y. Imry and K. Maki, 1983, Phys. Rev. A 27, 3345. Biller, P. and F. Petruccione, 1990, Phys. Rev. B 41, 2139. Bishop, AR., D.K. Campbell and St. Pnevmatikos, eds, 1989, Springer Proceedings in Physics, Vol. 39: Disorder and Nonlinearity (Springer, Berlin). Bourbonnais, R. and R. Maynard, 1990a, Phys. Rev. Lett. 64, 1397. Bourbonnais, R. and R. Maynard, 1990b, Inter. J. Mod. Phys. 1, 223. Braun, O.M. and Yu.S. Kivshar, 1991, Phys. Rev. B. 43, 1060. Bratus’, EM. and VS. Shumeiko, 1985, J. Low Temp. Phys. 60, 109. Bratus’, EM., S.A. Gredeskul, L.A. Pastur and VS. Shumeiko, 1988, Teor. Mat. Fiz. 76, 401 [Theor.Math. Phys. 76, 945]. Brekhovskikh, L.M., 1973, Waves in Layered Media (Moscow, Nauka) (in Russian). Burlakov, V.M., S.A. Kiselev and V.N. Pyrkov, 1990a, Solid State Commun. 74, 327. Burlakov, V.M., S.A. Kiselev and V.N. Pyrkov, 1990b, Phys. Rev. B 42, 4921. Campbell, D.K., J.F. Schonfeld and C.A. Wingate, 1983, Physica D 9, 1.
60
S.A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
Campbell, D.K., M. Peyrard and P. Sodano, 1986, Physica D 19, 165. Caputo, J.G., AC. Newell and M. Shelley, 1989, in: Lecture Notes in Physics No. 342: Integrable Systems and Applications. eds M. Balabane, P. Lochak and C. Sulem (Springer, Berlin) pp. 49—64. Courant, R. and D. Hilbert, 1962, Methods of Mathematical Physics, Vol. II (New York, Interscience). Currie, J.F., SE. Trullinger, AR. Bishop and iA. Krumhansl, 1977, Phys. Rev. B 15, 5567. Delyon F., YE. Levy and B. Souillard, 1986, Phys. Rev. Lett. 57, 2010. Devillard, P. and B. Souillard, 1986, J. Stat. Phys. 43, 423. Doucot, B. and R. Rammal, 1987a, Europhys. Lett. 3, 969. Doucot, B. and R. Rammal, 1987b, J. Physique 48, 509. Fei, Zhang, Yu.S. Kivshar and L. Vázquez, 1991a, Phys. Rev. A. in press. Fei, Zhang, Yu.S. Kivshar, BA. Malomed and L. Vázquez, 1991b, Phys. Leti. A 159, 318. Flytzanis, C., 1986, in: Nonlinear Phenomena in Solids, ed. Ml. Borissov (World Scientific, Singapore). Flytzanis, N., St. Pnevmatikos and M. Remoissenet, 1985, J. Phys. C 18, 4603. Fogel, MB., SE. Trullinger, AR. Bishop and iA. Krumhansl, 1976, Phys. Rev. Lett. 36, 1411. Fogel, MB., SE. Trullinger, AR. Bishop and l.A. Krumhansl, 1977, Phys. Rev. B 15, 1578. Fraggis, T., St. Pnevmatikos and EN. Economou, 1989, Phys. Lett. A 142, 361. Freylikher, V.D. and S.A. Gredeskul, 1990, J. Opt. Soc. Am. A 7, 868. Freylikher, V.D. and S.A. Gredeskul, 1991, Radio Sci. 26, 375. Freylikher, V.D. and S.A. Gredeskul, 1992, in: Progress in Optics, Vol. 30, ed. E. Wolf, to be published. Gerzenstein, N.E. and V.B. Vasilyev, 1959a, Radiotekh. Electron. 4, 611 (in Russian). Gerzenstein, N.E. and yB. Vasilyev, 1959b, Prob. Theory AppI. 4, 424 (in Russian). Gredeskul, S.A. and L.A. Pastur, 1985, 1. Stat. Phys. 38, 25. Gredeskul, S.A. and L.A. Pastur, 1991, in: Nonlinearity with Disorder, eds F.Kh. Abdullaev, AR. Bishop and St. Pnevmatikos (Springer. Berlin), in press. Gredeskul, S.A. 5 Yu.S. Kivshar and MV. Yanovskaya, 1990, Phys. Rev. A 41, 3994. Gredeskul, S.A., Yu.S. Kivshar, L.M. Maslov, A. Sanchez and L. Vázquez, 1992, Phys. Rev. A, in press. Grudinin, A.E., EM. Dianov, AM. Prokhorov and DV. Khaidarov, 1988, Pis’ma Zh. Tekh. Fiz. 14. 1010. Halding, I. and P.S. Lomdahl, 1988, Phys. Rev. A 37, 2608. Hasegawa, A., 1988, Solilons in Optical Fibers (Springer, Berlin). Hasegawa, A. and F.D. Tappert, 1973a, Appi. Phys. Lelt. 23. 142, Hasegawa, A. and F.D. Tappert, 1973b, AppI. Phys. Lett. 23, 171. Hochstrasser, D., H. Bi.itner, H. Desfontaines and M. Peyrard, 1988, Phys. Rev. A 36, 5332. Ishiwata, S., Y. Okada, S. Watanabe and H. Tanaka, 1990, J. Phys. Soc. Jpn 59, 3029. Karpman, VI. and EM. Maslov, 1977, Zh. Eksp. Teor. Fiz. 73, 537 [Soy. Phys. — JETP 46 (1978) 2811. Karpman, VI. and EM. Maslov, 1978, Zh. Eksp. Teor. Fiz. 75, 504 [Soy. Phys. — JETP 48 (1978) 252]. Kaup, Di., 1976, SIAM I. Appl. Math. 31, 121. Keller, J.B., G.C. Papanicolaou and I. Weilenmann, 1978, Commun. Pure AppI. Math. 32, 583. Kiselev, S.A. and VI. Rupasov, 1990, Phys. Lett. A 148, 355. Kivshar, Yu.S., 1991a, Phys. Rev. A 43, 3117. Kivshar, Yu.S., 1991b, in: Nonlinearity with Disorder, eds F.Kh. Abdullaev, AR. Bishop and St. Pnevmatikos (Springer, Berlin), in press. Kivshar, Yu.S., 1992, IEEE I. Quantum Electron., in press. Kivshar, Yu.S. and BA. Malomed, 1989, Rev. Mod. Phys. 61, 763. Kivshar, Yu.S. and M. Peyrard, 1992, Phys. Rev. A, in press. Kivshar, Yu.S., AM. Kosevich and O.A. Chubykalo, 1987a, Zh. Eksp. Teor. Fiz. 93. 968 [Soy. Phys. —JETP 66 (1987) 545]. Kivshar, Yu.S., A.M. Kosevich and O.A. Chubykalo, 1987b, Phys. Lett. A 125, 35. Kivshar, Yu.S., S.A. Gredeskul, A. Sanchez and L. Vázquez, 1990, Phys. Rev. Lett. 64, 1693. Kivshar, Yu.S., Zhang Fei and L. Vázquez, 1991a, Phys. Rev. Lett. 67, 1177. Kivshar, Yu.S.. A. Sanchez, O.A. Chubykalo, A.M. Kosevich and L. Vázquez, 1991b, Phys. Rev. B, submitted. Kivshar, Yu.S., S.A. Gredeskul, A. Sanchez and L. Vázquez, 1992, Waves in random media, to be published. Klyatskin, V.1., 1980, Stochastic Equations and Waves in Random Media (Moscow, Nauka) (in Russian). Klyatskin, VI., 1988, Imbedding Approach in the Theory of Wave Propagation (Moscow, Nauka) (in Russian). Knapp, R., G.P. Papanicolaou and B. White, 1989, in: Disorder and Nonlinearity, eds AR. Bishop, D.K. Campbell and St. Pnevmatikos, Springer Series in Solid State Science (Springer, Berlin). Knapp, R., G. Papanicolaou and B. White, 1991, 1. Stat, Phys. 63, 567. Kosevich, AM., 1990, Physica D 41, 253. Kotani, S., 1986, Contemp. Math. 50, 277. Krtikel, D., K. Halas, G. Guiliani and D. Grischkowsky, 1988, Phys. Rev. Lett. 60, 29. Landau, L.D. and EM. Lifshitz, 1977, Quantum Mechanics (Pergamon, Oxford).
S. A. Gredeskul, Yu.S. Kivshar, Propagation and scattering of nonlinear waves in disordered systems
61
Landauer, R., 1970, Philos. Mag. 21, 863. Li, 0., St. Pnevmatikos, EN. Economou and CM. Soukoulis, 1988a, Phys. Rev. B 37, 3534. Li, Q., CM. Soukoulis, St. Pnevmatikos and EN. Economou, 1988b, Phys. Rev. B 38, 11888. Lifshitz, tM. and V.Ya. Kirpichenkov, 1979, Zh. Eksp. Teor. Fiz. 77,989. Lifshitz, I.M., S.A. Gredeskul and L.A. Pastur, 1988, Introduction to the Theory of Disordered Systems (Wiley, New York). Malomed, BA., 1985, Physica D 15, 385. McKinstrie, C.J. and G.G. Luther, 1990, Phys. Scr. 30, 31. McLaughlin, D.W. and AC. Scott, 1978, Phys. Rev. A 18, 1652. Papanicolaou, G.C., 1971, 1. Appl. Math. 21, 13. Pascual, P.1. and L. Vazquez, 1985, Phys. Rev. B 32, 8305. Pascual, P.1., L. Vázquez, St. Pnevmatikos and AR. Bishop, 1989, in: Springer Proceedings in Physics, Vol. 39: Disorder and Nonlinearity, eds AR. Bishop, D.K. Campbell and St. Pnevmatikos (Springer, Berlin). Pastur, L.A. and E.P. Fel’dman, 1974, Zh. Eksp. Teor. Fiz. 67, 487. Payton, D.N., H. Rich and W.M. Visscher, 1967, Phys. Rev. 160, 706. Petruccione, F. and P. Biller, 1990, Phys. Rev. B 41, 2145. Peyrard, M. and AR. Bishop, 1989, in: Lecture Notes in Physics, No. 353: Nonlinear Coherent Structures, eds M. Barthes and 1. Leon (Springer, Berlin). Peyrard, M. and D.K. Campbell, 1983, Physica D 9,33. Peyrard, M., St. Pnevmatikos and N. Flytzanis, 1987, Phys. Rev. A 36, 903. Ping Sheng, ed., 1990, Scattering and Localization of Classical Waves in Random Media (world Scientific, Singapore). Pnevmatikos, St., 1988, Phys. Rev. Lett. 60, 1534. Remoissenet, M., 1986, Phys. Rev. B 33, 2386. Sanchez, A. and L. Vázquez, 1991, Inter. 1. Mod. Phys. B 5, 2825. Sievers, Al. and S. Takeno, 1988, Phys. Rev. Lett. 61, 970. Takeno, S. and K. Hon, 1990, 1. Phys. Soc. ipn 59, 3037. Takeno, S., K. Kisoda and A.i. Sievers, 1988, Prog. Theor. Phys. Suppl. No. 94, 242. Toda, M., 1967, 1. Phys. Soc. Jpn 22, 431; 23, 501. Tsuboi, T. and F.M. Toyama, 1991, Phys. Rev. A 44, 2691. Tsurui, A., 1973, 1. Phys. Soc. ipn 34, 1462. \Vadati, M., 1983, 1. Phys. Soc. Jpn 52, 2642. Wadati, M. and Y. Akutsu, 1984, 1. Phys. Soc. ipn 53, 3342. Weiner, AM., i.P. Heritage, RI. Hawkins, RN. Thurston, EM. Kirschner and WI. Tomlinson, 1988, Phys. Rev. Lett. 61, 2445. Zakharov, V.E. and A.B. Shabat, 1972, Zh. Eksp. Teor. Fiz. 61, 118 [Soy. Phys. — JETP 34, 62]. Zakharov, V.E. and A.B. Shabat, 1973, Zh. Eksp. Teor. Fiz. 64, 1627 [Soy. Phys. —JETP 37, 823]. Zakharov, V.E., S.V. Manakov, S.P. Novikov and L.P. Pitaevsky, 1980, Theory of Solitons (Moscow, Nauka) [EnglishTranslation by Consultants Bureau, New York, 1984]. Zolotaryuk, A.V., 1986, Teor. Mat. Fiz. 68, 415. Zolotaryuk, A.V., K.H. Spatschek and E.W. Ladke, 1984, Phys. Lett. A 101, 517.