text:
# 1 Introduction

## 1 Introduction

The great degree of miniaturization of today’s electronic components is such that the de Broglie wavelength of the charge carriers is frequently comparable to the dimensions of the system. Hence, quantum mechanical effects (e.g. tunneling) are expected to play a central role in the behavior of electronic components to be constructed in the next years. In order to describe these quantum effects, it is unlikely that classical transport models will be sufficient. Quantum transport equations, such as the Schrödinger-Poisson or the Wigner-Poisson systems , will therefore be a necessary tool in order to understand the basic properties of these physical systems.

In the present paper, we consider a one-dimensional quantum plasma, where the electrons are described by a statistical mixture of $`N`$ pure states, with each wavefunction $`\psi _i`$ obeying the Schrödinger-Poisson system

$`i\mathrm{}{\displaystyle \frac{\psi _i}{t}}`$ $`=`$ $`{\displaystyle \frac{\mathrm{}^2}{2M}}{\displaystyle \frac{^2\psi _i}{x^2}}e\varphi \psi _i,i=1,\mathrm{},N`$ (1)
$`{\displaystyle \frac{^2\varphi }{x^2}}`$ $`=`$ $`{\displaystyle \frac{e}{\epsilon _0}}({\displaystyle \underset{i=1}{\overset{N}{}}}|\psi _i|^2n_0),`$ (2)
where $`\varphi (x,t)`$ is the electrostatic potential. Electrons have mass $`M`$ and charge $`e`$, and are globally neutralized by a fixed ion background with density $`n_0`$. Finally, we assume periodic boundary conditions, with spatial period equal to $`L`$.

The system of Eqs. (1)–(2) takes into account diffraction, which is the most prominent quantum effect, but neglects dissipation, spin and relativistic corrections. These effects may be important in more realistic models for small semiconductor devices. Nevertheless, it is useful to consider simplified models that capture the main features of quantum plasmas. Indeed, Eqs. (1)–(2) are sufficiently rich to display a wide variety of behaviors, as it will be seen in the rest of this work. At the same time, this model is still amenable to analytic and numerical treatment.

A physically equivalent approach would consist in considering a Wigner function describing the same mixture. The Wigner function approach is a reformulation of quantum mechanics in the classical phase space language . The price to pay for this otherwise appealing formalism, is that Wigner functions can take negative values, and cannot therefore be regarded as true probability distributions. However, both for the analytical and the numerical treatment of the problems of interest in this paper, the Schrödinger-Poisson formalism is more appropriate. This is particularly true for the numerical simulations, since the Wigner formalism is cast into a two-dimensional phase space, whilst the Schrödinger-Poisson model only requires the discretization of the one-dimensional configuration space. Of course, if the number $`N`$ of streams is large, the numerical cost for the description of the system of Eqs. (1)-(2) is also considerable. Nevertheless, interesting physical phenomena (such as instabilities) can take place even with a few streams ($`N=1`$ or 2), as it will appear in the rest of this work. More subtle effects (Landau damping, for instance) would probably require a larger, although hopefully not prohibitive, number of streams ($`N2030`$).

For the analytical study, the hydrodynamic formulation of the Schrödinger-Poisson system is particularly convenient, since it makes direct use of macroscopic plasma quantities, such as density and average velocity. Moreover, it enables one to perform straightforward perturbation calculations in the same fashion as in the classical case. Let us introduce the amplitude $`A_i`$ and the phase $`S_i`$ associated to the pure state $`\psi _i`$ according to

$$\psi _i=A_i\mathrm{exp}(iS_i/\mathrm{}).$$
(3)
Both $`A_i`$ and $`S_i`$ are defined as real quantities. The density $`n_i`$ and the velocity $`u_i`$ of the $`i`$th stream of the plasma are given by

$$n_i=A_i^2,u_i=\frac{1}{M}\frac{S_i}{x}.$$
(4)
Introducing Eqs. (3)–(4) into Eqs. (1)–(2) and separating the real and imaginary parts of the equations, we find

$`{\displaystyle \frac{n_i}{t}}+{\displaystyle \frac{}{x}}(n_iu_i)`$ $`=`$ $`0,`$ (5)
$`{\displaystyle \frac{u_i}{t}}+u_i{\displaystyle \frac{u_i}{x}}`$ $`=`$ $`{\displaystyle \frac{e}{M}}{\displaystyle \frac{\varphi }{x}}+{\displaystyle \frac{\mathrm{}^2}{2M^2}}{\displaystyle \frac{}{x}}\left({\displaystyle \frac{^2(\sqrt{n_i})/x^2}{\sqrt{n_i}}}\right),`$ (6)
$`{\displaystyle \frac{^2\varphi }{x^2}}`$ $`=`$ $`{\displaystyle \frac{e}{\epsilon _0}}({\displaystyle \underset{i=1}{\overset{N}{}}}n_in_0).`$ (7)
Quantum effects are contained in the pressure-like, $`\mathrm{}`$-dependent term in Eq. (6). If we set $`\mathrm{}=0`$, we simply obtain the classical multistream model introduced by Dawson . Therefore, we shall refer to Eqs. (5)-(7) as the quantum multistream model.

Equations (5)–(7) constitute the mathematical model used in the rest of this work. We focus our attention on the one-stream (Sec. 2) and two-stream (Sec. 3) cases, the latter being related to the well-known two-stream instability. The relevant dispersion relations are derived, and the unstable branches are identified. We also investigate the properties of stationary solutions of the Schrödinger-Poisson system, which can be viewed as the quantum counterpart of the classical Bernstein-Greene-Kruskal (BGK) modes . The analytical calculations are checked against time-dependent numerical simulations, shown in Sec. 4. Our conclusions are presented in Sec. 5.

## 2 One-stream plasma

In order to fix the basic ideas, we first consider the one-stream case and take $`N=1`$, i.e. a single pure quantum state. For brevity, we write $`n_1n`$, $`u_1u`$. We obtain

$`{\displaystyle \frac{n}{t}}+{\displaystyle \frac{}{x}}(nu)`$ $`=`$ $`0,`$ (8)
$`{\displaystyle \frac{u}{t}}+u{\displaystyle \frac{u}{x}}`$ $`=`$ $`{\displaystyle \frac{e}{M}}{\displaystyle \frac{\varphi }{x}}+{\displaystyle \frac{\mathrm{}^2}{2M^2}}{\displaystyle \frac{}{x}}\left({\displaystyle \frac{^2(\sqrt{n})/x^2}{\sqrt{n}}}\right),`$ (9)
$`{\displaystyle \frac{^2\varphi }{x^2}}`$ $`=`$ $`{\displaystyle \frac{e}{\epsilon _0}}(nn_0).`$ (10)
A homogeneous zeroth-order solution for Eqs. (8)–(10) is provided by

$$n=n_0,u=u_0,\varphi =0,$$
(11)
where $`u_0`$ is a constant representing the equilibrium velocity of the stream. The linear stability of this solution is obtained by Fourier analyzing Eqs. (8)–(10),

$`n`$ $`=`$ $`n_0+n^{}\mathrm{exp}(i(kx\omega t)),`$ (12)
$`u`$ $`=`$ $`u_0+u^{}\mathrm{exp}(i(kx\omega t)),`$ (13)
$`\varphi `$ $`=`$ $`\varphi ^{}\mathrm{exp}(i(kx\omega t))`$ (14)
Retaining only terms up to first order in $`n^{}`$, $`u^{}`$ and $`\varphi ^{}`$ leads to the dispersion relation

$$(\omega ku_0)^2=\omega _p^2+\mathrm{}^2k^4/4M^2,$$
(15)
where $`\omega _p=(n_0e^2/M\epsilon _0)^{1/2}`$ is the plasma frequency. For $`u_0=0`$, the dispersion relation given in  is recovered, the term $`ku_0`$ merely representing a Doppler shift. As the frequency $`\omega `$ is always real, there can be neither instability nor damping of the wave.

The classical analog of this system is the “cold plasma” model, which is known to sustain nonlinear oscillations when the amplitude of the initial perturbation is smaller than a certain value. Beyond that value, the solution becomes singular in a finite time, which is a sign that the model is no longer valid. This phenomenon corresponds to the breaking of the plasma wave, due to particle overtaking in the phase space. Due to the pressure-like term in Eq. (9), the quantum solution nevers becomes singular, as it was shown by computer simulations .

Let us now turn our attention to the stationary regimes of the system. If all quantities only depend on position, then Eqs. (8)-(9) are reduced to

$`{\displaystyle \frac{d}{dx}}(nu)`$ $`=`$ $`0,`$ (16)
$`u{\displaystyle \frac{du}{dx}}`$ $`=`$ $`{\displaystyle \frac{e}{M}}{\displaystyle \frac{d\varphi }{dx}}+{\displaystyle \frac{\mathrm{}^2}{2M^2}}{\displaystyle \frac{d}{dx}}\left({\displaystyle \frac{d^2(\sqrt{n})/dx^2}{\sqrt{n}}}\right).`$ (17)
Equations (16)–(17) possess the first integrals

$`J`$ $`=`$ $`nu,`$ (18)
$`E`$ $`=`$ $`{\displaystyle \frac{Mu^2}{2}}e\varphi {\displaystyle \frac{\mathrm{}^2}{2M}}\left({\displaystyle \frac{d^2(\sqrt{n})/dx^2}{\sqrt{n}}}\right),`$ (19)
corresponding to charge and energy conservation. The constant $`E`$ can be eliminated by the global shift $`\varphi \varphi +E/e`$, and therefore we assume $`E=0`$ in the rest of this section. Eliminating $`u`$, introducing $`A=\sqrt{n}`$ and using Poisson’s equation, we obtain

$`\mathrm{}^2{\displaystyle \frac{d^2A}{dx^2}}`$ $`=`$ $`M\left({\displaystyle \frac{MJ^2}{A^3}}2eA\varphi \right),`$ (20)
$`{\displaystyle \frac{d^2\varphi }{dx^2}}`$ $`=`$ $`{\displaystyle \frac{e}{\epsilon _0}}(A^2n_0).`$ (21)
It can be easily verified that the $`J=0`$ case cannot sustain small-amplitude, periodic solutions. Hence, we assume $`J=n_0u_0`$ with $`u_00`$ and introduce the following rescaling

$`x^{}={\displaystyle \frac{\omega _px}{u_0}},A^{}={\displaystyle \frac{A}{\sqrt{n_0}}},`$
$`\varphi ^{}={\displaystyle \frac{e\varphi }{Mu_0^2}},H={\displaystyle \frac{\mathrm{}\omega _p}{Mu_0^2}}.`$ (22)
We obtain, in the transformed variables (omitting the stars for simplicity of notation),

$`H^2{\displaystyle \frac{d^2A}{dx^2}}`$ $`=`$ $`2\varphi A+{\displaystyle \frac{1}{A^3}},`$ (23)
$`{\displaystyle \frac{d^2\varphi }{dx^2}}`$ $`=`$ $`A^21,`$ (24)
a system that only depends on the rescaled parameter $`H`$, which is a measure of the importance of quantum effects.

Notice that the classical limit is singular, in the sense that Eq. (23) degenerates into an algebraic equation when $`H=0`$. The equation for the electrostatic potential then becomes

$$\frac{d^2\varphi }{dx^2}=\frac{1}{\sqrt{2\varphi }}1.$$
(25)
Equation (25) has a Hamiltonian character, and corresponds to the equation of motion of a “particle” moving in a potential $`V(\varphi )=\varphi \sqrt{2\varphi }`$. Using this analogy, one can see that Eq. (25) has periodic solutions around the equilibrium $`\varphi =1/2`$, provided that the initial condition satisfies $`0<\varphi (x=0)<2`$. The fact that no solution exists for sufficiently large values of the potential is easily understood. A large potential fluctuation induces a velocity fluctuation, which can drive $`u(x)`$ far from its nominal value $`u_0`$. If the potential is sufficiently strong, $`u(x)`$ can even vanish, but in that case the relation $`nu=J=`$ constant implies an infinite density. This is the well-known effect of particle overtaking that occurs in the cold plasma model.

