Symmetry-breaking instability in a prototypical driven granular gas

11 Pages·0.27 MB·

by  Evgeniy Khain

Symmetry-breaking instability in a prototypical driven granular gas
Evgeniy Khain and Baruch Meerson
Racah Institute of Physics, Hebrew University of Jerusalem, Jerusalem 91904, Israel
Symmetry-breakinginstabilityofalaterallyuniformgranularcluster(stripstate)inaprototypical
drivengranular gas isinvestigated. The system consists of smooth hard disksin atwo-dimensional
box, colliding inelastically with each other and driven, at zero gravity, by a “thermal” wall. The
2 limitofnearlyelasticparticlecollisionsisconsidered,andgranularhydrodynamicswiththeJenkins-
0 Richman constitutive relations is employed. The hydrodynamic problem is completely described
0 by two scaled parameters and the aspect ratio of the box. Marginal stability analysis predicts a
2
spontaneous symmetry breaking instability of the strip state, similar to that predicted recently for

n a different set of constitutive relations. If the system is big enough, the marginal stability curve
a becomesindependentofthedetailsoftheboundaryconditionatthedrivingwall. Inthisregime,the
J densityperturbationisexponentiallylocalized attheelastic walloppositetothethermalwall. The

1 short-andlong-wavelengthasymptoticsofthemarginalstabilitycurvesareobtainedanalyticallyin
3 thedilute limit. The physicsof thesymmetry-breakinginstability is discussed.

] PACSnumbers: 45.70.Qj
t
f
o
s I. INTRODUCTION ticle simulations [9] showed that, for typical initial con-
.
t ditions, the system evolvesto a state where the particles
a
m are separated into two groups. Almost all particles form
Granular materials play an important role in indus-
aclusterinasmallregionofspace,wheretheymovewith
- trialapplications,geophysicsandastrophysics. Theyare
d verysmallvelocities,whileaveryfewremainingparticles
alsoofagreatgeneralinteresttophysicists,aseachofthe
n movewithhighvelocities. Clearly,this steadystate can-
“phases”ofgranularmaterials: solid,liquidandgashave
o
not be described by granular hydrodynamics (actually,
c unusualpropertiesthatdistinguishthemfromtheirclas-
by any coarse-grained theory). Therefore, the results of
[ sic (atomic) counterparts [1, 2]. We will consider gran-
Kadanoff et al. [9] brought into question the validity of
1 ular gas (or rapid granular flow) and focus on a variant
granular hydrodynamics in general.
v of clustering instability, a striking tendency of granular
9 gases to form dense clusters. Clustering instability was This question was addressed in two subsequent the-
6 first discussedin the contextof a freely ”cooling”granu- oretical works [10, 11] that dealt with similar systems
5
largas[4,5,6]. Followingtheseworks,relatedclustering in two dimensions. Esipov and Po¨schel [11] investi-
1
phenomena wereinvestigatedindrivengranulargasesas gated an ensemble of nearly elastic hard disks in a cir-
0
2 well, both in experiments [7, 8] and in particle simula- cular box with the circumference serving as thermal
0 tions [9, 10, 11]. wall. Grossman et al. [10] considered a rectangular
/ box, one side of which served as thermal wall. Parti-
t Granular clustering results from energy losses by in-
a cle simulations [10, 11] showed granular clusters: dense
m elastic collisions, and it is a manifestation of thermal
and“cold”regionsofgranulatedevelopingawayfromthe
condensationinstability,alsoencounteredinotherfields,
- thermalwall. Intermsofthecoarse-grainedparticleden-
d for example in gases and plasmas that cool by their own
n radiation [12]. Since the discovery of the clustering in- sity, these steady-state clusters had simple shapes: az-
o imuthally uniform (circular state) [11] and laterally uni-
stability,thevalidityofgranularhydrodynamics(seeRef.
c form (strip state) [10]. Grossman et al. also showed
[13]forareview)hasbeenunderscrutiny[3]. Incontrast
:
v that, for nearly elastic collisions, the strip state is de-
totheclusteringinafreely”cooling”granulargas,where
i scribable by a steady-state solution of granular hydro-
X onedealswithacomplextime-dependentprocess,steady
dynamic equations. The empiric constitutive relations
r states are achievable in driven granular systems. One of
a the simplest settings ofthis type is driving the granulate suggestedby Grossmanet al. used simple interpolations
between the low-density limit, where the constitutive re-
byasidewallatazerogravity. Therefore,anensembleof
lationsarederivablesystematically[18],andhigh-density
inelasticallycollidinghardspheres,confinedinaboxand
limit where, close to the dense packing, free volume ar-
drivenbyoneortwo“thermal”wallshasservedasapro-
guments can be used.
totypicaldrivengranularsystem[7,9,10,14,15,16,17].
Steady states of this system have served as test beds The results of Refs. [10, 11] showed that the anomaly
for granular hydrodynamics and its violations. The first observedinthe one-dimensionalsetting [9] does not per-
analysisofthis systeminthe physicalliteraturewasper- sist in higher dimensions. Clustered states qualitatively
formed, in one dimension, by Kadanoff et al. [9]. The similar to those of Grossman et al. were observed in ex-
nearly elastic particles were constrained to move on a perimentofKudrolliet al. [7] whoinvestigatedasystem
straightlinewithenergyinputfromtheboundaries. Par- of spherical particles in a box, rolling on a smooth sur-
2
face and driven by a rapidly vibrating side wall. The dicted in Ref. [17] is sensitive to the constitutive rela-
number of particles served as the control parameter in tions. Livne et al. [17] employed the empiric relations
Ref. [7]. A dense cluster of the strip type was observed suggested by Grossman et al. [10]. Here we shall use
when the number of particles was big enough, in much the better known Jenkins-Richman (JR) relations [20].
the same way as in Ref. [10]. The basic physics of the WhiletherelationsofGrossmanetal. aremoreaccurate
strip state is simple and can be explained by the fol- forhighdensities(evenincludingthoseclosetothedense
lowinghydrodynamicargument. Becauseofthe inelastic packinglimit),theJRrelationsshouldworkbetteratlow
collisions the particle random motion slows down (that and intermediate densities. We shall see, however, that
is,the granulartemperaturedecreases)withthe increase the marginal stability curves, obtained with these two
of the distance from the driving wall. To maintain the sets of relations, are not much different from each other.
momentumbalance,thegranulardensityshouldincrease This implies that the symmetry-breaking instability is
withthisdistance. Whenthetotalnumberofparticlesis robust. Our second objective is to get more insight into
big enough (the rest of parameters being the same), the the marginal stability problem and, where possible, to
density contrast becomes large, and the enhanced den- obtainanalytic results. We shallshow that the marginal
sity regionaway from the driving wall is observedas the problem is equivalent to an eigenvalue problem of quan-
strip state. tum mechanics. An important finding here is a univer-
The prototypical system exhibits many interesting sal behavior of the marginal stability curves in the limit
phenomena of non-hydrodynamic nature. These include when the density perturbations are strongly localized at
inelasticcollapse[11,15],possiblelackofscaleseparation the elastic wall opposite to the thermal wall. In the di-
[10],non-Gaussianityintheparticlevelocitydistribution lute limit, we obtain analytically the short- and long-
[10, 16], normalstress difference and pressure nonunifor- wavelength asymptotics of the marginal stability curves
mities [14], etc. For nearly elastic collisions, however, anddensityeigenfunctions. We alsogivethe physicalin-
granular hydrodynamics was shown to yield an accurate terpretationto thesymmetry-breakinginstabilityandto
quantitative description in the dilute limit [14], and rea- the density borders of the instability region.
sonablyaccuratedescriptionformoderateandhighgran- The rest of the paper is organized as follows. In Sec.
ular densities [10]. Of course, the nearly elastic limit is 2 we formulate the model and briefly describe the strip
quiterestrictiveformostofrealisticgranularflows. Still, state: the simplest steady state of the prototypical sys-
this limit is conceptually important just because gran- tem. Section 3 presents marginal stability analysis of
ular hydrodynamics can be used there. Granular hy- thestripstateandcomparestheresultsobtainedfortwo
drodynamics has a great predictive power and helps to different sets of constitutive relations. More results on
identifyimportantcollectivephenomena(shearflowsand marginalstability, including some analytic results in the
vortices, shocks, different modes of clustering flows etc.) dilute limit, are presented in Sec. 4. Section 5 includes
thataredifficult,ifnotimpossible,toidentifyandpredict a discussion and summary.
in the language of individual particles. Once identified,
thesephenomenacanthenbe investigatedinexperiment
and simulations in more general (not necessarily hydro- II. PROTOTYPICAL SYSTEM AND STRIP
dynamic) formulations. STATE
Therefore, granular hydrodynamics provide a leading-
orderapproachtoabigensembleofnearly-elasticallycol- .
liding grains. This approach has been taken recently by Theprototypicaldrivengranularsystemintwodimen-
Livne et al. [17] who employed granular hydrodynamics sions include inelastically colliding hard disks of diame-
for astability analysisofthe stripstate describedabove. ter d and mass m=1, moving in a box with dimensions
The analysis revealed a spontaneous symmetry-breaking L L . The gravity force is zero. Collisions of disks
x y
×
instabilityofthestripstatewithrespecttoperturbations with the walls x = 0,y = 0 and y = L are assumed
y
along the strip. Well within the instability region, the elastic. The wall x = L is ”thermal” wall: upon colli-
x
clustered states with broken symmetry, found by a nu- sion a particle is assigned a random velocity taken from
merical solution of the steady state hydrodynamic equa- aGaussianensemblewithtemperatureT0. Energyinput
tion, are strongly localized in the lateral direction: most atthethermalwallbalancestheenergydissipationdueto
oftheparticlesarelocatedindense“islands”[17]. These inter-particlecollisions,sothe systemcanreacha steady
results indicate that the prototypical system can show a state. Weshallparameterizetheinelasticityoftheparti-
non-trivial behavior even in the leading-order, hydrody- cle collisions by a constant normal restitution coefficient
namiclimit. Indeed,thissystemscanbeputintothelist r and work in the nearly elastic limit: 1 r2 1. In
− ≪
of pattern-forming systems far from equilibrium [19]. this limit, the Navier-Stokes granular hydrodynamics is
The present work focuses on a more detailed stability expected to be sufficiently accurate in a system with a
analysis of this system. Our first objective is to check bignumberofparticlesandsmallKnudsennumber. The
to what extent the symmetry-breaking instability pre- possible steady states ofthe systemare describedby the
3
steady state versions of the momentum and energy bal- where f = n /n is the area fraction of the grains and
c
h i
ance equations: n = N/(L L ) is the average number density of the
x y
h i
grains.
p=const and (κ T)=I, (1)
∇· ∇ Equations (2)-(4) and the boundary conditions make
a complete set. One can see that the governing param-
where p is the granular pressure, T is the granular tem-
eters of this system are the scaled parameters η, f and
perature, κ is the thermal conductivity and I is the
∆. If the system is infinite in the y-direction, only two
rate of energy losses by collisions. We assume that the
governing parameters: η and f remain. Notice that the
number density n is not too big: n/n 0.5, where
c
n =2/(√3d2) is the (hexagonal)dense pa≤cking density. steady-state density distributions are independent of the
c
This assumption enables us to employ the constitutive wall temperature T0, in contrast to the similar problem
with gravity, where the gravity acceleration, combined
relations derived by Jenkins and Richman [20]. For the
steady state problem, the required constitutive relations with T0 and the (finite) system size in the direction of
gravity, would make an additional governing parameter.
include the equation of state p = p(n,T) and relations
for κ and I in terms of n and T. The laterally uniform steady state (strip state) corre-
Let us introduce scaledcoordinates: r/L r. In the sponds to the one-dimensional (y-independent) solution
x
newcoordinatestheboxdimensionsare1 ∆,→where∆= z =Z(x). It is described by the equations
×
L /L is the aspect ratio. Introducing the normalized
y x
inverse density z(x,y)= nc/n(x,y), one can rewrite the (FZ′)′ =ηQ, Z′|x=0 =0,
1
energy balance equation in Eq. (1) in terms of z(x,y): 1
and Z (x)dx=f, (5)
−
Z0
(F(z) z)=ηQ(z), (2)
∇· ∇
where primes stand for the x derivatives. For the strip
where F(z)=A(z)B(z),
state, the boundary condition at the wall x = 1 drops
G 1+ 9π 1+ 2 2 out. This implies, in particular, that the density profile
A(z) = hz1/2(116+(cid:0)2G)35G/2(cid:1) i, othfethderisvtirnigp wstaaltle(tisheirnmdeapleonrdevnibtroaftinthgeweaxlal)ct[2n1a]t.urTehoisf
π z(z+ 16π√3) degeneracy of the strip state is caused by the character
B(z) = 1+2G+ √3 (z π )3 , of particle interaction: the hard-core potential does not
− 2√3 introduce any characteristicenergy [11]. Notice that, in-
6 z1/2G stead of prescribing the grain area fraction f, one can
Q(z) = ,
π(1+2G)3/2 prescribe the inverse density Z = Z0 at x = 0. This
G(z) = π z− 327√π3 , (3) dcoenfidnietsioan,Ccaoumchbyinpedrowblietmh tfhoer nZo(-xfl)u.xScoolvnidnigtiothneaCtaxu=chy0
2√3(z π )2
− 2√3 problem, one can then compute, from the last equation
inEq. (5),therespectivevalueoff. Atfixedη,thereisa
and η = (2π/3)(1 r)(L /d)2. Notice that, for an arbi-
− x one-to-onecorrespondencebetweenZ0 andf. Therefore,
trarysmallbutfiniteinelasticity1 r,the dimensionless
− analternativeparameterizationofthestripstateisgiven
parameter η can be made arbitrary large, if the system
bythescalednumbersη andZ0. Weshallseebelowthat
sizeL islargeenough. Theparameterη differsfromthe
x the same property keeps (and can be conveniently used)
parameter used by Livne et al. [17] only by a numer-
L in the marginalstability problem. Figure 1 shows a typ-
ical factor of order unity. Of most interest are regimes
ical example of the scaled density profile n(x)/n of the
where η 1, see below. c
≫ strip state obtained by solving Eqs. (5) numerically.
The boundary conditions for Eq. (2) are determined
bythepropertiesoftheparticle-wallinteractions. Atthe
elastic walls x=0, y =0 and y =∆ we shouldprescribe
III. INSTABILITY OF THE STRIP STATE:
azeronormalcomponentoftheheatflux. Intermsofthe
MARGINAL STABILITY CURVES
inverse density z we have z = 0 at these three walls.
n
∇
Here index n denotes the gradientcomponent normal to
In general, the strip state is only one of the possible
the wall. The constant temperature at the “thermal”
solutions of Eq. (2). Because of its nonlinearity, Eq. (2)
wall x = 1 yields the condition ∂z(x= 1,y)/∂y = 0. To
mayhaveadditionalsolutionssatisfyingthesamebound-
maketheformulationoftheproblemcomplete,onemore
ary conditions. When exist, these additional solutions
condition is needed. In experiment or particle simula-
are truly two-dimensional: the translational symmetry
tions, the number of particles N is fixed. This yields a
along y is broken. An important class of these solutions
normalization condition:
bifurcate supercritically fromthe stripstate [17]. There-
1 ∆
1 dxdy fore, close to the bifurcation point, these solutions can
=f, (4)
∆Z0 Z0 z(x,y) be found by linearizing Eq. (2) around the strip state.
4
0.4
marginal stability curves: the curves k = k(f) at differ-
ent values of η, computed numerically. In these compu-
0.3 tations, the parameter η was taken large enough. The
strip state is unstable below the respective curve and
stable above the curve. Notice that, at fixed η, the in-
c
n 0.2 stability is possible only within a finite interval of f:
/
n f1(η) < f < f2(η). The same property was reported
in Ref. [17] for another set of constitutive relations. We
0.1 shall give a physical explanation to this finding in Sec.
4. Notice (see also Ref. [17]) that, at large η, the high-
density stability border f2 is quite small. The curves in
0 Fig. 2 are actually plotted in scaled coordinates: kη−1/2
0 0.2 0.4 x 0.6 0.8 1 versus fη1/2. It can be seen that, in the scaled coor-
dinates, all the curves exit from the same point of the
FIG. 1: An example of the scaled density profile n(x)/n of horizontal axis fη1/2. In addition, the maxima of all
c
thestripstatefortheJRconstitutiverelations. Thegoverning the curves are equal. These observations will be also ex-
parameters are η=104 and f =0.0342. plained in Sec. 4.
0.06
A similar analysis was performed in Ref. [17] for the
constitutive relations of Grossman et al. [10]. In the
framework of a time-dependent hydrodynamic formula-
tion,thisanalysiscorrespondstomarginal stability anal- 0.04
ysisofthestripstatewithrespecttosmallperturbations 2
/
along the strip [17]. 1
−
Substituting z(x,y) = Z(x) + ψ (x)cosky and lin- η
k k
earizing Eq. (2) with respect to the small correction
ψ (x)cosky, we obtain: 0.02
k
ηQ
Z 2
φ +k φ=0. (6)
′′
−(cid:18) F (cid:19)
Hereφ(x)=F(x)ψk(x),functionsF andQareevaluated 0
at z = Z(x), and subscript Z means the z derivative 1 2 fη1/2 3 4
evaluated at z =Z(x). The boundary conditions are
φ(x=0)=0 and φ(x=1)=0. (7) FIG. 2: Marginal stability curves for different values of η,
′
plotted in scaled coordinates: kη−1/2 versus fη1/2. For a
Equation(6)coincideswiththe Schr¨odingerequationfor fixedηthestripstateisstableabovetherespectivecurveand
an even wave function φ(x) of a particle in the potential unstablebelowthecurve. Thevaluesofηare: 104(solidline),
2.5·104 (dashedline), 5·104 (dottedline)and 105 (dash-dot
well
line).
ηQZ if x <1,
U(x)=(cid:26)+F oth|er|wise. (8) Ifthesystemisinfiniteinthelateraldirection,∆= ,
∞ while η and f are fixed, a whole continuum spectrum∞of
2
The quantity k serves as the energy eigenvalue. The wavenumbersbetweenk =0andk =k(η,f)isunstable.
−
energylevelsinthepotential(8)arealwaysdiscrete,and Both in experiment, and in numerical simulations ∆ is
there is an infinite number of them. However, as the finite. In this case k becomes discrete because of the
wave number k should be real, only negative or zero en- boundary conditions: k = mπ/∆, where m = 1,2,....
ergy levels are admissible. At fixed values of η and f, For each m we can find the critical value of the aspect
the potential (8) admits at most one such energy level. ratio ∆ (let us call it ∆ ) such that for ∆ > ∆ the
m m
Theabsenceofnegativeenergylevelsimpliesthat,inthe strip state looses stability with respect to the m-mode.
vicinity of the strip state, there are no steady states dif- Obviously, ∆m = m∆1, and ∆1 is the lowest critical
ferent from it. The presence of a negative energy level value for the symmetry-breaking instability. Figure 3
correspondsto a “weaklytwo-dimensional”steady state, shows ∆1 as a function of f for different values of η.
bifurcating from the strip state. We shall exploit the For fixed η and f, the strip state is unstable above the
quantum-mechanical analogy more fully in Sec. 4. Here curve. It is seen from Fig. 3 that, in order to observe
we report some numerical results. Figure 2 shows the the symmetry-breaking instability, one does not need to
5
workwithverylargeaspectratios∆: itissufficientifthe 4
system is big enough, so that parameter η is sufficiently
large.
3
∆
3
2
2 1
∆
0
0.01 0.02 0.03 0.04
1 f
FIG. 4: The critical aspect ratio ∆1(f) for the symmetry-
0 breakinginstabilityascomputedfortheconstitutiverelations
0.01 f 0.03 0.04 of Grossman et al. [10] (dash-dotted line) and of JR [20]
(dotted line). Parameter η=11,094.
FIG. 3: The critical aspect ratio ∆1 for the symmetry-
breaking instability as a function of f for different values of coordinate: x¯=xη1/2. Interms ofthe original,physical
η. For a fixed η, the strip state is stable below the respec-
coordinate x the new rescaling is independent of the
tivecurveandunstableabovethecurve. Theparametersare: ph
η1 =104 (solid line), η2 =2.5·104 (dashed line), η3 =5·104 system size:
(dotted line) and η4 =105 (dash-dot line).
2π 1/2 (1 r)1/2x
ph
x¯= − .
To what extent is the symmetry-breaking instability (cid:18) 3 (cid:19) d
sensitivetothepreciseformoftheconstitutiverelations?
Equations (5) for the strip state become
We compared the marginal stability curves ∆ = ∆1(f)
for different values of η with the respective curves [17]
(FZ′)′ =Q, Z(x¯=0)=Z0
found for the constitutive relations of Grossman et al.
[10]. A typical example of this comparison is shown in and Z′(x¯=0)=0, (9)
Fig. 4. One can see that, qualitatively, the results are
while the eigenvalue problem (6) and (7) reads
thesame: thebothcurvesdescribeasymmetry-breaking
instability at a critical value of the aspect ratio that Q
φ Z +k¯2 φ=0, φ(x¯=0)=0
depends on the area fraction. In both cases, there are ′′ ′
−(cid:18) F (cid:19)
sharplow-andhigh-densitybordersofinstabilityregion.
and φ(x¯=η1/2)=0. (10)
Therefore, we can conclude that the instability is robust
and does not require a very special form of the consti-
Now the primes denote the derivatives with respect to
tutive relations. On the other hand, there is a notice- x¯, while k¯ = kη 1/2 is the new scaled wave number.
−
able (about 15%) difference in the exact positions of the Like the coordinate x¯, the new scaled wave number k¯ is
marginal stability curves, so the instability provides a
independent of L :
x
good quantitative test for constitutive relations of gran-
ular hydrodynamics. 1/2
3 k d
k¯= ph ,
(cid:18)2π(cid:19) (1 r)1/2
−
IV. MORE RESULTS ON MARGINAL where k is the physical wave number. The problem
ph
STABILITY
(9) and (10) is determined by two parameters: Z0 and
η. However, η enters the rescaled equations only in one
In this Section we investigate the marginal stability place: in the last boundary condition in Eq. (10) where
problem in more detail and obtain some analytic results it determines the scaled system size. If the wave func-
in the dilute limit. tion φ(x¯) is strongly localized in the potential well U(x¯)
(correspondingly,the negative energy level is sufficiently
deep),theresultsfork¯andφ(x¯)becomeindependentofη
A. Localization and universality atsufficiently largeη. Indeed, in this case one can safely
move the boundary x¯ = η1/2 to infinity. It is important
Letuscharacterizethestripstatebythescaledparam- that, in this case, the exact form of the boundary con-
eters η and Z0 and introduce a different rescaling of the dition at the driving wall becomes insignificant, leading
6
0.06
onlytoexponentiallysmallcorrections[21]. Thisuniver- 0.03
sal“localizationregime”wasdiscoveredinRef. [17]that
2
employedadifferentsetofconstitutiverelations. Thecri- −1/ 0.02
terion for localization can be obtained from the require- η
ment that the localizationlength(which is of orderk¯−1) 0.04 k
be much smaller than the (scaled) system size η1/2. In /2 0.01
1
thephysicalunitsitcorrespondstotheshort-wavelength −
η
tlihmisitcoriftetrhieonb,ifk¯urschaotuinldg bsoelufatironen:okupghhLxfro≫m1t.heTboofrudlefirlsl k 02 2.3 2.6
of the instability interval, where k¯ vanishes. In the next 0.02 Z
0
subsection we will work in the dilute limit and rewrite
thiscriterionintermsofthe governingparametersofthe
problem.
Figure5showsthemarginalstabilitycurvesk¯=k¯(Z0) 0
fordifferentvaluesofη, obtainednumerically. Insteadof 0 50 Z 100 150
the borders f1(η) and f2(η) of the instability interval in 0
terms of parameter f, the respective borders in terms of
parameter Z0 appear. One can see that, for large values FIG. 5: Marginal stability curves for different values of η,
of η, the marginal stability curves coincide in a wide re- plottedincoordinatesk¯=kη−1/2versusZ0. Forafixedηthe
strip state is stable above the respective curve and unstable
gion of Z0 not too close to the borders of the instability
below the curve. The inset shows the splitting of the curves
interval. This region corresponds to strong localization.
near the high-density stability border. The values of η are:
Figures 6-8 show the form of the potential (8) and the 104 (solid line), 2.5·104 (dashed line), 5·104 (dotted line)
negative energy level k2 in three characteristic cases and 105 (dash-dot line).
−
(in these figures we returned to the rescaling of the co-
ordinates and wavenumber by the system size L ). Fig-
x
0
ure 6 corresponds to the region of parameters where the
energy level is deep and eigenfunction is localized. Fig-
ures 7 and 8 correspond to the parameter regions close −50
to the low- and high-density borders of the instability,
respectively. There is no localization here. Notice the −100
U
qualitative change in the form of the potential near the
high-density stability border. Figure 9shows the respec- −150
tive eigenfunctions in these three cases.
−200
−250
−1 −0.5 0 0.5 1
B. Marginal stability borders and physics of the x
instability
FIG. 6: An example of the negative energy level E =−k2 ≃
−83.01 (dashed line) in the potential U(x) (solid line) in the
The low- and high-density stability borders f1(η) and regime of localization. The parameters are η = 2.5·104 and
f2(η) (or respective borders in terms of Z0) are deter-
minedbythezero-eigenvalue(k¯ =0)solutionofEq.(10). Z0 =6.
This solution can be found if we know the strip solution
Z(x¯;Z0)ofEq.(9). Indeed, itiseasyto checkthatfunc-
scaled pressure
tion φ0(x¯) = F ∂Z/∂Z0 is a solution of Eq. (10) with
k¯ =0, satisfying the boundary condition φ(x¯=0). Em-
′
ploying the second boundary condition φ(x¯ = η1/2) = 0 p 1+2G T
P = = .
we obtain ∂Z1/∂Z0 = 0, where Z1 = Z(x¯ = η1/2;Z0). ncT0 Z T0
For a given η, this equation is an algebraic equation for
Z0. Our numerical results imply that this equation has As P is independent of the coordinates,we can compute
only two solutions corresponding to the low- and high- itatthethermalwallx¯=η1/2. HereT =T0andZ =Z1,
density instability borders. The instability borders have so we arrive at
aclearphysicalmeaningwhichshedslightonthephysics
of the instability. Let us consider the granular pressure 1+2G(Z1)
P = =P(η,f).
p = nT(1+2G) [20] of the strip state, and introduce a Z1
7
0
1.5
−1
1
−2
φ
U
−3
0.5
−4
0
−5 0 0.2 0.4 0.6 0.8 1
−1 −0.5 0 0.5 1 x
x
FIG. 7: An example of the negative energy level E =−k2 ≃ FIG.9: Eigenfunctionsφ(x)correspondingtotheeigenvalues
shown in Figs. 6 (strong localization, solid line), 7 (dash-dot
−0.155 (dashed line) in the potential U(x) (solid line) in the
absence of localization. The parameters η = 2.5·104 and line)and8(dashedline). Theeigenfunctionsare(arbitrarily)
normalized so that φ(0)=1.
Z0 = 75 correspond to the region close to the low-density
stability border.
the interval. The physical interpretation of these results
1200 200 is clear. The presence of the anomalous (falling) part
of the P(f) curve indicates instability, and it is caused
by the destabilizing role of collisional heat losses [5, 12].
800 U0
We can say that, on the interval (f1,f2), the granulate
U has negative lateral compressibility. At f < f1 the heat
−200
400 losses are too small to cause instability. The presence of
−0.1 0 0.1 the high-densityborderf2 is causedbythe finite-density
x
corrections to the constitutive relations (that is, by the
0
finite size of the particles). This is in contrast to radia-
tive condensations in gases and plasmas, where such a
−400 stabilizing effect would be absent. Now consider a small
−1 −0.5 x0 0.5 1 density modulation of the strip state with a very long
wavelength: k 0. For this perturbation, the stabiliz-
→
ingeffectofthelateralheatconductionvanishes,andthe
FIG. 8: Another example of the negative energy level E =
−k2 ≃−4.1(dashedline)inthepotentialU(x)(solidline)in negativecompressibilitymakesthestripstateunstableon
theabsenceof localization. Theparametersη=2.5·104 and the interval (f1,f2). For a non-zero k, the lateral heat
Z0 = 2.3 correspond to the region close to the high-density conduction has a stabilizing effect. Therefore, a density
stability border. modulation with too a short lateral wavelength should
be stable, as Fig. 2 indeed shows.
Now let us compute the derivative ∂P/∂f at a constant
η:
C. Dilute limit
∂P ∂P ∂Z1 ∂Z0
= .
∂f ∂Z1 ∂Z0 ∂f Inthissubsectionweshallworkinthe dilute limitand
deriveseveralanalyticresults. Weshallseethat,atlarge
One can easily check that the first and third multipliers
η, the dilute limit faithfully reproduces the low-density
intherighthandsideofthisrelationarealwaysnegative.
parts of the marginal stability curves.
Therefore, the sign of ∂P/∂f is determined by the sign
of∂Z1/∂Z0. Aswehaveseen,themarginalstabilitybor-
ders are determined by equation ∂Z1/∂Z0 = 0. There-
fore, the steady-state pressure has its extremum points 1. Strip state and marginal stability problem
points exactly at the points f1 and f2. Figure 10 shows
an example of the dependence P =P(f) at a constantη In the dilute limit, Z 1, Eq. (9) for the strip state
≫
for the strip state, found numerically. One can see, that becomes
P is a decreasing function of f within the instability in-
terval (f1,f2), and an increasing function of f outside (Z3/2)′′ =3Z−1/2, (11)
8
6.2
eigenvalue problem includes a single parameter ξ, while
κ serves as the eigenvalue.
6 We havebeenunableto solveEq. (15)analyticallyfor
a general κ. Still, several important asymptotics can be
obtained.
5.8
P
5.6
2. Zero-energy state and stability border f1
5.4 For k¯ =0 Eq. (15) can be solved analytically:
0 0.02 0.04 φ(ζ,k¯ =0) φ0(ζ)= ζ ζ 1 arccosh ζ, (16)
f ≡ − −
p p p
In other words, we impose a zero eigenvalue k¯ = 0 and
FIG. 10: The scaled steady-state granular pressure P versus find the low-density stability border f1 =f1(η) from the
the grain area fraction f for the strip state. The two circles
boundary condition φ0(ζ1) = 0. We obtain an algebraic
cefforercetsipveonladtetroatlhceommparregsinsiablilsittyaboilfittyhebogradseirssnfe1gaantidvef2o.nTthhee equation coth(ξ1) = ξ1 for ξ1 = f1η1/2. Its solution is
interval (f1,f2). Parameter η=104. ξ1 = 1.19968.... This result explains why all marginal
stabilitycurvesshowninFig.2depart(almost)fromthe
same point at the low-density side. Figure 11 compares
wheretheprimesdenotethex¯derivatives. Theboundary the scaled quantity f1η1/2 at different η, found numeri-
conditions are Z′(x¯ = 0)=0 and Z(x¯ = 0)= Z0, where cally from Eq. (6), with this analytic prediction (a con-
Z0 is related to f and η by the normalization condition stant). Theagreementisverygoodforlargeη. Asηgoes
√ηZ 1dx¯ = fη1/2. The solution of this problem is down, f1 increases and the dilute approximation starts
0 −
eRlementary: to deteriorate.
1.4
x¯= Z0 arccosh ζ+ ζ2 ζ , (12)
2 (cid:16) p p − (cid:17)
where ζ =Z/Z0 and 2
/
1
η
2η1/2 1
Z0 = fη1/2+ 1 sinh(2fη1/2). (13) f
2
1.2
(Returningforamomenttotheoldrescalingofthecoor-
dinate, x /L x, one can see that the density profile
ph x
(12) is determin→ed by a single parameter: ξ = fη1/2.)
Eq. (10)for φ(x¯)=(√3/2)Z1/2(x¯)ψ (x¯)takesthe form
k
1
φ κ φ=0, (14)
′′−(cid:18) − ζ2(cid:19)
1
where κ = k¯2Z02, and ζ = ζ(x¯) is given, in an implicit 100 200 η1/2 300 400
form,byEq. (12). TheboundaryconditionsforEq. (14)
remain the same as in Eq. (10). As ζ(x¯) is a monotonic
function ofx¯, we canchange the independent variablein FIG. 11: The scaled low-density stability boundary ξ1 =
Eq. (14) from x¯ to ζ. The resulting equation for φ(ζ) is f1η1/2atdifferentη1/2asfoundnumericallyfromEq. (6)(cir-
cles). Solid line shows the analytical result ξ1 = 1.19968...
2 obtained in thedilute limit.
4(ζ 1)ζφ′′+2φ′+(1 κζ )φ=0, (15)
− −
where the primes now denote the ζ derivatives. The Figure 12 compares the analytic result for the zero-
function φ(ζ) is defined on the interval 1 ζ ζ1, energyeigenfunction φ0(x), givenby Eqs. (16)and(12),
2 ≤ ≤ 4
where ζ1 = Z1/Z0 = cosh ξ. One boundary condition withanumericalsolutionofEq. (6)forη =10 . Theco-
is φ(ζ = ζ1) = 0 to which we may add an arbitrary ordinatexinFig.12is rescaledbyLx. The analyticand
normalization condition φ(ζ = 1) = 1. An additional numeric results are obtained for slightly different values
boundarycondition, φ′(ζ =1)=(κ 1)/2,canbe found of f1 (see Fig. 11). One can see that the agreement is
−
from Eq. (15) itself, after substituting there ζ =1. This excellent.
9
1
Returning to the parametersη andf,andto the wave
number k = k L , we can rewrite the asymptotics k¯ =
0.8 ph x
A/Z0 as
0.6 A 1
0 k = fη1/2+ sinh(2fη1/2) . (17)
φ 2 (cid:20) 2 (cid:21)
0.4
This asymptotics is shown in Fig. 14.
20
0.2
0
0 0.2 0.4 x 0.6 0.8 1
15
FIG. 12: The zero-energy eigenfunction φ0(x) computed nu-
merically from Eq. (6) (circles) and given analytically by
Eqs.(16) and (12) (solid line). The parameters are η = 104, k 10
f1,num =0.0124 (circles) and f1,anal=0.0120 (solid line).
3. Short-wavelength limit: localization and universality
5
In the short-wavelength limit the system boundary
2
ζ = ζ1 = cosh ξ can be moved to infinity. This re-
quires a strong inequality ξ = fη1/2 1. In this limit, 0
theeigenvalueproblem(15)doesnoti≫ncludeanyparam- 1 2 fη1/2 3 4
eter. The eigenvalue κ should therefore be a number of
order of unity, hence k¯ = A/Z0, with constant A of or- FIG. 14: Marginal stability curves for different values of
der of unity. The constant can be found numerically: η, plotted in coordinates k versus ξ = fη1/2. For a fixed
A 0.525. This simple result represents the low-density η, the strip state is stable above the respective curve and
lim≃it of the “universal” marginal stability curve, corre- unstable below the curve. The values of η are: 104 (solid
sponding to strong localization. Figure 13 shows this line), 2.5 ·104 (dashed line), 5· 104 (dotted line) and 105
5 (dash-dot line). Also shown is the dilute-limit asymptotics
asymptoticsforη =10 . Onecanseeexcellentagreement
(17) of the universal marginal stability curve.
for large enough Z0, but not too close to the higher-Z0
(low-density)instabilityborder. Neartheinstabilitybor-
The asymptotics (17) is valid when ξ = fη1/2 1.
der ξ becomes of order of unity, and localization breaks
≫
It is easy to check that this criterion coincides, in the
down.
dilute limit, with the localization criterion k¯ 1 η1/2
0.06 − ≪
discussedinSec. 4A.Ontheotherhand,theparameterξ
shouldnotbetoolarge,sothatthedilutelimitcondition
Z0 1 is still satisfied, see Eq. (13). These two criteria
≫
can be rewritten as a strong double inequality for f:
0.04
2
/
−1 1 log(8η1/2)
η f ,
k η1/2 ≪ ≪ 2η1/2
0.02
that can be satisfied only for extremely large η.
4. Long-wavelength limit: perturbation theory
0
0 50 Z 100 150
0
Close to the low-density stability border, f f1
FIG. 13: The marginal stability curve for η = 105, plotted f1 = 1.19968...η−1/2 we can assume that κ −1 an≪d
in coordinates k¯ = kη−1/2 versus Z0 (solid line) and two solve Eq. (15) perturbatively. In the physic≪al units,
dilute-limit asymptotics: the dilute-limit part k¯ = A/Z0 of thisstronginequalitycorrespondstothelong-wavelength
the(“universal”)short-wavelengthcurve(dottedline)andthe limit: k L 1. In its turn, the dilute limit requires
ph x
long-wavelength asymptotics (24) (dashed line). ≪
η 1. We substitute in Eq. (15) φ(ζ)=φ0(ζ)+κΩ(ζ),
≫
10
where φ0 is givenby Eq. (16). Neglecting the κ2-termin V. SUMMARY AND DISCUSSION
Eq. (15), we obtain
2 We determined the criteria for the spontaneous
4(ζ 1)ζΩ′′+2Ω′+Ω=ζ φ0(ζ), (18)
− symmetry-breaking instability of the laterally uniform
The normalizationandboundaryconditionsatζ =1are granular cluster (strip state) in a prototypical driven
Ω(ζ = 1) = 0 and Ω(ζ = 1) = 1/2, respectively. The granulargas. Workinginthelimitofnearlyelasticparti-
′
latter condition follows from Eq. (18) itself. Equation clecollisionsandlowormoderatedensities,weemployed
(18) can be solved analytically. With the account of the granular hydrodynamics with the Jenkins-Richman con-
two boundary conditions we obtain stitutive relations [20]. The instability of the strip state
can be interpreted in terms of negative compressibility
1
Ω(ζ)= Φ2(ζ)+Φ2(ζ)I1(ζ) Φ1(ζ)I2(ζ), (19) of the granulate in the lateral direction. An important
−8 −
limit is found, when the marginalstability curvesare in-
where Φ1(ζ) = (ζ 1)1/2 and Φ2(ζ) = 2φ0(ζ). I1 and dependent of the details of the boundary condition at
I2 are indefinite in−tegrals: − the driving wall. In this regime the density perturba-
tionisexponentiallylocalizedattheelasticwallopposite
ζ Φ1G1 ζ Φ2G1 to the driving wall. Working in the dilute limit, we ob-
I1 = dζ and I2 = dζ, (20)
Z G2W Z G2W tained some analytic asymptotics of the marginal stabil-
ity curves.
where G1 = ζ2φ0/4, G2 = ζ(ζ −1), W = Φ1Φ′2−Φ′1Φ2 The results of this work show that the symmetry-
and the primes denote the ζ derivatives. Integrals I1
breaking instability predicted in Ref. [17] is robust and
andI2 canbe evaluatedinelementaryfunctions, butthe
does not require very special constitutive relations. The
results are too cumbersome to be presented here. The
marginalstability curves obtained in this work are quite
additionalboundarycondition[atζ =Z1/Z0 ζ1]reads
≡ similarto those obtainedearlier[17]for a different setof
φ(ζ1)=φ0(ζ1)+κΩ(ζ1)=0 which yields
constitutive relations (see Fig. 4). There are some quan-
κ(ξ)= φ0(ζ1)/Ω(ζ1), (21) titative differences, however. Therefore, the instability
− provides a sensitive test to the accuracy of constitutive
2
where ζ1 = cosh ξ. At the low-density stability border relations.
f =f1 we haveξ =ξ1 =1.19968....Inthe perturbative This work was focused on the criteria of instability
treatment,one shouldexpandκ(ξ) nearξ =ξ1 upto the of the strip state. In systems sufficiently long in the
lineartermξ ξ1. The zero-ordertermvanishes,andwe lateral direction, instability occurs in a whole range of
−
obtain
wave numbers k (below the respective marginalstability
Φ′2(ξ1) curve). Correspondingly, multiple steady state solutions
κ(ξ)= (ξ ξ1) . with different k are possible. In a laterally infinite sys-
−2Φ1(ξ1)I2(ξ1) −
tem,thesesolutionsareperiodicinthelateralcoordinate.
In the physical variables we have A finite system selects a finite number of wavelengths
[17]. An important issue that was not addressed in this
kphLx = η1/2 Φ′2(ξ1) 1/2(ξ ξ1)1/2. workisselection: whatisthewavelengthofthe resulting
Z0(ξ1)(cid:18)−2Φ1(ξ1)I2(ξ1)(cid:19) − symmetry-broken cluster in an infinite, or long enough,
(22) system? Theselectionhasdynamicalnature;thisimpor-
As Z0 is proportional to η1/2, the right hand side of Eq. tant issue is addressed elsewhere [22].
(22)isactuallyindependentofη. Evaluatingtheintegral
Recently, the predicted symmetry-breaking instability
I2(ξ1),weobtainI2(ξ1)= 0.883381.... Thefinalresult
− has been observed in particle simulations [23]. We hope
is
it will be investigated in experiment, too. The experi-
kphLx =2.5115...(ξ ξ1)1/2. (23) mental setting can be of the type used by Kudrolliet al.
− [7, 16]: a system of steel spheres, rolling on a smooth
Alternatively, we obtain surface and driven by a rapidly vibrating side wall. The
present work (see also Ref. [17]) provides the region of
3/2 1/2
k¯(Z0) = ξ1 Z0−Z0(ξ1) parameterswheretheinstabilitycanbeobserved. Anim-
(cid:18)Z0(ξ1)(cid:19) (cid:18) 2I2(ξ1) (cid:19) portant issue is to eliminate the static friction between
= 0.000485...(Z0(ξ1) Z0)1/2. (24) theparticlesandsurfacethatoccursfarenoughfromthe
− driving wall. In experiment, this is achieved by slightly
Theasymptotics(24)isdepictedinFig. 13. Closetothe incliningthesystem,sothataverysmallgravityappears
higher-Z0 (low-density) stability border it shows good [7,16]. Astheresult,thestripstatemovesdown,toward
agreement with the marginal stability curve found nu- the driving wall [7]. The model problem investigated in
merically. the present work does not include gravity. We expect,