Going back to the quantum mechanical case, we have shown that Eqs. (23)–(24) describe the inhomogeneous QBGK (quantum-BGK) equilibria of the one-component quantum plasma. We are faced with the mathematical problem of understanding the qualitative properties of the solutions of a coupled, nonlinear system of two second order differential equations depending on a parameter. Only few analytical results can be obtained, shown in the rest of this section.

Equations (23)–(24) can be put into Hamiltonian form by using the variables

$$\overline{A}=iA,\overline{\varphi }=\varphi /H.$$
(26)
Notice that the rescaled amplitude $`\overline{A}`$ is a purely imaginary quantity. We have

$$d^2\overline{A}/dx^2=U/\overline{A},d^2\overline{\varphi }/dx^2=U/\overline{\varphi },$$
(27)
where $`UU(\overline{A},\overline{\varphi })`$ is the pseudo-potential

$$U(\overline{A},\overline{\varphi })=\frac{1}{H}(1+\overline{A}^2)\overline{\varphi }+\frac{1}{2H^2\overline{A}^2}.$$
(28)
Since the equations of motion are autonomous with respect to the independent variable $`x`$, the Hamiltonian formulation immediately gives the first integral (subscripts denote differentiation)

$$I=\frac{1}{2}(\overline{A}_x^2+\overline{\varphi }_x^2)+U(\overline{A},\overline{\varphi }),$$
(29)
which is the Hamiltonian function in transformed coordinates. Transforming back to the original variables, one obtains the first integral for Eqs. (23)–(24)

$$I=\frac{1}{2}(A_x^2+H^2\varphi _x^2)+\frac{1}{H^2}(1A^2)\varphi \frac{1}{2H^2A^2}.$$
(30)
According to Liouville-Arnold theorem , an autonomous two-degree of freedom Hamiltonian system is completely integrable if it possesses two first integrals in involution and with compact level surfaces. Even if $`I`$ has not compact level surfaces, a second constant of motion would be a strong indication of integrability of the QBGK spatial dynamics. We have tried to find a second constant of motion for Eqs. (23)–(24) by a variety of methods (geometrical Noether and Lie symmetries, for instance), but without success. However, numerical integrations of Eqs. (23)–(24) for a wide range of values of $`H`$, and different initial conditions, strongly suggest that, when bounded solutions exist, they are always regular. An additional first integral must therefore exist, although its actual expression may be difficult to guess.

It is interesting to perform a linear stability analysis in order to see in what conditions the system supports small amplitude spatially periodic solutions. Writing

$$A=1+A^{}\mathrm{exp}(ikx),\varphi =1/2+\varphi ^{}\mathrm{exp}(ikx),$$
(31)
and retaining in Eqs. (23)–(24) only terms up to first order in the primed variables, we obtain the relation

$$H^2k^44k^2+4=0.$$
(32)
Again, we point out the singular character of the classical limit: for $`H=0`$, Eq. (32) degenerates into a quadratic equation, with solutions $`k=\pm 1`$. The wavenumbers being always real, this corresponds to spatially periodic solutions. When $`H0`$, we obtain

$$k^2=\frac{2\pm 2\sqrt{1H^2}}{H^2}.$$
(33)
For $`H<1`$ (semiclassical regime), both wavenumbers are real, and therefore the system can sustain spatially periodic oscillations. For $`H>1`$ (strong quantum effects), the solutions are spatially unstable, and grow exponentially. For $`H=1`$, the spectrum is degenerate, with associated secular terms. The corresponding solution is also spatially unstable, growing linearly with $`x`$. We conclude that small-amplitude stationary solutions of the one-stream Schrödinger-Poisson system can only exist in the semiclassical regime, $`H<1`$.

## 3 Two-stream plasma

### 3.1 Two-stream instability

We now turn to the more interesting case of two streams, which classically can give rise to instability. For this, we consider Eqs. (5)–(7) with $`N=2`$. We first linearize around the equilibrium solution $`n_1=n_2=n_0/2`$, $`u_1=u_2=u_0`$, $`\varphi =0`$, where $`u_00`$ is a nonzero reference velocity, and then Fourier transform both in space and time variables. In terms of the dimensionless variables

$$\mathrm{\Omega }=\omega /\omega _p,K=u_0k/\omega _p,H=\mathrm{}\omega _p/Mu_0^2,$$
(34)
the dispersion relation becomes

$$\mathrm{\Omega }^4\left(1+2K^2+\frac{H^2K^4}{2}\right)\mathrm{\Omega }^2K^2\left(1\frac{H^2K^2}{4}\right)\left(1K^2+\frac{H^2K^4}{4}\right)=0.$$
(35)
Consequently, one obtains

$$\mathrm{\Omega }^2=\frac{1}{2}+K^2+\frac{H^2K^4}{4}\pm \frac{1}{2}(1+8K^2+4H^2K^6)^{1/2}.$$
(36)
The solution for $`\mathrm{\Omega }^2`$ has two branches, one of which is always positive and gives stable oscillations. The other solution is negative ($`\mathrm{\Omega }^2<0`$) provided that

$$(H^2K^24)(H^2K^44K^2+4)<0,$$
(37)
Notice that instability, when it occurs, always arises through the marginal mode ($`\mathrm{\Omega }=0`$). In the classical case, Eq. (37) yields $`K^2<1`$, which is the classical criterion for the occurrence of the two-stream instability. In the quantum case, Eq. (37) bifurcates for $`H=1`$. If $`H>1`$, the second factor is always positive, and the plasma is unstable if $`HK<2`$. If $`H<1`$, there is instability if either

$$0<H^2K^2<22\sqrt{1H^2},$$
(38)
or

$$2+2\sqrt{1H^2}<H^2K^2<4.$$
(39)
This yields the stability diagram plotted on Fig. 1. The lower instability zone is semiclassical, as it represents an extension of the classical instability criterion. The upper instability zone, on the other hand, has no classical analog. The two zones merge for $`H=1`$. We call $`K_A,K_B`$ and $`K_C`$ the wavenumbers at which the growth rate vanishes. These are defined by

$`H^2K_{A,B}^2`$ $`=`$ $`2\pm 2\sqrt{1H^2},`$ (40)
$`H^2K_C^2`$ $`=`$ $`4.`$ (41)
It is easy to verify the following property

$$K_A^2+K_B^2=K_C^2.$$
(42)
We shall see that these wavelengths are related to the stationary solutions of the Schrödinger-Poisson system (quantum-BGK, or QBGK, modes). ¿From the previous analysis, it appears that quantum mechanics has a destabilizing effect in the semiclassical regime, where more modes are unstable compared to the classical case. On the other hand, when $`H>2`$, fewer modes turn out to be unstable than in the classical regime. This is the result obtained by Suh, Feix and Bertrand , who found that quantum effects are stabilizing for large enough $`H`$.

It would be interesting to know whether the quantum instability is stronger or weaker than the classical one. For this, one should search for the maximum (over all wavenumbers) growth rate for a fixed value of $`H`$, and compare it to the maximum classical growth rate. Although the algebra is rather involved, direct inspection of the function $`\gamma (K)`$ for different values of $`H`$ indicates that the maximum classical growth rate is always larger than the maximum quantum growth rate (Fig. 2) \[on this figure, the intersections with the $`K`$ axis correspond to wavenumbers $`K_A`$, $`K_B`$ and $`K_C`$ as defined in Eqs. (40)–(41), for which the growth rate vanishes\]. Of course, for some wavenumbers, the quantum growth rate may be larger than the classical one, as can be seen from Fig. 2. Note that the secondary maximum (between wavenumbers $`K_B`$ and $`K_C`$) is considerably smaller than the one between $`K=0`$ and $`K_A`$ , and goes to zero in the classical limit $`H0`$.

In the previous discussion, we have ignored the fact that, due to the periodic boundary conditions, the momentum variable $`p=\mathrm{}k`$ is discrete, and can only be a multiple of $`\mathrm{}k_0`$, where $`k_0=2\pi /L`$ is the fundamental wavenumber. For the previous stability calculations, we have assumed an initial equilibrium solution for which the density $`n_0`$ is spatially uniform and the velocity is equal to $`\pm u_0`$. This corresponds to wavefunctions of the type

$$\psi _\pm (x,0)=\sqrt{n_0/2}\mathrm{exp}(\pm iMu_0x/\mathrm{}).$$
(43)
In order to satisfy the periodicity, one must have

$$Mu_0=n\mathrm{}k_0,n=1,2,3\mathrm{}$$
(44)
Using the dimensionless variables introduced earlier, this condition becomes

$$HK_0=1/n,$$
(45)
which implies $`HK_01`$. This result sets an upper bound on the fundamental wavenumber $`K_0`$. On the instability diagram of Fig. 1, this means that the upper instability region cannot be accessed for $`K_0`$, although of course it can for some of its harmonics. Computing the intersection of the curve $`HK=1`$ with the lower curve on the diagram ($`K=K_A`$), we obtain $`H^2=3/4`$. Therefore, in the region $`H^2>3/4`$, the fundamental wavenumber is always unstable.

The previous discussion raises the question of the physical meaning of the upper instability region in Fig. 1. In particular, we ask whether it is possible to locate the fundamental wavenumber in the stable region (between $`K_A`$ and $`K_B`$), and a higher harmonic $`K_m=mK_0`$ (where $`m`$ is an integer) in the unstable upper zone (between $`K_B`$ and $`K_C`$). The constraints to be satisfied are

$`K_B`$ $`<`$ $`K_m<K_C,`$ (46)
$`K_0`$ $`>`$ $`K_A.`$ (47)
Taking into account Eq. (45), we can write $`K_m=m/(nH)`$, and the previous inequalities become

$`2+2\sqrt{1H^2}`$ $`<`$ $`{\displaystyle \frac{m^2}{n^2}}<4,`$ (48)
$`22\sqrt{1H^2}`$ $`<`$ $`{\displaystyle \frac{1}{n^2}}.`$ (49)
The square root can be eliminated by summing Eqs. (48)-(49). We obtain

$$4n^21<m^2<4n^2,$$
(50)
which cannot be satisfied for any pair of integer numbers $`(n,m)`$. In summary, it is not possible to excite a harmonic in the unstable upper zone of Fig. 1, without also exciting the fundamental mode in the lower unstable region. Therefore, at least in the semiclassical regime $`H<1`$, we cannot expect to observe a “purely quantum” instability. However, as remarked earlier, for $`H^2>3/4`$ the fundamental mode is unstable, as a result of quantum effects. We stress that the above restrictions are a result of the periodic boundary conditions, and do not apply to a truly infinite plasma, for which momentum space is continuous. Still, periodic conditions can be relevant, for instance, to solid state plasmas, where the periodicity is induced by the underlying ion lattice.

### 3.2 Stationary solutions – QBGK modes

Let us now consider the stationary states of the two-stream Schrödinger-Poisson system (QBGK modes). If all quantities are dependent only on position, Eqs. (5)–(6) for $`N=2`$ possess the first integrals

$`J_1`$ $`=`$ $`n_1u_1,J_2=n_2u_2,`$ (51)
$`E_1`$ $`=`$ $`{\displaystyle \frac{Mu_1^2}{2}}e\varphi {\displaystyle \frac{\mathrm{}^2}{2M}}{\displaystyle \frac{d^2(\sqrt{n_1})/dx^2}{\sqrt{n_1}}},`$ (52)
$`E_2`$ $`=`$ $`{\displaystyle \frac{Mu_2^2}{2}}e\varphi {\displaystyle \frac{\mathrm{}^2}{2M}}{\displaystyle \frac{d^2(\sqrt{n_2})/dx^2}{\sqrt{n_2}}}.`$ (53)
We are particularly interested in the case of two symmetric streams, each carrying the same current (with opposite sign) and same kinetic energy. Therefore we write

$$J_1=J_2=n_0u_0/2,E_1=E_2=Mu_0^2/2$$
(54)
where $`u_00`$ is a nonzero reference velocity. Let $`n_1A_1^2`$ and $`n_2A_2^2`$, and transform to the dimensionless variables

$`x^{}=\omega _px/u_0`$ , $`A_{1,2}^{}=A_{1,2}/\sqrt{n_0},`$ (55)
$`\varphi ^{}=(e\varphi +E_1)/Mu_0^2`$ , $`H=\mathrm{}\omega _p/Mu_0^2.`$ (56)
In these variables, the two conservation laws, Eqs. (52)-(53), and Poisson’s equation (7) take the form of a six-dimensional dynamical system (we omit the stars for simplicity of notation)

$`H^2{\displaystyle \frac{d^2A_1}{dx^2}}`$ $`=`$ $`{\displaystyle \frac{1}{4A_1^3}}2\varphi A_1,`$ (57)
$`H^2{\displaystyle \frac{d^2A_2}{dx^2}}`$ $`=`$ $`{\displaystyle \frac{1}{4A_2^3}}2\varphi A_2,`$ (58)
$`{\displaystyle \frac{d^2\varphi }{dx^2}}`$ $`=`$ $`A_1^2+A_2^21.`$ (59)
Equations (57)–(59) constitute a coupled, nonlinear system of three second-order ordinary differential equations, depending on the control parameter $`H`$. They can be put into a Hamiltonian form, using a procedure similar to the one employed for the one-stream QBGK equations (23)–(24). This immediately provides a first integral, which is the Hamiltonian function itself (its actual expression is rather involved, and not particularly illuminating). Just as in the one-stream case, we could not find any additional constants of the motion, so that we cannot prove rigorously that the system is integrable. However, numerical integrations of Eqs. (57)–(59) (see end of this section) strongly suggest that the system is indeed integrable, with quasi-periodic solutions.

For the study of small amplitude oscillations, it suffices to expand Eqs. (57)–(59) in the vicinity of the spatially homogeneous equilibrium $`A_1=1/\sqrt{2}`$, $`A_2=1/\sqrt{2}`$, $`\varphi =1/2`$. After Fourier transforming, the following system is obtained, for the perturbed quantities $`A_i^{},\varphi ^{}`$

$`(4H^2K^2)A_i^{}+\sqrt{2}\varphi ^{}`$ $`=`$ $`0,i=1,2`$ (60)
$`\sqrt{2}(A_1^{}+A_2^{})+K^2\varphi ^{}`$ $`=`$ $`0,`$ (61)
where $`K=ku_0/\omega _p`$ is the dimensionless wavenumber. By searching for nontrivial solutions, one obtains the relation

$$(H^2K^24)(H^2K^44K^2+4)=0,$$
(62)
Notice that this is the same (with an equality sign) as the previously obtained Eq. (37). Solutions of Eq. (62) represent wavenumbers for which both the real and the imaginary part of the frequency vanish, and can be considered as the homogeneous limit of generally inhomogeneous stationary states (QBGK modes). If $`H<1`$ there are three such solutions, which are the wavenumbers $`K_A,K_B`$ and $`K_C`$ defined in Eqs. (40)–(41). If $`H>1`$, only the solution $`K_C`$ survives. The other two solutions become complex, so that spatially periodic QBGK modes can no longer exist.

It is also interesting to look for the eigenvectors corresponding to the eigenvalues $`K_A,K_B`$ and $`K_C`$. From Eqs. (60)–(61) we can express the perturbed amplitudes $`A_i^{}`$ in terms of the potential. Two cases are possible : (1) if $`K=K_C=2/H`$, then one must have $`\varphi ^{}=0`$ and $`A_1^{}=A_2^{}`$ ; (2) otherwise, if $`K=K_A`$ or $`K_B`$, we have $`A_1^{}=A_2^{}=\sqrt{2}\varphi ^{}/(H^2K^24)`$. The first case is particularly interesting. It means that the mode characterized by wavenumber $`K_C`$ is a “quasi-neutral” mode, in the sense that the associated electrostatic potential is zero to first order. Indeed, one can see that no plasma parameters (such as the plasma frequency) enter the definition of $`K_C2/H`$. These modes can be actually accessed, for example by choosing the fundamental wavenumber $`K_0=1/H`$, which is the largest admissible value for $`K_0`$ \[see Eq. (45)\]. In this case, the harmonic $`2K_0=2/H`$ corresponds to the quasi-neutral mode.

Numerical integration of Eqs. (57)–(59) confirms the previous results. For instance, it was verified that periodic solutions only exist for $`H<1`$. We take $`H=0.7`$ and initialize the amplitudes and the potential (at $`x=0`$) with their equilibrium value, plus a small perturbation $`ϵ`$, i.e. $`\varphi (0)=1/2+ϵ_\varphi ,A_i(0)=(1+ϵ_i)/\sqrt{2}`$. In agreement with the discussion of the previous paragraph, if we choose $`ϵ_\varphi =0`$ and $`ϵ_1=ϵ_2`$, the wavenumber $`K_C2.857`$ is linearly excited and thus dominates (Fig. 3), while the potential remains very small. On the other hand, if $`ϵ_1=ϵ_2`$ and $`ϵ_\varphi `$ arbitrary, modes $`K_A1.08`$ and $`K_B2.645`$ are linearly excited (Fig. 4). For generic perturbations, all three wavenumbers are excited. Of course these results are strictly valid only for infinitesimally small perturbations. For moderate values, other modes appear (visible on Figs. 3-4, for which $`ϵ=0.02`$), although the linear wavenumbers are still dominant. For even larger perturbations, bounded solutions no longer exist.

## 4 Time-dependent numerical simulations

A standard numerical technique has been employed in order to integrate the time-dependent Schrödinger-Poisson system. Let us write the Schrödinger equation as

$$i\mathrm{}\frac{\psi }{t}=K\psi +\mathrm{\Phi }\psi ,$$
(63)
where $`K`$ is the kinetic part of the Hamiltonian, and $`\mathrm{\Phi }(x,t)`$ is the potential. The Hamiltonian is split into these two parts, and each is treated separately. For the potential part, the solution is trivial

$$\psi ^{n+1}=\psi ^n\mathrm{exp}(i\mathrm{\Phi }\mathrm{\Delta }t/\mathrm{}),$$
(64)
where $`\psi ^n\psi (n\mathrm{\Delta }t)`$, and $`\mathrm{\Delta }t`$ is the timestep. For the kinetic part, we use a Crank-Nicolson scheme, which is exact to second order in $`\mathrm{\Delta }t`$,

$$i\mathrm{}\frac{\psi ^{n+1}\psi ^n}{\mathrm{\Delta }t}=\frac{1}{2}(K\psi )^{n+1}+\frac{1}{2}(K\psi )^n.$$
(65)
The kinetic operator $`K`$ is spatially discretized by using the standard centered differences formula. The time evolution is obtained by subsequently applying the potential and kinetic steps described above. Poisson’s equation is solved with a Fast Fourier Transform technique just before the potential step (notice that the kinetic step does not alter the spatial density, and therefore the potential). The resulting numerical scheme is unconditionally stable and second order accurate in both space and time variables. Another crucial property of the scheme is that it conserves exactly the integral $`|\psi |^2𝑑x`$ . A typical resolution used in the simulation is $`N_x=512`$ points, and timestep $`\mathrm{\Delta }t=0.02`$.

The initial condition for the simulations is obtained by applying a sudden sinusoidal potential to the equilibrium wavefunctions given in Eq. (43). In dimensionless variables, we have

$$\psi _\pm (x,0)=2^{1/2}\mathrm{exp}(\pm inK_0x)\mathrm{exp}(iϵH^1\mathrm{cos}(K_mx)),$$
(66)
where $`ϵ`$ and $`K_m=mK_0`$ are the amplitude and the wavenumber of the perturbation, and $`n,m`$ are integer numbers. We remind that one must have $`HK_0=1/n`$.

Several simulations have been run in order to compare with the analytical results obtained in the linear regime, with excellent agreement between the two. As an example, we use the parameters $`H=0.25`$ and $`n=5`$ ($`K_0=0.8`$), and perturb the fundamental mode ($`m=1`$). Figure 5 shows the evolution of two modes of the electrostatic potential. The straight line corresponds to the linear growth rate for $`K_0`$, as computed from Eq. (36), and closely matches the measured growth rate. At saturation, several modes are present.

The total momentum distribution $`F(p)`$ is given by the sum of the square modulus of the Fourier transform of each wavefunction, with $`p=\mathrm{}k`$. As pointed out earlier, momentum space is discrete, with $`\mathrm{\Delta }p=\mathrm{}k_0`$ (= 0.2 in the above case). For our simulations, the total momentum distribution of the unperturbed wavefunctions is simply

$$F(p)=\frac{1}{2}\delta (pMu_0)+\frac{1}{2}\delta (p+Mu_0),$$
(67)
where $`\delta `$ is the Dirac delta function. During the linear phase, the momentum distribution remains virtually unchanged, whereas at saturation we observe a significant spreading in momentum space, which extends to $`p\pm 2.5Mu_0`$ (Fig. 6). This is similar to the behavior of the classical two-stream instability.

As was shown in the preceding section, one cannot excite a mode in the unstable upper region of Fig. 1, without also exciting an unstable mode in the lower zone. However, if only the larger wavenumber is initially perturbed, one can hope to see it grow with the correct rate before other unstable modes can be excited. In order to do so, we take $`H=0.9`$ and $`n=5`$. The fundamental wavenumber is thus $`K_0=2/9`$, and is of course unstable since $`K_A=1.18`$ in this case. We only perturb the harmonic $`9K_0=2`$, which falls within the upper unstable zone, since $`K_B=1.883`$ and $`K_A=2.222`$. Of course, several other modes are also unstable (namely, those from $`K_0`$ to $`5K_0`$), but they are not initially perturbed. Figure 7 shows the evolution of mode $`9K_0`$, which closely agrees with the result of the linear calculation for the growth rate. Several other linearly unstable modes appear at a later time, so that at saturation the spectrum is large. This is reflected in the momentum distribution, plotted on Fig. 8.

The results for the fully quantum regime, and particularly for $`H1`$, are more surprising. As an example, let us take $`H=1`$, $`K_0=1`$, and perturb the fundamental mode itself. Instability occurs as expected with the correct growth rate (see Fig. 9). However, instead of saturating at a certain level, the system appears to decay with the same rate, and then to grow again. These nonlinear periodic oscillations do not damp even for very long times, as it has been checked numerically. Higher harmonics, which are linearly stable, are driven by the fundamental mode, and show a similar pattern, although at a lower level. The period of the oscillations is not universal, and depends on the amplitude of the initial perturbation, which confirms that this is indeed a nonlinear effect. A similar, although less pronounced behavior, is also observed for $`H>1`$. The momentum distribution is virtually unchanged over the entire duration of the simulation. It must be noted that, for parameters such that $`HK_0=1`$, the minimum nonzero value of the momentum is $`p_{\mathrm{min}}=Mu_0`$, and is therefore equal to the momentum of the unperturbed streams. In other words, the streams occupy the lowest possible level in momentum space. Although we do not have a detailed explanation for this phenomenon, it appears to be an example of a completely reversible quantum system, in which the initial condition is almost perfectly reconstructed after one period . This is to be compared with the inherently irreversible classical dynamics, for which returning to the initial state after saturation is virtually impossible.

Finally, we have studied the evolution of perturbed QBGK equilibria. It is particularly instructive to consider the case where only the $`K_C=2/H`$ mode is present. This, as detailed in the previous section, is a “quasi-neutral” mode, in the sense that the associated electrostatic potential is zero to first order in the perturbation parameter. We construct a weakly inhomogeneous equilibrium by using the form of Eq. (3) for the wavefunctions $`\psi _{1,2}`$. Both the amplitude and the phase should be sinusoidal functions. For the amplitudes, we have

$$A_{1,2}=\frac{1}{\sqrt{2}}\left(1\pm ϵ\mathrm{cos}(K_Cx)\right),$$
(68)
where $`K_C=mK_0`$ is the wavenumber of the QBGK mode; $`m`$ therefore represents the number of density oscillations. The phases $`S_{1,2}`$ are obtained from the definition of Eq. (4), and remembering that, in dimensionless units, $`u_in_i=\pm 1/2`$. One obtains, to first order in $`ϵ`$

$$S_{1,2}=\pm x2ϵK_C^1\mathrm{sin}(K_Cx).$$
(69)
By virtue of the relations $`nHK_0=1`$ and $`K_C=mK_0=2/H`$, we obtain that $`m=2n`$ must be an even number.

As an example, we take $`H=0.5,K_C=4`$. With $`n=1`$, we obtain $`K_0=2`$, and the number of spatial oscillations is $`m=2n=2`$. Notice that, in this case, the fundamental wavenumber is stable. If we had chosen $`n=2`$, we would have had $`K_0=1`$, which is unstable. With this value of $`H`$, it is therefore possible to construct a stable QBGK mode displaying at most two spatial oscillations. In order to have more oscillations, a smaller value of $`H`$ should be used. For small values of the perturbation parameter $`ϵ`$, the simulations confirm the linear analytical results, and virtually no evolution is observed; besides, the potential fluctuations stay small compared to the density ones. For larger values of the perturbation ($`ϵ=0.1`$), some temporal variations are observed (Fig. 10), since the wavefunctions given by Eqs. (68)–(69) no longer represent an exact stationary state. However, the periodic structure is not destroyed, and the potential fluctuations remain an order of magnitude smaller than the density fluctuations. These results show that a quantum plasma can support almost stationary, quasi-neutral, periodic solutions. These display significant density fluctuations for each stream, but small potential fluctuations. They have no analog in a classical two-stream plasma.

## 5 Conclusion

In this work, we have introduced a quantum multistream model to describe some physical phenomena arising in quantum plasmas. The quantum multistream model may be considered as a discrete version of the Wigner-Poisson model, in the same sense as the classical multistream model is a discrete form of the Vlasov-Poisson system. Indeed, it is well known  that the Wigner-Poisson system is formally equivalent to an infinite set of Schrödinger equations, coupled by a scalar potential obeying Poisson’s equations. However, it is often more appropriate to work with the hydrodynamic formulation of quantum mechanics, since it makes direct use of the same physical quantities that are employed in classical physics (density, velocity, pressure). Moreover, the stability analysis and perturbation calculations become straightforward in the hydrodynamic formulation. On the other hand, we have used the Schrödinger representation for the time-dependent simulations, since accurate numerical techniques for this equation are well-known from the computational literature.

For the case of the two-stream instability, it has been shown that the dispersion relation possesses three branches, one of a semiclassical character and two of a purely quantum nature. It is interesting to observe that even the classical branch reveals some unexpected features: for small $`H`$, quantum effects tend to enhance the two-stream instability. More precisely, some classically stable wavenumbers are destabilized for sufficiently large values of $`H<1`$. On the other hand, strong quantum effect can yield the opposite result: for $`H>2`$, some classically unstable wavenumbers become stable. The purely quantum region of the dispersion relation (the upper unstable zone of Fig. 1) cannot however be excited without also exciting some wavenumber in the lower (semiclassical) region. This means that a purely quantum instability cannot be observed for $`H<1`$. Extensive numerical simulations, run for different values of the relevant parameters, wholly support the analytical results obtained from linear theory. In the fully quantum case, we have observed a surprising, yet unexplained, regime of undamped nonlinear oscillations. This is a purely quantum effect, which is probably linked to quantum recurrences and echoes .

We have also considered the stationary states of the Schrödinger-Poisson system, which can be viewed as the quantum analog of the classical BGK modes. Such QBGK modes are described by a nonlinear system of coupled second-order differential equations, parametrized by the dimensionless Planck’s constant $`H`$. Such system provides an adequate framework for the analysis of QBGK modes, a task which would be rather difficult in the Wigner-Poisson formalism. In particular, numerical simulations have shown that quasi-neutral, spatially periodic, stationary states can be created in the two-stream plasma, and can survive over long times.

Some interesting questions remain to be addressed. In the present work, we have considered in detail only the one and two-stream cases. Further investigations are needed to explain the properties of the quantum multistream model when the number of streams is large. In this case, the Wigner-Poisson system could be a more appropriate model, despite its intrinsic mathematical difficulties. However, an intermediate number of streams might be sufficient to describe the main physical phenomena. Linear and nonlinear Landau damping, for example, should be good candidates to test these ideas, both analytically and numerically.

Acknowledgments

We would like to thank Pierre Bertrand for his valuable comments and suggestions. One of us (F. H.) thanks the Laboratoire de Physique des Milieux Ionisés for the hospitality while this work was carried out, and the Brazilian agency Conselho Nacional de Desenvolvimento Cientifico e Tecnólogico (CNPq) for financial support.

Figure Captions

Figure 1 : Stability diagram for the two-stream plasma. The filled area is unstable. The dashed line corresponds to $`HK=1`$. Lower and middle solid curves: $`K_A^2`$ and $`K_B^2`$ as defined in Eq. (40). Upper solid curve: $`K_C^2`$ as defined in Eq. (41).

Figure 2 : Plot of the squared growth rate $`\gamma `$ as a function of the dimensionless wavenumber $`K`$, for different values of Planck’s constant. $`H=0.5`$, solid line; $`H=1`$, dashed line; $`H=2`$, dotted line. The intersections of these curves with the $`K`$ axis correspond to wavenumbers $`K_A`$, $`K_B`$ and $`K_C`$ as defined in Eqs. (40)–(41). For $`H=0.5`$, only the intersection at $`K_A1.035`$ is shown; the intersections at $`K_B3.864`$ and $`K_C=4`$ are outside the $`K`$ axis range. For $`H=1`$, the intersections are at $`K_A=K_B=\sqrt{2}`$ and $`K_C=2`$. For $`H=2`$, there is only one intersection at $`K_C=1`$.

Figure 3 : Stationary solution of the two-stream Schrödinger-Poisson plasma (QBGK mode), with $`H=0.7`$, $`ϵ_\varphi =0,ϵ_1=ϵ_2=0.02`$. (a) Spatial variation of the density fluctuations $`A_1^{}`$ (solid line), $`A_2^{}`$ (dashed line), and potential fluctuations $`\varphi ^{}`$ (dotted line). Notice that the potential remains small. (b) Fourier transform of $`A_1^{}`$: the linear wavenumber $`K_C2.857`$ is dominant.

Figure 4 : Stationary solution of the two-stream Schrödinger-Poisson plasma (QBGK mode), with $`H=0.7`$, $`ϵ_\varphi =0,ϵ_1=ϵ_2=0.02`$. (a) Spatial variation of the density fluctuations $`A_1^{}`$ (solid line), $`A_2^{}`$ (dashed line), and potential fluctuations $`\varphi ^{}`$ (dotted line). Note that the solid and dashed lines are superposed, since $`A_1^{}A_2^{}`$. (b) Fourier transform of $`A_1^{}`$: the linear wavenumbers $`K_A1.080`$ and $`K_B2.645`$ are dominant.

Figure 5 : Two-stream instability – evolution of the fundamental mode $`K_0=0.8`$ (solid line), and first harmonic $`2K_0`$ (dashed line) of the electrostatic potential, for $`H=0.25`$. The straight line corresponds to the growth rate of $`K_0`$ computed from linear theory, $`\gamma =0.3116`$.

Figure 6 : Momentum distribution for the same case of Fig. 5, at times $`t=0`$ (dotted line) and $`t=80`$ (solid line). The final distribution has been magnified by a factor two. Momentum space is discrete with $`\mathrm{\Delta }p=\mathrm{}k_0=0.2`$.

Figure 7 : Two-stream instability – evolution of mode $`9K_0`$ of the electrostatic potential, for $`H=0.9`$ and $`K_0=2/9`$. The straight line corresponds to the growth rate computed from linear theory, $`\gamma =0.1085`$.

Figure 8 : Momentum distribution for the same case of Fig. 7, at times $`t=0`$ (dotted line) and $`t=160`$ (solid line). The final distribution has been magnified by a factor two. Momentum space is discrete with $`\mathrm{\Delta }p=\mathrm{}k_0=0.2`$.

Figure 9 : Two-stream instability – evolution of the fundamental mode $`K_0=1`$ (solid line), and first harmonic $`2K_0`$ (dashed line) of the electrostatic potential, for $`H=1`$. The straight line corresponds to the growth rate of mode $`K_0`$ computed from linear theory, $`\gamma =0.2297`$.

Figure 10 : Evolution of a strongly nonlinear ($`ϵ=0.1`$) quasi-neutral ($`K=K_C=4`$) QBGK equilibrium, for the set of parameters: $`H=0.5,K_0=2`$. The number of spatial oscillations is $`m=2`$. The density fluctuations $`n_1^{}`$ (solid line) and $`n_2^{}`$ (dashed line), and potential fluctuations $`\varphi ^{}`$ (magnified by a factor 100, dotted line), are shown at different times.


format: text

text:
# Untitled Document

ON NODAL LINES OF NEUMANN EIGENFUNCTIONS

Research partially supported by NSF grant DMS-0071486.

Krzysztof Burdzy

Abstract.  AMS mathematics subject classification (2000): 35P99, 60J65  Key words or phrases: nodal line, reflected Brownian motion, mirror coupling, eigenfunction, Neumann problem. We present a new method for locating the nodal line of the second eigenfunction for the Neumann problem in a planar domain.

1. Introduction. This note is inspired by recent progress on the “hot spots” conjecture of J. Rauch originally proposed in 1974. The conjecture states that the maximum of the second Neumann eigenfunction in a Euclidean domain is attained on the boundary. This formulation is somewhat vague but the counterexamples given by Burdzy and Werner (1999) and Bass and Burdzy (2000) show that even the weakest version of the conjecture fails in some planar domains. The first positive results on the “hot spots” conjecture appeared in Kawohl (1985) and Bañuelos and Burdzy (1999). The last article contains a proof that the conjecture holds in all triangles with an obtuse angle. The theorem was generalized to “lip domains”, i.e., planar domains between two graphs of Lipschitz functions with the Lipschitz constant 1; this more general result was hinted at in Bañuelos and Burdzy (1999) and proved in Atar and Burdzy (2002). See Atar (2001), Jerison and Nadirashvili (2000) and Pascu (2002) for some other positive results.

There is currently a significant gap between the best positive results and the strongest counterexamples to the “hot spots” conjecture. It is widely believed that the conjecture holds in all convex domains in $`𝐑^n`$ for all $`n2`$. But the problem remains open even for triangles whose all angles are acute.

One of the key assumptions in some theorems of Bañuelos and Burdzy (1999), and all results of Jerison and Nadirashvili (2000) and Pascu (2002), is the symmetry of the domain with respect to a straight line. The assumption is only needed to show that the nodal line of the second Neumann eigenfunction lies on the line of symmetry. Then one can replace the original Neumann problem on the whole domain with the mixed Dirichlet-Neumann problem on a nodal subdomain. As it often happens in analysis, the Dirichlet problem (in this case, the mixed Dirichlet-Neumann problem) turns out to be considerably easier to deal with than the original Neumann problem. Hence, one would like to be able to find the location of the nodal line for the second Neumann eigenfunction in domains which are not necessarily symmetric. If this goal is achieved, some further progress on the “hot spots” conjecture may be expected.

Very few methods for finding the nodal line are known; the following short list is probably complete.

(a) If a domain has a line of symmetry, the nodal line lies on that line, subject to some natural extra assumptions.

(b) There are a handful of classes of domains, such as rectangles and ellipses, for which explicit formulae for eigenfunctions are known.

(c) One can find an approximate location of the nodal line in domains which are “long and thin” in the asymptotic sense; see Jerison (2000).

The purpose of this note is to add the following item to this list.

(d) In some domains, one can use the probabilistic method of “mirror couplings” to delineate a region which the nodal line must intersect.

It should be emphasized that one cannot expect to find the location of the nodal line of the second Neumann eigenfunction in the sense of an explicit formula, except in some trivial cases. What one can realistically hope to achieve is to obtain sufficiently accurate information about the nodal line so that this information can in turn be used to prove some other results of interest; the combination of Lemmas 2 and 3 of Burdzy and Werner (1999) is an example of such a result on nodal lines.

Finally, we would like to mention the results of Melas (1992) on the nodal lines in the Dirichlet case.

2. Mirror couplings and nodal lines. Suppose that $`D`$ is a planar domain (open and connected set) with a piecewise smooth boundary. Informally speaking, a mirror coupling is a pair $`(X_t,Y_t)`$ of reflecting Brownian motions in $`\overline{D}`$, such that the line of symmetry $`K_t`$ between $`X_t`$ and $`Y_t`$ does not change on any interval $`(s,u)`$ such that $`X_tD`$ and $`Y_tD`$ for all $`t(s,u)`$. More formally, $`(X_t,Y_t)`$ satisfy the following system of stochastic differential equations,

$$dX=dW+dL,dY=dZ+dM,dZ=dW2mmdW,m=\frac{YX}{|YX|},$$

for times less than $`\zeta =inf\{s:X_s=Y_s\}`$. We set $`Y_t=X_t`$ for $`t\zeta `$. Here $`W`$ is a planar Brownian motion, and $`Z`$ is another Brownian motion for which the increments are mirror images of those of $`W`$, the mirror being the line with respect to which $`X`$ and $`Y`$ are symmetric. The processes $`X`$ and $`Y`$ are reflecting Brownian motions with the singular drifts $`L`$ and $`M`$ on the boundary, i.e., $`_0^{\mathrm{}}\mathrm{𝟏}_D(X_s)|dL_s|=0`$, and similarly for $`M`$. A rigorous construction of the mirror coupling was given in Atar and Burdzy (2002) although mirror couplings had been informally used in the past; see Atar and Burdzy (2002) for a short review of the history of this concept.

Recall that the first eigenvalue for the Laplacian in $`D`$ with Neumann boundary conditions is equal to 0 and the corresponding eigenfunction is constant. The second eigenvalue need not be simple but its multiplicity can be only 1 or 2 (Nadirashvili (1986, 1988)). The nodal set is the set of points in $`D`$ where the second eigenfunction vanishes. If $`D`$ is not simply connected, the nodal set need not be a connected curve. The nodal set of any second Neumann eigenfunction divides the domain $`D`$ into exactly two nodal domains (connected and open sets).

Theorem 1. Suppose that $`A\overline{D}`$ is closed, $`xDA`$, and $`D_1`$ is the connected component of $`DA`$ which contains $`x`$. Assume that for some $`yDD_1`$, the mirror $`K_t`$ for the coupling starting from $`(X_0,Y_0)=(x,y)`$ satisfies $`K_tDA`$ for all $`t[0,\zeta ]`$ a.s. Then for any second Neumann eigenfunction in $`D`$, none if its nodal domains can have a closure which is a subset of $`\overline{D}A`$ containing $`x`$.

Proof. Suppose that $`B`$ is one of the nodal domains for a second eigenfunction and $`x\overline{B}\overline{D}A`$. Let $`C_1`$ be a non-empty open disc centered at a boundary point of $`B`$, so small that $`B_{}=BC_1`$ satisfies $`x\overline{B}_{}\overline{D}A`$, and there exists a non-empty open disc $`C_2DB_{}`$. Assume that for some $`yDD_1`$, if the mirror coupling starts from $`(X_0,Y_0)=(x,y)`$ then $`K_tDA`$ for all $`t[0,\zeta ]`$ a.s. Since the stationary measure for $`X_t`$ is the uniform distribution in $`D`$, $`X_t`$ will hit $`C_2`$ with probability one, and so it will leave $`B_{}`$ with probability 1. Recall that there are exactly two nodal domains—$`B`$ is one of them; let $`B_1`$ denote the other one. Let $`\tau _X(B_{})=inf\{t0:X_tB_{}\}`$ and $`\tau _Y(B_1)=inf\{t0:Y_tB_1\}`$. We will argue that $`\tau _X(B_{})<\tau _Y(B_1)`$ a.s. Suppose otherwise. By the continuity of $`tY_t`$, $`Y_{\tau _Y(B_1)}`$ belongs to $`\overline{B_1D}`$, and so $`Y_{\tau _Y(B_1)}\overline{B}\overline{B}_{}`$. We have assumed that $`\tau _X(B_{})\tau _Y(B_1)`$, so $`X_{\tau _Y(B_1)}\overline{B}_{}`$. Since the points $`X_{\tau _Y(B_1)}`$ and $`Y_{\tau _Y(B_1)}`$ belong to the connected set $`\overline{B}_{}`$, their line of symmetry $`K_{\tau _Y(B_1)}`$ must also intersect $`\overline{B}_{}`$. We have $`K_tDA`$ for all $`t[0,\zeta ]`$ so $`A\overline{B}_{}\mathrm{}`$, which is a contradiction. We conclude that $`\tau _X(B_{})<\tau _Y(B_1)`$ a.s.

Let $`\mu _2>0`$ denote the second Neumann eigenvalue. Then $`\mu _2`$ is the first eigenvalue for the mixed problem in the nodal domain $`B`$, with the Neumann boundary conditions on $`D`$ and the Dirichlet boundary conditions on the nodal line. The fact that $`B_{}`$ is strictly larger than $`B`$ easily implies that $`\mu _2>\mu _2^{}`$, where $`\mu _2^{}`$ is the analogous mixed eigenvalue for $`B_{}`$. Using the well known identification of Brownian motion density with the heat equation solution, we obtain from Proposition 2.1 of Bañuelos and Burdzy (1999) that

$$\underset{t\mathrm{}}{lim}P(\tau _Y(B_1)>t)e^{\mu _2t}(0,\mathrm{}).$$

Since $`\tau _X(B_{})<\tau _Y(B_1)`$ a.s.,

$$\underset{t\mathrm{}}{lim}P(\tau _X(B_{})>t)e^{\mu _2t}<\mathrm{},$$

and so $`\mu _2^{}\mu _2`$, but this contradicts the fact that $`\mu _2>\mu _2^{}`$. Our initial assumption that $`x\overline{B}\overline{D}A`$ for some nodal domain $`B`$ must be false.

We will illustrate our main result with two examples. Before we do that, we recall a few crucial facts about mirror couplings from Burdzy and Kendall (2000) and Bañuelos and Burdzy (1999).

Suppose that $`D`$ is a polygonal domain and $`I`$ is a line segment contained in its boundary. Let $`J`$ denote the straight line containing $`I`$ and recall that $`K_t`$ denotes the mirror line. Let $`H_t`$ denote the “hinge,” i.e., the intersection of $`K_t`$ and $`J`$. Note that $`H_t`$ need not belong to $`D`$. Suppose that for all $`t`$ in $`[t_1,t_2]`$, the reflected Brownian motions $`X_t`$ and $`Y_t`$ do not reflect from any part of $`D`$ except $`I`$. Let $`\alpha _t`$ denote the smaller of the two angles formed by $`K_t`$ and $`J`$. Then all possible movements of $`K_t`$ have to satisfy the following condition.

(M) The hinge does not move within the time interval $`[t_1,t_2]`$, i.e., $`H_t=H_{t_1}`$ for all $`t[t_1,t_2]`$. The angle $`\alpha _t`$ is a non-decreasing function of $`t`$ on $`[t_1,t_2]`$.

A domain $`D`$ with piecewise smooth boundary can be approximated by polynomial domains $`D_n`$. Mirror couplings in the approximating domains $`D_n`$ converge weakly to a mirror coupling in $`D`$. One can deduce which motions of the mirror $`K_t`$ in $`D`$ are possible by analyzing all allowed movements of mirrors in $`D_n`$’s and then passing to the limit. We will not present the details of this limit theorem here. The readers who are concerned about full rigor should add an extra assumption to Example 2 that $`D`$ is a polygonal domain.

Example 1. Suppose that $`D`$ is an obtuse triangle with vertices $`C_1,C_2`$ and $`C_3`$. Let $`C_3`$ be the vertex with an angle greater than $`\pi /2`$. Let $`\overline{C_jC_k}`$ denote the line segment with endpoints $`C_j`$ and $`C_k`$, and let $`d(C_j,C_k)`$ denote the distance between these points. Let the points $`C_4,C_{10}\overline{C_1C_3}`$, $`C_6,C_{11}\overline{C_2C_3}`$, and $`C_5,C_7,C_8,C_9\overline{C_1C_2}`$ be such that $`d(C_1,C_4)=d(C_4,C_3)=d(C_1,C_5)`$, $`d(C_2,C_6)=d(C_6,C_3)=d(C_2,C_7)`$, and the following pairs of line segments are perpendicular: $`\overline{C_4C_8}`$ and $`\overline{C_2C_3}`$, $`\overline{C_6C_9}`$ and $`\overline{C_1C_3}`$, $`\overline{C_5C_{10}}`$ and $`\overline{C_1C_3}`$, $`\overline{C_7C_{11}}`$ and $`\overline{C_2C_3}`$. See Figure 1.

Figure 1.

Let $`A`$ be the closed subset of $`\overline{D}`$ whose boundary consists of line segments $`\overline{C_3C_4}`$, $`\overline{C_3C_6}`$ and $`\overline{C_5C_7}`$ and circular arcs $`\widehat{C_4C_5}`$ and $`\widehat{C_6C_7}`$. The arcs $`\widehat{C_4C_5}`$ and $`\widehat{C_6C_7}`$ are centered at $`C_1`$ and $`C_2`$, respectively. Let $`A_1`$ be the closure of the union of pentagons with vertices $`C_3C_4C_8C_9C_6`$ and $`C_3C_{10}C_5C_7C_{11}`$.

The second Neumann eigenvalue is simple in obtuse triangles by a theorem of Atar and Burdzy (2002).

We will show that

(i) The nodal line for the second Neumann eigenfunction must intersect $`A`$.

(ii) The nodal line lies within $`A_1`$.

Let $`\alpha _t`$ denote the angle formed by the mirror $`K_t`$ with $`\overline{C_1C_2}`$. Let $`\beta _1`$ and $`\beta _2`$ denote the angles such that if $`\alpha _t=\beta _1`$ then $`K_t`$ is perpendicular to $`\overline{C_2C_3}`$ and when $`\alpha _t=\beta _2`$ then $`K_t`$ is perpendicular to $`\overline{C_1C_3}`$. An argument totally analogous to the proof of Theorem 3.1 of Bañuelos and Burdzy (1999) shows that if $`\alpha _0[\beta _1,\beta _2]`$ then $`\alpha _t[\beta _1,\beta _2]`$ for all $`t0`$ a.s.

Consider any point $`xDA`$ to the left of the arc $`\widehat{C_4C_5}`$. Let $`y\overline{D}`$ be such that $`C_1,x`$ and $`y`$ lie on the same line and $`d(C_1,y)=d(C_1,C_3)`$. Consider a mirror coupling $`(X_t,Y_t)`$ with $`X_0=x`$ and $`Y_0=y`$. We will show that $`K_tDA`$ for all $`t0`$.

We have $`K_0DA`$, by the definition of $`y`$. Assume that it is not true that $`K_tDA`$ for all $`t`$. We will show that this assumption leads to a contradiction. If follows from this assumption that for some $`\epsilon >0`$ there exists $`t`$ with $`d(C_1,K_t)d(C_1,C_3)/2\epsilon `$ or $`d(C_2,K_t)d(C_2,C_3)/2\epsilon `$. Let $`T_1`$ be the infimum of such $`t`$. Without loss of generality assume that $`d(C_1,K_{T_1})d(C_1,C_3)/2\epsilon `$. Let $`T_2`$ be the last time $`t`$ before $`T_1`$ when $`d(C_1,K_t)d(C_1,C_3)/2\epsilon /2`$. Note that $`d(C_1,K_t)d(C_1,C_3)/2\epsilon /2`$ for all $`t[T_2,T_1]`$. This implies that $`d(C_1,X_t)d(C_1,C_3)`$ and $`d(C_1,Y_t)d(C_1,C_3)`$ for all $`t[T_2,T_1]`$. Hence, the processes $`X_t`$ and $`Y_t`$ can reflect only on $`\overline{C_1C_2}`$ and $`\overline{C_1C_3}`$ on the time interval $`[T_2,T_1]`$. According to (M), the distance from $`C_1`$ to $`K_t`$ cannot decrease between times $`T_2`$ and $`T_1`$, which contradicts the definition of these times.

Claim (i) now follows from Theorem 1. The second claim is a consequence of the first one and the bounds on the direction of the gradient of the second eigenfunction proved in Theorem 3.1 of Bañuelos and Burdzy (1999).

We will argue heuristically that it is impossible to obtain much sharper estimates for the nodal line location without imposing some extra assumptions on the obtuse triangle. To see this, consider an obtuse triangle with two angles very close to $`\pi /2`$—one slightly less than the right angle and another one slightly larger than that (see Figure 2). Such a triangle has a shape very close to a thin circular sector. Let us assume that the diameter of the triangle is equal to 1. The nodal line in a circular sector is an arc with center at its vertex $`C`$. The nodal line distance from $`C`$ is equal to $`a_0/a_10.63`$, where $`a_0`$ is the first positive zero of the Bessel function of order 0 and $`a_1`$ is the first positive zero of its derivative. This follows from known results on eigenfunctions in discs and estimates for Bessel function zeroes, see Bandle (1980), p. 92. Our methods yield 0.5 as the lower bound for the distance of the nodal line from $`C`$; we see that this estimate cannot be improved beyond 0.63.

Figure 2.

We end this example with a conjecture which, if proved, would considerably strengthen claims (i) and (ii).

Conjecture 1. Suppose that $`D`$ is a triangle with vertices $`C_1,C_2`$ and $`C_3`$ and that the angle at $`C_3`$ is greater than $`\pi /2`$. Let $`C_{12}\overline{C_1C_2}`$ be the unique point such that $`\overline{C_{12}C_3}`$ is perpendicular to $`\overline{C_1C_2}`$, and let $`D_1`$ be the larger of the two connected components of $`D\overline{C_{12}C_3}`$. Then the nodal line for the second Neumann eigenfunction is contained in $`D_1`$.

Little is known about eigenfunctions in triangles different from equilateral; see Pinsky (1980, 1985) for that special case.

Example 2. Consider a convex domain $`D`$ in the plane whose boundary contains line segments $`\{(x_1,x_2):ax_1a,x_2=0\}`$ and $`\{(x_1,x_2):ax_1a,x_2=1\}`$ for some $`a>0`$. Suppose that the remaining parts of $`D`$ are contained in $`\{(x_1,x_2):abx_1a,0x_21\}`$ and $`\{(x_1,x_2):ax_1a+b,0x_21\}`$. We will assume that $`b<a1/(4a)`$. This means that we have to have $`a>1/2`$. Let

$$\begin{array}{cc}\hfill A& =\{(x_1,x_2)\overline{D}:b1/(4a)x_1b+1/(4a)\},\hfill \\ \hfill A_1& =\{(x_1,x_2)\overline{D}:b1/(4a)1x_1b+1/(4a)+1\}.\hfill \end{array}$$

Figure 3 shows an example of $`D`$ and the corresponding rectangle $`A`$.

Figure 3.

We will prove that

(i) The nodal line for any second Neumann eigenfunction must intersect $`A`$.

We will say that $`D`$ is a “lip domain” if its boundary consists of graphs of two Lipschitz functions with the Lipschitz constant 1. The domain depicted in Figure 3 is a lip domain. The second eigenvalue is simple in every lip domain, by a result of Atar and Burdzy (2002).

Our second claim is

(ii) If $`D`$ is a lip domain then the nodal line for the Neumann eigenfunction lies within $`A_1`$.

Suppose that $`xDA`$ lies to the left of $`A`$. If $`x=(x_1,x_2)`$ for some $`x_1a`$ then let $`y=(x_1,x_2)`$. Otherwise let $`y`$ be any point with the first coordinate greater than or equal to $`a`$. Let $`(X_t,Y_t)`$ be a mirror coupling with $`X_0=x`$ and $`Y_0=y`$. Our goal is to show that $`K_tDA`$ for all $`t0`$.

It is elementary, and so it is left to the reader, to check that for any points $`x=(x_1,x_2)`$ and $`y=(y_1,y_2)`$ in $`D`$ with $`x_1a`$ and $`y_1a`$, the line of symmetry for $`x`$ and $`y`$ intersects $`D`$ within $`A`$. Hence $`K_0DA`$, by the choice of $`y`$. Next we follow the same line of thought as in Example 1. Suppose that the condition $`K_tDA`$ fails for some $`t0`$. We will show that this assumption leads to a contradiction. If follows from our assumption that for some $`\epsilon >0`$ there exists $`t`$ and $`(z_1,z_2)K_t\overline{D}`$ with $`z_1b1/(4a)\epsilon `$ or $`z_1b+1/(4a)+\epsilon `$. We can assume without loss of generality that $`z_1b1/(4a)\epsilon `$. Let $`T_1`$ be the infimum of $`t`$ such that $`K_t\overline{D}`$ intersects the line $`\{(z_1,z_2):z_1=b1/(4a)\epsilon \}`$. Let $`T_2`$ be the supremum of $`tT_1`$ such that $`K_t\overline{D}`$ does not intersect $`\{(z_1,z_2):z_1=b1/(4a)\epsilon /2\}`$. Write $`X_t=(X_t^1,X_t^2)`$ and $`Y_t=(Y_t^1,Y_t^2)`$. The following argument applies when $`\epsilon >0`$ is sufficiently small—obviously, we can make such an assumption on $`\epsilon `$. Note that $`X_t^1a`$ and $`Y_t^1a`$ for $`t[T_2,T_1]`$. During this interval, $`Y_t`$ can only reflect on the horizontal parts of the boundary of $`D`$, so according to (M), this reflection cannot decrease $`inf\{z_1:(z_1,z_2)K_t\overline{D}\}`$. The same holds true if $`X_t`$ reflects on the horizontal parts of $`D`$. Suppose that $`X_t`$ reflects on the part of $`D`$ to the left of $`\{(z_1,z_2):z_1=a\}`$ which is not horizontal. Then the hinge $`H_t`$ lies outside $`\overline{D}`$ and it easily follows from (M) that both intersection points of $`K_t`$ with $`D`$ move to the right. Hence, $`inf\{z_1:(z_1,z_2)K_t\overline{D}\}`$ is a non-decreasing function of $`t`$ for $`t[T_2,T_1]`$. This contradicts the definition of $`T_1`$ and $`T_2`$.

The above argument and Theorem 1 yield claim (i). The second claim follows from the first one and the bounds on the direction of the gradient of the second eigenfunction given in Example 3.1 of Bañuelos and Burdzy (1999).

REFERENCES

 R. Atar (2001). Invariant wedges for a two-point reflecting Brownian motion and the “hot spots” problem. Elect. J. of Probab. 6, paper 18, 1–19.

 R. Atar and K. Burdzy (2002), On Neumann eigenfunctions in lip domains (preprint).

 C. Bandle Isoperimetric Inequalities and Applications. Monographs and Studies in Mathematics, 7. Pitman, Boston, Mass.-London, 1980.

 R. Bañuelos and K. Burdzy (1999), On the “hot spots” conjecture of J. Rauch. J. Func. Anal. 164, 1–33.

 R. Bass and K. Burdzy (2000) Fiber Brownian motion and the ‘hot spots’ problem Duke Math. J. 105, 25–58.

 K. Burdzy and W. Kendall (2000) Efficient Markovian couplings: examples and counterexamples Ann. Appl. Probab. 10, 362-409.

 K. Burdzy and W. Werner (1999) A counterexample to the ”hot spots” conjecture. Ann. Math. 149, 309–317.

 D. Jerison (2000) Locating the first nodal line in the Neumann problem. Trans. Amer. Math. Soc. 352, 2301-2317.

 D. Jerison and N. Nadirashvili (2000) The “hot spots” conjecture for domains with two axes of symmetry. J. Amer. Math. Soc. 13, 741–772.

 B. Kawohl, Rearrangements and Convexity of Level Sets in PDE, Lecture Notes in Mathematics 1150, Springer, Berlin, 1985.

 A. Melas (1992) On the nodal line of the second eigenfunction of the Laplacian in $`R^2`$. J. Differential Geom. 35, 255–263.

 N.S. Nadirashvili (1986) On the multiplicity of the eigenvalues of the Neumann problem, Soviet Mathematics, Doklady, 33, 281–282.

 N.S. Nadirashvili (1988) Multiple eigenvalues of the Laplace operator, Mathematics of the USSR, Sbornik, 133-134, 225–238.

 M. Pascu (2002) Scaling coupling of reflecting Brownian motions and the hot spots problem Trans. Amer. Math. Soc. (to appear)

 M. Pinsky (1980) The eigenvalues of an equilateral triangle. SIAM J. Math. Anal. 11, 819–827.

 M. Pinsky (1985) Completeness of the eigenfunctions of the equilateral triangle. SIAM J. Math. Anal. 16, 848–851.

Department of Mathematics
University of Washington
Box 354350
Seattle, WA 98195-4350, USA
Email: burdzy@math.washington.edu
http://www.math.washington.edu/~burdzy/


text:
# On the size of the smallest scales in cosmic string networks

## I Introduction

Topological defects are a generic prediction of grand unified theories. They can be formed in phase transitions when the topology of the vacuum manifold of the low energy theory is non-trivial . For a review see .

The result of a phase transition that produces strings is a network of long strings that stretch across the horizon and a collection of closed loops. If the phase transition occurs at energy $`\eta _s`$, then the mass per unit length in the strings is $`\mu \eta _s^2`$. Immediately after formation, strings in the early universe undergo an epoch of heavy damping until a time $`t_{}`$, which depends on the string formation scale $`\eta _s`$ as well as the specifics of the particle physics model that produces them. As a result of this damping all sub-horizon features, such as loops and perturbations on long strings, are wiped out and the long strings and loops are Brownian with a persistence length $`dt_{}`$.

After that time, numerical simulations show that in an expanding universe there is an attractor solution in which the network evolves into a “scaling regime” (see  and references therein), where the energy density of the string network is a small constant fraction of the radiation or matter density and the statistical properties of the system, such as the correlation lengths of long strings and average sizes of loops, scale with the cosmic time $`t`$. This solution is possible because of intercommutations that produce cosmic string loops which in turn decay by radiating gravitationally. Simulations also have found that most loops and the perturbations on long strings have the smallest possible size, the simulation resolution, which does not scale. The prevailing opinion on this issue is that the size of small-scale structure in fact also scales with the cosmic time $`t`$ and its value is given by the gravitational back-reaction scale $`\mathrm{\Gamma }G\mu t`$, where $`\mathrm{\Gamma }`$ is a number of order 100 and $`G`$ is Newton’s constant. This possibility was first pointed out in .

Cosmic strings are good candidates for a variety of interesting cosmological phenomena such as gamma ray bursts , gravitational waves  and ultra high energy cosmic rays . Some of the predictions of these models, however, depend sensitively on the so far unresolved question of the size of the small-scale structure.

## II Back-reaction model

Hindmarsh  (see also ) showed that the power per unit length radiated into gravitational waves from two colliding arbitrary small perturbations of long repeat length $`L`$ on an infinite string is

$$\frac{dP}{dl}\pi G\mu ^2\underset{n,m}{}(\kappa _n+\kappa _m)ϵ_n^2ϵ_m^2,$$
(1)
where $`\kappa _n=2\pi n/L`$ and $`\kappa _m=2\pi m/L`$ are the wavenumbers of the right- and left-moving Fourier modes that make up the arbitrary perturbations and $`ϵ_n`$ and $`ϵ_m`$ their (small) amplitude to wavelength ratios.

We can use this expression to construct a simple back-reaction model. If we split the sum of the wavenumbers in Eq. (1)

$$\frac{dP}{dl}=\underset{n}{}\frac{dP_n}{dl}+\underset{m}{}\frac{dP_m}{dl},$$
(2)
we can identify individual modes with the power they radiate

$$\frac{dP_n}{dl}\pi G\mu ^2ϵ_n^2\kappa _n\underset{m}{}ϵ_m^2,$$
(3)
with a similar expression for the $`m`$ modes travelling in the opposite direction. Each of the modes contributes to the effective mass per unit length of the string by an amount

$$\delta \mu _n\mu ϵ_n^2$$
(4)
and the power radiated into gravitational waves by each of the modes decreases this contribution

$$\frac{dP_n}{dl}=\frac{d}{dt}(\delta \mu _n).$$
(5)
By putting Eq. (5) together with Eqs. (4) and (3) we arrive at an expression for the time evolution of the amplitude to wavelength ratio of each of the modes

$$\dot{ϵ}_n\pi G\mu \kappa _nϵ_n\underset{m}{}ϵ_m^2$$
(6)
with an analogous expression for the left-moving $`m`$ modes. Each mode therefore loses amplitude at a rate proportional to its own frequency, to its amplitude to wavelength ratio and to a sum which represents the interaction of the mode with every other mode moving in the opposite direction.

There remains, however, some puzzling behavior. In the case of just two colliding modes we can write Eq. (1) as

$$\frac{dP}{dl}\pi G\mu ^2(\kappa _a+\kappa _b)ϵ_a^2ϵ_b^2$$
(7)
where the $`a`$ and $`b`$ subscripts differentiate between the right- and left-moving modes. If we conformally stretch (say) the $`a`$ mode such that $`ϵ_a`$ remains constant while $`\kappa _a0`$, Eq. (7) approaches the constant

$$\frac{dP}{dl}\pi G\mu ^2\kappa _bϵ_a^2ϵ_b^2.$$
(8)
This result is in contradiction to the fact that the power vanishes in the case of perturbations travelling in only one direction on a string : By conformally stretching one of the modes we are making the string on which the other mode is travelling straighter and we expect the power radiated to approach zero in that limit.

This problem was investigated in . Here we will summarise the results. It turns out that in order to be in the Hindmarsh regime (the regime where Eq. (7) is valid) it is not merely sufficient for the amplitudes of each of the colliding modes to be smaller than their corresponding wavelengths: It is also necessary for both amplitudes to be small compared to the geometric mean of the wavelengths, namely it is necessary to have

$$ϵ_a^2\frac{\kappa _b}{\kappa _a}=\frac{A_a^2}{\lambda _a\lambda _b}1$$
(9)
and

$$ϵ_b^2\frac{\kappa _a}{\kappa _b}=\frac{A_b^2}{\lambda _a\lambda _b}1,$$
(10)
where $`A_a`$ and $`A_b`$ are the amplitudes and $`\lambda _a`$ and $`\lambda _b`$ are the wavelengths of the two colliding modes. If instead we are in the regime where one of these two Lorentz invariant quantities, say,

$$ϵ_a^2\frac{\kappa _b}{\kappa _a}=\frac{A_a^2}{\lambda _a\lambda _b}1$$
(11)
the power radiated is exponentially suppressed. This is precisely the regime where the power given by Eq. (7) tends to the constant Eq. (8). It was also found that when the quantity $`ϵ_a^2\kappa _b/8\kappa _a`$ increases past $`1`$ the power drops discontinuously.

We can approximate this behaviour by introducing a cutoff in Eq. (6) such that the sum excludes interactions between modes that satisfy

$$ϵ_a^2\frac{\kappa _b}{\kappa _a}>8,$$
(12)
as follows

$$\dot{ϵ}_n\pi G\mu \kappa _nϵ_n\underset{m}{}ϵ_m^2\theta (8ϵ_n^2\frac{\kappa _m}{\kappa _n})\theta (8ϵ_m^2\frac{\kappa _n}{\kappa _m})$$
(13)
with a similar equation for the $`m`$-modes traveling in the opposite direction. This approximation is justified because of the discontinuous drop in the power that takes place when $`ϵ_a^2\kappa _b/\kappa _a=8`$ and the exponential suppression that sets in when Eq. (11) is satisfied. In fact we cannot easily associate the power radiated with a given mode outside of the Hindmarsh regime and it is fortunate that in this case the power radiated is exponentially supressed and therefore negligible. It should be noted that the calculation in  was performed for just two modes travelling in opposite directions. Here we have assumed (reasonably, we believe) that when a mode interacts with a collection of modes all of which satisfy Eq. (11) the power is also negligible.

Eq. (13) can be used to find the evolution of a set of amplitude to wavelength ratios from some initial conditions. There are two distinct cases of interest. On a long string, one can consider Eq. (13) to give the damping of the average perturbation with a certain wavelength. By symmetry the left-moving and right-moving perturbations will have the same average amplitudes. However, in the case of a loop, a statistical fluctuation may produce an excess of power in a particular direction over a broad range of wavelengths. In that case, since a perturbation passes the same oppositely-directed perturbations over and over again, it is essential to treat the right- and left-movers separately. Here we will consider the case of long strings only.

## III The size of the smallest scales

When the time-scale over which modes on a long string are significantly stretched by the expansion of the universe, which is about $`t`$, the cosmic time, is large compared to the gravitational damping time, Eq. (13) can be used to estimate the size of the smallest scales on long strings: The $`n`$th mode will be significantly damped in amplitude provided $`\tau _g<t`$ with $`\tau _g`$ given by

$$\frac{\dot{ϵ}_n}{ϵ_n}\frac{1}{\tau _g}.$$
(14)
We expect that a combination of stretching by the expansion of the universe and self-intersections on strings will result in an average spectrum of pertubations where the small scales are suppressed relative to the large scales. This seems reasonable even when kinks (produced at intercommutation events) are present since the amplitude to wavelength ratio of Fourier components of kinks as a function of mode number decays like $`1/n`$. The presence of a flat or decreasing spectrum means the sum in Eq. (13) is going to be bounded from below and above: A given mode does not see every other mode on the string but rather a range of modes that lie around it. If we assume a power law spectrum

$$ϵ_nn^\beta ,$$
(15)
one can use Eq. (12) to show that the lower bound on $`m`$ for a given mode $`n`$ is

$$m_{\text{Min}}\left(\frac{n}{8}\right)^{1/(1+2\beta )},$$
(16)
which corresponds to the largest wavelength mode that the $`n`$th mode interacts with while still in the Hindmarsh regime. The smallest wavelength mode that it interacts with can also be found from Eq. (12) and is given by

$$m_{\text{Max}}8n^{1+2\beta }.$$
(17)
We can therefore write Eq. (13) as

$$\frac{\dot{ϵ}_n}{ϵ_n}\pi G\mu \kappa _n\underset{m=m_{\text{Min}}}{\overset{m_{\text{Max}}}{}}ϵ_m^2$$
(18)
with $`m_{\text{Min}}`$ and $`m_{\text{Max}}`$ given by Eqs. (16) and (17) above.

Previously, it was assumed that all modes interact with all other modes and therefore for spectra where the small scales are sufficiently suppressed, namely those with $`\beta >1/2`$,

$$\underset{m}{}ϵ_m^21.$$
(19)
This means that we would expect the survival of modes whose wavelengths

$$\lambda >2\pi ^2G\mu t,$$
(20)
independently of the specific spectrum of perturbations present on the string.

In fact we need to evaluate the sum between the two limits. We can approximate the sum in Eq. (18) by an integral

$$\underset{m=m_{\text{Min}}}{\overset{m_{\text{Max}}}{}}ϵ_m^2_{m_{\text{Min}}}^{m_{\text{Max}}}𝑑mm^{2\beta }$$
(21)
and take the dominant contribution which when $`\beta >1/2`$ is given by

$$_{m_{\text{Min}}}^{m_{\text{Max}}}𝑑mm^{2\beta }\frac{1}{12\beta }m_{\text{Min}}^{}{}_{}{}^{12\beta }.$$
(22)
The dominant contribution therefore comes from the largest wavelength mode that the $`n`$th mode can interact with.

Using Eq. (22) we can write Eq. (18) as

$$\frac{\dot{ϵ}_n}{ϵ_n}\frac{1}{12\beta }\pi G\mu \kappa _nm_{\text{Min}}^{}{}_{}{}^{12\beta }.$$
(23)
Using Eqs. (14) and (16) we can see that

$$\frac{1}{\tau _g}\frac{\mathrm{\Gamma }G\mu }{\lambda }n^{(12\beta )/(1+2\beta )},$$
(24)
where

$$\mathrm{\Gamma }=\frac{8^{(2\beta 1)/(2\beta +1)}}{2\beta 1}2\pi ^2.$$
(25)
Since we want the gravitational lifetime of the mode to be larger than the cosmic time and using $`nt/\lambda `$,

$$\frac{\mathrm{\Gamma }G\mu }{\lambda }\left(\frac{t}{\lambda }\right)^{(12\beta )/(1+2\beta )}<\frac{1}{t},$$
(26)
and therefore only modes with wavelengths

$$\lambda >\left(\mathrm{\Gamma }G\mu \right)^{(1+2\beta )/2}t,$$
(27)
have a significant amplitude at time $`t`$.

## IV Simple Examples

When the amplitude is small compared to the wavelength, one can find the effect of cosmological expansion on the spectrum of perturbations on a string . If the wavelength of a mode on a long string is larger than the horizon, then it is stretched conformally; its amplitude and wavelength grow with the scale factor. If, on the other hand, the wavelength of the mode is smaller then the horizon, only its wavelength is stretched. The net result of this cosmological processing is that the amplitude to wavelength ratio of modes larger than the horizon as a function of the wavenumber remains constant and for modes inside the horizon it becomes a power law.

It is not hard to see why this is so. If a mode enters the horizon at a time $`t_0`$, then its wavelength and amplitude are both $`t_0`$. While inside the horizon the amplitude remains fixed and the wavelength is stretched so that at sometime $`t`$ later

$$\lambda _0(t)t_0\left(\frac{t}{t_0}\right)^\alpha =t_0^{1\alpha }t^\alpha ,$$
(28)
where $`\alpha =1/2`$ or $`2/3`$ depending on whether we are in the radiation or matter era. We can therefore write the amplitude to wavelength ratio of that mode as

$$ϵ_0(t)\frac{t_0}{t_0^{1\alpha }t^\alpha }=\frac{t_0^\alpha }{t^\alpha }.$$
(29)
If we imagine that our mode of wavelength $`\lambda _0(t)`$ is some fraction of the largest mode, $`\lambda _0(t)=\lambda /n`$ with $`\lambda t`$, it is easy to see that

$$t_0=tn^{1/(1\alpha )}$$
(30)
which when substituted into Eq. (29) yields the spectrum

$$ϵ_nn^{\alpha /(1\alpha )}=\{\begin{array}{cc}n^1& \text{Radiation Era}\\ n^2& \text{Matter Era.}\end{array}$$
(31)
Feeding the spectrum Eq. (31) into Eq. (27) yields

$$\lambda >\left(\mathrm{\Gamma }G\mu \right)^{(1+\alpha )/2(1\alpha )}t=\{\begin{array}{cc}\left(\mathrm{\Gamma }G\mu \right)^{3/2}t& \text{Radiation Era}\\ \left(\mathrm{\Gamma }G\mu \right)^{5/2}t& \text{Matter Era}\end{array}$$
(32)
as the small scale structure cutoffs in the radiation and matter eras. However, this analysis neglects intercommutations, which produce kinks, and the possibility that small perturbations propagating on curved horizon-sized strings may not be efficiently stretched.

## V Conclusions

We have presented a method for the calculation of gravitational back-reaction on cosmic strings taking into account that not all Fourier modes that make up the perturbations on strings interact with all other modes with the same efficiency.

In particular, modes of a given wavelength only interact significantly with a narrow range of other modes whose wavelengths are comparable to it. This range depends on the spectrum of perturbations on the string at late times (when gravitational effects become important) which in turn has the effect of making the small scale structure cutoff sensitive to this spectrum.

Using this method and assuming a power law spectrum we have arrived at an expression for the small scale structure cutoff. We have further applied it to the spectrum resulting from the stretching of small amplitude waves in a radiation- or matter-dominated FRW universe.

It is unclear to what extent the results obtained for these two simple cases apply to a realistic network of cosmic strings because we have ignored the effect of intercommutations, which leads to kinks, as well as the possibility that small perturbations propagating on curved horizon-sized strings may not be significantly stretched by the expansion of the universe. Both of these effects have the potential to change the spectrum of perturbations and therefore the value of the small scale structure cutoff.

However, it is clear that because of the restricted range of interaction affecting every mode, the small scale structure cutoff will be smaller than the one given by the usual calculation where an unrestricted range of interaction is assumed.

###### Acknowledgements.

We would like to thank Jose Juan Blanco-Pillado for useful discussions. The work of KDO and AV was partially funded by the NSF.


text:
# Rolling as a “continuing collision”

11institutetext:  Moscow State University, Physics Department, Moscow 119899, Russia
Humboldt-Universität, Institut für Physik, Invalidenstr. 110, D-10115 Berlin, Germany
http://summa.physik.hu-berlin.de/$``$kies/ (Received: / Revised version: )

## Abstract

We show that two basic mechanical processes, the collision of particles and rolling motion of a sphere on a plane, are intimately related. According to our recent findings, the restitution coefficient for colliding spherical particles $`ϵ`$, which characterizes the energy loss upon collision, is directly related to the rolling friction coefficient $`\mu _{\mathrm{roll}}`$ for a viscous sphere on a hard plane. We quantify both coefficients in terms of material constants which allow to determine either of them provided the other is known. This relation between the coefficients may give rise to a novel experimental technique to determine alternatively the coefficient of restitution or the coefficient of rolling friction.

###### Keywords:

Granular materials, restitution coefficient, rolling friction

Is the collision of two spheres related to rolling motion of a sphere on a plane? Probably the answer will be “No”, if the idealised version of these two basic processes (i.e. absolutely elastic collision and ideal rolling without any resistance) is assumed. However, in reality neither ideal collision nor ideal rolling occur since in both processes mechanical energy is lost according to dissipative material deformation. Quantifying the losses a profound similarity of these two processes has been revealed. It results in a relation between quantities describing energy loss due to collision and rolling and provides a novel view on rolling which may be considered as a continuing collision.

It is well known since Newton’s times that two balls colliding with velocities $`v_1`$ and $`v_2`$ have smaller after-collisional velocities $`v_1^{}`$ and $`v_2^{}`$ as compared with the perfectly elastic collision. The decrease of the particle velocity or in other words the loss of mechanical energy is described by the coefficient of restitution $`ϵ`$. For a pair of identical particles colliding on a line one has:

$`v_1^{}=v_1{\displaystyle \frac{1+ϵ}{2}}g`$ (1)
$`v_2^{}=v_2+{\displaystyle \frac{1+ϵ}{2}}g`$
where $`g=v_1v_2`$ is the initial relative velocity. The ideal elastic collision is described by $`ϵ=1`$, so that $`v_1^{}=v_2`$ and $`v_2^{}=v_1`$. The relative velocity $`g=v_1v_2`$ changes at the collision as

$$g^{}=ϵg$$
(2)
We want to remark that for the more general case of a collision Eqs. (1,2) hold true but the pre- and after-collisional velocities are vectors $`\stackrel{}{v}_{1/2}`$, $`\stackrel{}{v}_{1/2}^{}`$ while $`\stackrel{}{g}`$ and $`\stackrel{}{g}^{}`$ are the pre- and after-collisional relative normal velocities

$`\stackrel{}{g}`$ $`=`$ $`\left[\left(\stackrel{}{v}_1\stackrel{}{v}_2\right)\stackrel{}{e}\right]\stackrel{}{e}`$ (3)
$`\stackrel{}{g}^{}`$ $`=`$ $`\left[\left(\stackrel{}{v}_1^{}\stackrel{}{v}_2^{}\right)\stackrel{}{e}\right]\stackrel{}{e}`$
where $`\stackrel{}{e}=(\stackrel{}{r}_1\stackrel{}{r}_2)/\left|\stackrel{}{r}_1\stackrel{}{r}_2\right|`$ with $`\stackrel{}{r}_{1/2}`$ being the particle positions at the instant of collision.

The non-ideal character of the collision originates from the dissipative force

$$F_{\mathrm{diss}}=\frac{3}{2}A\rho \sqrt{\xi }\dot{\xi }$$
(4)
acting on the colliding spheres which has been derived recently BSHP ; KuwabaraKono ; MorgadoOppenheim . Here, $`\xi `$ is the time-dependent compression of the particles during the collision

$$\xi =2R\left|\stackrel{}{r}_1\stackrel{}{r}_2\right|$$
(5)
with $`R`$ being the radius of the spheres, see Fig.1. The material constant $`\rho `$ depends on the Young modulus $`Y`$ and Poisson ratio $`\nu `$ of the particle material as

$$\rho =\frac{2Y}{3(1\nu ^2)}\sqrt{R/2}$$
(6)
The constant $`A`$ is expressed in terms of the elastic and viscous material constants BSHP . The dissipative force always acts against the relative velocity $`\dot{\xi }`$ of the particles, so that elastic energy stored during the collision is not completely reconverted into kinetic energy after the collision. In linear approximation with respect to the dissipative parameter $`A`$ the solution of the collision problem which accounts for the dissipative force yields for the restitution coefficient TomThor

$$1ϵ=C_1\left(\frac{3A}{2}\right)\left(\frac{2\rho }{m}\right)^{2/5}g^{1/5}\pm \mathrm{}$$
(7)
with TomThor ; Rozaetal

$$C_1=\frac{\mathrm{\Gamma }(3/5)\sqrt{\pi }}{2^{1/5}5^{2/5}\mathrm{\Gamma }(21/10)}=1.15344.$$
(8)
$`\mathrm{\Gamma }(x)`$ is the Gamma-function and $`m`$ is the mass of the (identical) particles. Hence, the restitution coefficient for viscoelastic collisions depends sensitively on the normal component of the relative velocity $`g`$ as described by Eq. (7).

Each collision of particles certainly terminates after some time if attractive forces are excluded, but one can ask, whether some “collision” can proceed permanently.

In spite of being far from apparent one can show that the mechanics of a rolling sphere on a hard plane is intrinsically similar to that of a collision BTEurLet  provided the main part of the dissipation originates from viscoelastic bulk deformations of the sphere 6 ; 7 ; 8 . Indeed, when a soft sphere rolls we notice the same sequence of compression and subsequent decompression as in the collision process. Thus, in rolling motion, which one can treat as a “continuing” collision, we observe a steady dissipation due to incomplete retransformation of elastic energy during decompression. This results in a rolling friction moment $`M`$, which acts against the motion:

$$M=\mu _{\mathrm{roll}}F^N$$
(9)
(see Fig.1). Here, $`F^N`$ is the normal force exerted by the plane onto the sphere caused by the sphere’s own weight. Calculations performed for a soft sphere rolling on a hard (undeformable) plane reveal that

$$\mu _{\mathrm{roll}}=AV,$$
(10)
where $`V`$ is the sphere’s linear velocity and $`A`$ is exactly the same material constant as in the law of collision BTEurLet . The relation (10) has been obtained for the case when the deformation of the rolling sphere is small compared to its radius, the velocity $`V`$ is much less then the speed of sound in the material, and when the relaxation time of the rolling process, estimated as the ratio of the sphere deformation and the velocity $`V`$ is much larger then the dissipative relaxation times of the viscoelastic material BTEurLet .

Thus,

$$\frac{1ϵ}{b\left(\rho /m\right)^{2/5}g^{1/5}}=\frac{\mu _{\mathrm{roll}}}{V}$$
(11)
with

$$b=3C_1/2^{3/5}=2.28296$$
(12)
relates the rolling friction coefficient $`\mu _{\mathrm{roll}}`$ to the coefficient of normal restitution $`ϵ`$ for identical particles colliding with relative velocity $`g`$.

The linear Eq.(11) refers to the case of of small impact velocities of the colliding spheres. This situation is the most favourable for performing the experimental studies. One can, however, apply a more general nonlinear relation which accounts for high-order corrections of the dissipative parameter $`A`$. This follows from the high-order expansion for the restitution coefficient Rozaetal :

$`1ϵ`$ $`=`$ $`C_1\left({\displaystyle \frac{3A}{2}}\right)\left({\displaystyle \frac{2\rho }{m}}\right)^{2/5}g^{1/5}`$ (13)
$``$ $`C_2\left({\displaystyle \frac{3A}{2}}\right)^2\left({\displaystyle \frac{2\rho }{m}}\right)^{4/5}g^{2/5}`$
$`+`$ $`C_3\left({\displaystyle \frac{3A}{2}}\right)^3\left({\displaystyle \frac{2\rho }{m}}\right)^{6/5}g^{3/5}`$
$``$ $`C_4\left({\displaystyle \frac{3A}{2}}\right)^4\left({\displaystyle \frac{2\rho }{m}}\right)^{8/5}g^{4/5}\pm \mathrm{}`$
where Rozaetal

$`C_2`$ $`=`$ $`3/5C_1^2`$
$`C_3`$ $`=`$ $`0.315119C_1^3`$
$`C_4`$ $`=`$ $`0.161167C_1^4`$
Introducing then

$$Zb\left(\frac{\rho }{m}\right)^{2/5}g^{1/5}\frac{\mu _{\mathrm{roll}}}{V}$$
(14)
one can write the generalization of the linear Eq.(11) as

$$\frac{1ϵ}{Z}=10.6Z+0.315119Z^20.161167Z^3\pm \mathrm{}$$
(15)
Physically the relations (11) and (15) are based on the intrinsic mechanical similarity of the collision and rolling processes. In practice Eq.(11) \[or (15)\] may be used to determine either of coefficients $`ϵ`$ or $`\mu _{\mathrm{roll}}`$ provided the other one is known.

The measurement of restitution coefficients for spheres is a complicated experimental problem, in particular if one is interested in collisions at very low impact velocity, e.g. Bridges . These values are of great importance in several problems, e.g. for the description of the kinetics of planetary ring material where the particles typically collide with velocities of the range $`\left(10^2\mathrm{}10^3\right)`$ m/sec. The derived relation between the coefficient of rolling friction and the coefficient of normal restitution allows to determine the latter value by measuring the resistance of a sphere against rolling on a hard plane which might be experimentally less complicated.

###### Acknowledgements.

The authors want to thank Michel Louge and Thomas Schwager for helpful discussion