Dissipationless collapses in MOND

0.39 MB·English

by  C. Nipoti

Accepted, January15,2007
PreprinttypesetusingLATEXstyleemulateapjv.10/09/06
DISSIPATIONLESS COLLAPSES IN MOND
Carlo Nipoti
Dept. ofAstronomy,UniversityofBologna,I-40127Bologna,Italy
Pasquale Londrillo
INAF-BolognaAstronomicalObservatory,I-40127Bologna,Italy
Luca Ciotti
7 Dept. ofAstronomy,UniversityofBologna,I-40127Bologna,Italy
0 Accepted, January 15, 2007
0
2 ABSTRACT

n Dissipationless collapses in Modified Newtonian Dynamics (MOND) are studied by using a new
a particle-mesh N-body code based on our numerical MOND potential solver. We found that low
J
surface-density end-products have shallower inner density profile, flatter radial velocity-dispersion
5
profile, and more radially anisotropic orbital distribution than high surface-density end-products.
1
The projected density profiles of the final virialized systems are well described by Sersic profiles with

1 index m<4, down to m 2 for a deep-MOND collapse. Consistently with observations of elliptical
v galaxies,∼the MOND end-∼products, if interpreted in the context of Newtonian gravity, would appear
8
to have little or no dark matter within the effective radius. However, we found impossible (under
1
4 the assumption of constant mass-to-light ratio) to simultaneously place the resulting systems on the
1 observedKormendy,Faber-JacksonandFundamentalPlanerelationsofellipticalgalaxies. Finally,the
0
simulations provide strong evidence that phase mixing is less effective in MOND than in Newtonian
7
gravity.
0
/ Subject headings: gravitation — stellar dynamics — galaxies: kinematics and dynamics — galaxies:
h
p elliptical and lenticular, cD – galaxies: formation — methods: numerical
-
o
r 1. INTRODUCTION throughout this paper we use
t
s y
In Bekenstein & Milgrom’s (1984, hereafter BM) La-
a µ(y)= . (5)
: grangianformulationofMilgrom’s(1983)ModifiedNew- 1+y2
v
Xi tonian Dynamics (MOND), the Poisson equation In the so-called ‘deep pMOND regime’ (hereafter
ar ∇2φN =4πGρ (1) adM),OµN(yD)),=dyesacnridbisnogelqowua-aticocnel(e2r)atsiiomnpslyifisteesmtos (k∇φk≪
0
for the Newtonian gravitationalpotential φ is replaced
N
( φ φ)=4πGa ρ. (6)
by the field equation 0
∇· k∇ k∇
φ The source term in equation (2) can be eliminated by
µ k∇ k φ =4πGρ, (2) using equation (1), giving
∇· a ∇
(cid:20) (cid:18) 0 (cid:19) (cid:21)
where a0 ≃ 1.2× 10−10ms−2 is a characteristic accel- µ k∇aφk ∇φ=∇φN+S, (7)
eration, ... is the standard Euclidean norm, φ is the (cid:18) 0 (cid:19)
k k
MOND gravitational potential produced by the density where S = curlh is a solenoidal field dependent on ρ
distribution ρ, and in finite mass systems φ 0 for and in general different from zero. When S = 0 equa-
∇ →
x . The MOND gravitationalfield g experienced tion(7)reducestoMilgrom’s(1983)formulationandcan
k k→∞
by a test particle is besolvedexplicitly. Suchreductionispossibleforconfig-
urations with spherical, cylindrical or planar symmetry,
g= φ, (3)
−∇ whicharespecialcasesofa moregeneralfamily ofstrat-
and the function µ is such that ifications (BM; Brada & Milgrom 1995). Though the
solenoidal field S has been shown to be small for some
y for y 1,
µ(y) ≪ (4) configurations(Brada& Milgrom1995;Ciotti, Londrillo
∼ 1 for y 1;
(cid:26) ≫ &Nipoti2006,hereafterCLN),neglectingitwhensimu-
Electronicaddress: [email protected] latingtime-dependentdynamicalprocesseshasdramatic
2 Nipoti, Londrillo & Ciotti
effects such as non-conservation of total linear momen- tion 4 and discussed in Section 5.
tum (e.g. Felten 1984; see also Section 3.1). 2. THEN-BODYCODE
Nowadays severalastronomical observationaldata ap-
While most N-body codes for simulations in New-
pear consistent with the MOND hypothesis (see, e.g.,
tonian gravity are based on the gridless multipole ex-
Milgrom 2002; Sanders & McGaugh 2002). In addition,
pansion treecode scheme (Barnes & Hut 1986; see also
Bekenstein(2004)recentlyproposedarelativisticversion
Dehnen2002),thenon-linearityoftheMONDfieldequa-
of MOND (Tensor-Vector-Scalar theory, TeVeS), mak-
tion(2)forcesonetoresorttoothermethods,suchasthe
ing it an interesting alternative to the cold dark matter
particle-meshtechnique(seeHockney&Eastwood1988).
paradigm. However,dynamicalprocessesinMONDhave
Inthisapproach,particlesaremovedundertheactionof
beeninvestigatedverylittlesofar,mainlyduetodifficul-
a gravitational field which is computed on a grid, with
ties posed by the non-linearity of equation (2). Here we
particle-mesh interpolation providing the link between
recall the spherically symmetric simulations (in which
the two representations. In our MOND particle-mesh
S = 0) of gaseous collapse in MOND by Stachniewicz
N-body code, we adopt a spherical grid of coordinates
&Kutschera(2005)andNusser&Pointecouteau(2006).
(r, ϑ, ϕ), made of N N N points, on which the
Theonlygenuinethree-dimensionalMONDN-bodysim- r × ϑ × ϕ
MONDfieldequationissolvedasinCLN.Particle-mesh
ulations(inwhichequation[2]issolvedexactly)arethose
interpolations are obtained with a quadratic spline in
by Brada & Milgrom (1999, 2000), who studied the sta-
eachcoordinate,while time stepping is givenby aclassi-
bility of disk galaxies and the external field effect, and
cal leap-frog scheme (Hockney & Eastwood 1988). The
those of Tiret & Combes (2007). Other attempts to
time-step ∆t is the same for all particles and is allowed
study MOND dynamical processes have been conducted
to vary adaptively in time. In particular, according to
usingthree-dimensionalN-body codesbyarbitrarilyset-
the stability criterion for the leap-frog time integration,
ting S = 0: Christodoulou (1991) investigated disk sta-
weadopt∆t=η/ max 2φ,whereη<0.3isadimen-
bility, while Nusser (2002) and Knebe & Gibson (2004) sionless parameter. We f|o∇und|that η =∼0.1 assures good
explored cosmological structure formation1. p
conservation of the total energy in the Newtonian cases
In this paper we present results of N-body simula-
(see Section 3.1). In the present version of the code, all
tions of dissipationless collapse in MOND. The simula-
the computations on the particles and the particle-mesh
tions were performed with an originalthree-dimensional
interpolations can be split among different processors,
particle-mesh N-body code, based on the numerical
while the computations relative to the potential solver
MOND potential solver presented in CLN, which solves
are not performed in parallel. The solution of equa-
equation (2) exactly. These numerical experiments are
tion (2) over the grid is then the bottleneck of the sim-
interesting both from a purely dynamical point of view,
ulations: however, the iterative procedure on which the
allowing for the first time to explore the relaxation pro-
potential solver is based (see CLN) allows to adopt as
cesses in MOND, and in the context of elliptical galaxy
seed solution at each time step the potential previously
formation. In fact, the ability of dissipationless collapse
determined.
at producing systems strikingly similar to real ellipti-
TheMONDpotentialsolvercanalsosolvethe Poisson
calsisaremarkablesuccessofNewtoniandynamics(e.g.,
equation(obtained by imposing µ=1 in equation 2), so
van Albada 1982; Aguilar & Merritt 1990; Londrillo,
Newtonian simulations can be run with the same code.
Messina & Stiavelli 1991; Udry 1993; Trenti, Bertin &
We exploited this property to test the code by running
vanAlbada 2005;Nipoti, Londrillo,& Ciotti 2006,here-
severalNewtoniansimulationsofbothequilibriumdistri-
after NLC06), while there have been no indications so
butions and collapses, comparing the results with those
far that MOND can work as well in this respect. Here
ofsimulations(startingfromthe sameinitialconditions)
we study the structural and kinematical properties of
performed with the FVFPS treecode (Londrillo, Nipoti
the end-products of MOND simulations, and we com-
& Ciotti 2003; Nipoti, Londrillo & Ciotti 2003). One of
pare them with the observed scaling relations of ellipti-
these tests is described in Section 4.1.2.
cal galaxies: the Faber–Jackson (FJ) relation (Faber &
We also verified that the code reproduces the New-
Jackson 1976), the Kormendy (1977) relation, and the
tonian and MOND conservation laws (see Section 3.1):
Fundamental Plane (FP) relation (Djorgovski & Davis
note that the conservation laws in MOND present some
1987, Dressler et al. 1987).
peculiarities with respect to the Newtonian case, so we
The paper is organized as follows. The main features
give here a brief discussion of the subject. As already
ofthenewN-bodycodearepresentedinSection2,while
stressed by BM, equation (2) is obtained from a vari-
Section 3 describes the set-up and the analysis of the
ational principle applied to a Lagrangian with all the
numerical simulations. The results are presented in Sec-
required symmetries, so energy, linear and angular mo-
1CosmologicalN-bodysimulationsinthecontextofarelativistic mentum are conserved. Unfortunately, as also shown by
MONDtheorysuchasTeVeShavenotbeenperformedsofar. BM, the total energy diverges even for finite mass sys-
Dissipationless collapses in MOND 3
in practice 2 systems of arbitrarymass and size. Each
TABLE 1 ∞
Time, velocity,andenergyunitsforNewtonian and of them is obtained by assigning specific values to the
MOND(subscript n),anddMOND (subscriptd) N-body length and mass units, r∗ and M∗, in which the initial
simulations. conditions are expressed. Also dMOND gravity is scale
free, because a appears only as a multiplicative factor
0
t∗n=r∗3/2(GM∗)−1/2 t∗d=r∗(GM∗a0)−1/4 in equation (6), and so a simulation in dMOND gravity
v∗n=(GM∗)1/2r∗−1/2 v∗d=(GM∗a0)1/4 representssystems with arbitrarymass and size (though
E∗n=GM∗2r∗−1 E∗d=(Ga0)1/2M∗3/2 inprinciple the results applyonly to systems with accel-
tems,thusposingacomputationalchallengetocodeval- erations much smaller than a ). MOND simulations can
0
idation. Wesolvedthisproblembycheckingthevolume- also be rescaled, but, due to the presence of the charac-
limited energy balance equation teristicaccelerationa inthe non-linearfunctionµ, each
0
d k+ρφ+ a20 ||∇φ|| d3x= sMim∗uclaantinoontdbeescchriobseesnoinnldyep∞en1desynsttlyemosf,eabcehcaoutsheerr.∗ and
dt 8πGF a
ZV0(cid:20) (cid:18) 0 (cid:19)(cid:21) Onthebasisoftheabovediscussion,wefixthephysical
1 ∂φ
µ < φ,nˆ >da, (8) units as follows (see Appendix C for a detailed descrip-
4πG ∂t ∇
Z∂V0 tion of the scaling procedure). Let the initial density
which is derived in Appendix A. In equation (8) V0 is distribution be characterized by a total mass M∗ and a
an arbitrary (but fixed) volume enclosing all the system characteristicradiusr∗. Werescalethefieldequationsso
mass, k is the kinetic energy per unit volume, and that the dimensionless source term is the same in New-
tonian, MOND and dMOND simulations. We also re-
y
(y) 2 µ(ξ)ξdξ, (9) quire that the Second Law of Dynamics, when cast in
F ≡
Zy0 dimensionless form, is independent of the specific force
where y0 is an arbitrary constant; note that only finite law considered, and this leads to fix the time unit. As a
quantities are involved. Another important relation be- result,NewtonianandMONDsimulationshavethesame
tween global quantities for a system at equilibrium (in time unit t∗n = r∗3/2(GM∗)−1/2, while the natural time
MOND as in Newtonian gravity) is the virial theorem unit in dMOND simulations is t∗d = r∗(GM∗a0)−1/4.
Note that MOND simulations are characterized by the
2K+W =0, (10)
dimensionless parameter κ = GM∗/r∗2a0, and scaling of
where K is the total kinetic energy and W = TrWij is aspecific simulationisallowedprovidedthe value ofκ is
the trace of the Chandrasekhar potential energy tensor maintained constant. So, simulations with lower κ val-
∂φ(x) ues describe lower surface-density, weaker acceleration
W ρ(x)x d3x (11)
ij i systems; dMOND simulations represent the limit case
≡− ∂x
Z j κ 1, while Newtonian ones describe the regime with
(e.g., Binney & Tremaine 1987). Note that in MOND ≪
κ 1. With the time units fixed, the corresponding ve-
K + W is not the total energy, and is not conserved. ≫
locity and energy units are v∗n r∗/t∗n, v∗d r∗/t∗d,
However, W is conserved in the limit of dMOND, be- E∗n =M∗v∗2n, and E∗d =M∗v∗2d≡(see Table 1 fo≡r a sum-
ing W = (2/3) Ga M3 for all systems of finite total
0 ∗ mary).
−
mass M∗ (see Appendix B for the proof). As a con-
p
sequence, in dMOND the virial theorem writes simply 3.1. Initial conditions and analysis of the simulations
σV4 = 4GM∗a0/9, where σV ≡ 2K/M∗ is the system We performed a set of five dissipationless-collapse N-
virial velocity dispersion (this relation was proved for
p body simulations, starting from the same phase-space
dMOND spherical systems by Gerhard & Spergel 1992;
configuration: theinitialparticledistributionfollowsthe
seealsoMilgrom1984). Inoursimulationswealsotested
Plummer (1911) spherically symmetric density distribu-
thatequation(10)issatisfiedatequilibrium,andthatW
tion
is conserved in the dMOND case (see Sections 3.1 and
4). ρ(r)= 3M∗r∗2 , (12)
4π(r2+r2)5/2
∗
3. NUMERICALSIMULATIONS where M∗ is the total mass and r∗ a characteristic ra-
The choice of appropriate scaling physical units is an dius. The choice of a Plummer sphere as initial con-
important aspect of N-body simulations. This is espe- dition is quite artificial, and not necessarily the most
cially true in the presentcase,in whichwe wantto com- realistic to reproduce initial conditions in the cosmolog-
pareMONDandNewtoniansimulationshavingthesame ical context (e.g., Gunn & Gott 1972). We adopt such a
initial conditions. As well known, due to the scale-free distribution to adhere to other papers dealing with col-
nature of Newtonian gravity, a Newtonian N-body sim- lisionless collapse (e.g., Londrillo et al. 1991; NLC06;
ulation starting from a given initial condition describes see also Section 5, in which we present the results of a
4 Nipoti, Londrillo & Ciotti
Fig. 1.— Time evolution of 2K/|W|, K, W, and K+W for simulations D, M1, and N. K, W, and K+W are in units of E∗d (left
column),andE∗n (central andrightcolumns). Forclarity,thetimeaxes arezoomed-inbetween 0and10.
set of simulations starting from different initial condi- and a grid with N = 64, N = 16 and N = 32) are
r ϑ ϕ
tions). Theparticlesareatrest,sotheinitialvirialratio evolved up t = 150t . In all cases the modulus of the
dyn
2K/W = 0. What is different in each simulation is center of mass position oscillates around zero with r.m.s
| |
the adopted gravitationalpotential, which is Newtonian <0.1r∗; similarly, the modulus of the total angular mo-
in simulation N, dMOND in simulation D, and MOND m∼entum oscillates around zero2 with r.m.s. <0.02, in
with acceleration ratio κ in simulations Mκ (κ=1, 2, 4). unitsofr∗M∗v∗n(simulationsMκandN)ando∼fr∗M∗v∗d
For each simulation we define the dynamical time t (simulationD).K+W intheNewtoniansimulationand
dyn
as the time at which the virial ratio 2K/W reaches its W inthedMONDsimulationareconservedtowithin2%
| |
maximum value. In particular, we find tdyn 2t∗d in and 0.6%, respectively. The volume-limited energy bal-
∼
simulation D, and tdyn 2t∗n in simulations N, M1, M2 anceequation(8)isconservedwithanaccuracyof1%in
∼
andM4. Wenotethattdyn GM∗5/2(2K+W )−3/2 in MONDsimulations,independentlyoftheadoptedV0. To
∼ | |
simulation N. estimate possible numerical effects, we reran one of the
Following NLC06, the particles are spatially dis- MOND collapse simulations (M1) using N = 2 106,
×
tributed according to equation (12) and then randomly N = 80, N = 24, and N = 48: we found that the
r ϑ ϕ
shifted in position (up to r∗/5 in modulus). This arti- end-products of these two simulations do not differ sig-
ficial, small-scale ”noise” is introduced to enhance the
2 As an experiment we also ran a simulation, with the same
phase mixing at the beginning of the collapse, because
initial conditions and parameter κ as M1, in which the force was
the numerical noise is small, and the velocity dispersion calculated from equation (7) imposing S = 0. In this simulation
is zero (see alsoSection 4.2). As such, these fluctuations thelinearandangular momentum arestronglynotconserved: for
are not intended to reproduce any physical clumpiness. instance, the center of mass is already displaced by ∼ 7r∗ after
All the simulations (realized with N = 106 particles, ∼30tdyn.
Dissipationless collapses in MOND 5
nificantly,asfarasthe propertiesrelevanttothe present generaloverviewofthetimeevolutionofthevirialquan-
work are concerned. titiesinoursimulations,postponingtoSection4.2amore
The intrinsic and projected properties of the collapse detaileddescriptionofthephase-spaceevolution. Inpar-
end-products are determined as in NLC06. In partic- ticular, in Fig. 1 we show the time evolution of 2K/W ,
| |
ular, the position of the center of the system is deter- K,W, and K+W for simulations D, M1, and N. In the
mined using the iterative technique described by Power diagramstime isnormalizedto t ,soplots referringto
dyn
et al. (2003). Following Nipoti et al. (2002),we measure different simulations are directly comparable (the values
theaxisratiosc/aandb/aofthe inertiaellipsoid(where of t in time units for the five simulations are given
dyn
a,bandcarethemajor,intermediateandminoraxis)of in Section 3.1). In simulation N (right column) we find
the final density distributions, their angle-averagedpro- thewellknownbehaviorofNewtoniandissipationlesscol-
file and half-mass radius r . We fitted the final angle- lapses: 2K/W has a peak, then oscillates, and eventu-
h
| |
averageddensityprofileswiththeγ-model(Dehnen1993; ally converges to the equilibrium value 2K/W =1; the
| |
Tremaine et al. 1994) total energy K +W is nicely conserved during the col-
lapse, though it presents a secular drift, a well known
ρ r4
ρ(r)= 0 c , (13) feature of time integration in N-body codes. The time
rγ(r +r)4−γ
c evolution of the same quantities is significantly different
where the inner slope γ and the break radius r are inadMONDsimulation(leftcolumn). Inparticular,the
c
free parameters, and the reference density ρ is fixed virial ratio 2K/W quickly becomes close to one, but is
0
| |
by the total mass M∗. The fitting radial range is still oscillating at very late times because of the oscilla-
0.06 < r/r < 10. Inordertoestimatetheimportanceof tionsofK,whileW isconstantasexpected. Asweshow
h
proje∼ctioneff∼ects,foreachend-productweconsiderthree inSection4.2, these oscillationsare relatedto a peculiar
orthogonalprojectionsalongtheprincipalaxesofthein- behavior of the system in phase space. Finally, simula-
ertia tensor, measuring the ellipticity ǫ = 1 b /a , the tionM1(centralcolumn)representsanintermediatecase
e e
−
circularizedprojecteddensityprofileandthecircularized betweenmodels N andD: the systemstarts as dMOND,
effectiveradiusR √a b (wherea andb arethema- but soon its core becomes concentrated enough to enter
e e e e e
≡
jorandminorsemi-axisoftheeffectiveisodensityellipse). the Newtonian regime. After the initial phases of the
We fit (over the radial range 0.1 < R/R < 10) the cir- collapse, Newtonian gravity acts effectively in damping
e
cularized projected density profil∼es of the∼end-products the oscillationsofthe virialratio. Overall,it is apparent
with the R1/m Sersic (1968) law: how the system is in a “mixed” state, neither Newto-
nian (K +W is not conserved) nor dMOND (W is not
1/m
R
constant).
I(R)=I exp b(m) 1 , (14)
e
(− "(cid:18)Re(cid:19) − #)
4.1. Properties of the collapse end-products
whereI I(R )andb(m) 2m 1/3+4/405m(Ciotti
e e
≡ ≃ − 4.1.1. Spatial and projected density profiles
&Bertin1999). Inthefittingproceduremistheonlyfree
parameter, because Re and Ie are determined by their We found that all the simulated systems, once viri-
measuredvaluesobtainedbyparticlecount. Inaddition, alized, are not spherically symmetric. However, while
we measure the central velocity dispersion σ0, obtained the dMOND collapse end-product is triaxial (c/a 0.2,
∼
by averaging the projected velocity dispersion over the b/a 0.4), MOND and Newtonian end-products are
∼
circularized surface density profile within an aperture of oblate (c/a c/b 0.5). The ellipticity ǫ of the pro-
∼ ∼
Re/8. Some of these structural parameters are reported jected density distributions (measured for each of the
in Table 2 for the five simulations described above, as principal projections) is found in the range 0.5 0.8 in
−
wellasforthreeadditionalsimulations,whichstartfrom D, and 0 0.5 in M1, M2, M4 and N. These values are
−
different initial conditions (see Section 5). consistentwiththoseobservedinrealellipticals,withthe
exception of ǫ in model D (see Table 2), which would
b
4. RESULTS correspond - if taken at face value - to an E8 galaxy.
In Newtonian gravity, collisionless systems reach viri- These result could be just due to the procedure adopted
alization through violent relaxation in few dynamical to measure the ellipticity (see Section 3.1), however we
times, as predicted by the theory (Lynden-Bell 1967) find interesting that dMOND gravity could be able to
and confirmed by numerical simulations (e.g. van Al- produce some system that would be unstable in New-
bada 1982). Onthe other hand, due to the non linearity tonian gravity. We remark that a similar result, in the
of the theory and the lack of numerical simulations, the different context of disk stability in MOND, has been
detailsofrelaxationprocessesandvirializationinMOND obtained by Brada & Milgrom (1999).
are much less known. Thus, before discussing the spe- Inordertodescribethe radialmassdistributionofthe
cificpropertiesofthe collapseend-productswepresenta final virialized systems, we fitted their angle-averaged
6 Nipoti, Londrillo & Ciotti
TABLE 2
End-product properties.
κ c/a b/a γ rc/rh ma mb mc ǫa ǫb ǫc
D - 0.21 0.41 0.17+0.39 0.27+0.09 2.87±0.01 2.50±0.03 2.16±0.02 0.48 0.80 0.58
−0.17 −0.02
M1 1 0.47 0.85 1.24+0.40 0.44+0.20 3.20±0.07 3.00±0.09 3.07±0.13 0.42 0.51 0.17
−0.36 −0.12
M2 2 0.48 0.92 1.45+0.26 0.53+0.19 3.38±0.08 3.24±0.08 3.28±0.12 0.49 0.45 0.08
−0.34 −0.13
M4 4 0.47 0.90 1.54+0.26 0.58+0.20 3.55±0.10 3.40±0.11 3.34±0.15 0.51 0.45 0.10
−0.32 −0.15
N - 0.45 0.91 1.69+0.13 0.74+0.14 4.21±0.07 4.35±0.08 3.96±0.13 0.48 0.55 0.12
−0.15 −0.12
D′ - 0.25 0.45 0.72+0.36 0.37+0.12 3.06±0.06 2.90±0.04 2.71±0.08 0.44 0.76 0.56
−0.50 −0.11
M′ 20 0.42 0.83 1.26+0.44 0.47+0.25 3.41±0.09 3.36±0.06 3.20±0.13 0.49 0.57 0.16
−0.40 −0.15
N′ - 0.45 0.93 1.78+0.15 0.78+0.24 4.29±0.10 4.56±0.15 4.19±0.22 0.51 0.55 0.09
−0.18 −0.18
Firstcolumn: nameofthesimulation. κ=GM∗/r∗2a0: acceleration ratio. c/aandb/a: minor-to-majorandintermediate-to-majoraxis
ratios. γ,rc: best-fitγ-modelparameters. ma,mb,mc andǫa,ǫb,ǫc: best-fitSersicindicesandellipticitiesforprojectionsalongthe
principalaxes.
Fig. 2.—Angle-averageddensity,radialvelocity-dispersionandanisotropy-parameterprofiles(frombottomtotop)fortheend-products
of simulations D,M1, and N. Dotted lines inthe bottom panels represent ρ∝r−1 profiles, whichareshown for reference. Empty circles
intherightcolumnshowthecorrespondingprofilesobtained withtheFVFPStreecodefromthesameinitialconditions.
Dissipationless collapses in MOND 7
Fig. 3.— Line-of-sightvelocity-dispersionprofiles(top),circularizedprojecteddensityprofilesandresidualsoftheSersicfit(bottom)of
theend-products ofsimulationsD,M1,andN(squares;1-σ errorbarsarealwayssmallerthanthesymbolsize). Thedotted linesarethe
best-fittingSersicmodels.
8 Nipoti, Londrillo & Ciotti
density profiles with the γ-model (13) over the radial anisotropic for r>r . For each model projection we
h
range 0.06 < r/r < 10. The best-fit γ and r for the computed the line∼-of-sight velocity dispersion σ , con-
h c los
final distrib∼ution of∼each simulation are reported in Ta- sidering particles in a strip of width R /4 centered on
e
ble 2 together with their 1σ uncertainties (calculated the semi-major axis of the isophotal ellipse. The line-of-
from ∆χ2 = 2.30 contours in the space γ r ). As sightvelocity-dispersionprofiles (for the major-axispro-
c
−
also apparent from Fig. 2 (bottom), the Newtonian col- jection)areplottedinthetoppanelsofFig.3. TheNew-
lapseproducedthesystemwiththesteepestinnerprofile tonian profile is very steep within R , while MOND and
e
(γ 1.7), the dMOND end-product has inner logarith- dMONDprofilesaresignificantlyflatterthere. Aswellas
∼
mic slope close to zero, while MOND collapses led to σ ,σ decreasesfordecreasingradiusintheinnerregion
r los
intermediate cases, with γ ranging from 1.2 (κ = 1) ofmodelD.ThekinematicalpropertiesofM2andM4are
∼
to 1.5 (κ = 4). We also note that the ratio r /r intermediate between those of M1 and of N: overall we
c h
∼
(indicating the position of the knee in the density pro- find only weak dynamical non-homology among MOND
file)increasessystematicallyfromdMONDtoNewtonian end-products. The empty symbols in Fig. 2 (right col-
simulations. umn) refer to a test Newtonian simulation run with the
The circularized projected density profiles of the end- FVFPS treecode (with 4 105 particles). The struc-
×
products are analyzed as described in Section 3.1. The tural and kinematical properties of the end-product of
best-fit Sersic indices m , m and m (for projections this simulation are clearly in good agreementwith those
a b c
along the axes a, b, and c, respectively) are reported in of the end-product of simulation N (solid lines), which
Table2,togetherwiththe1σuncertaintiescorresponding started from the same initial conditions.
to∆χ2 =1;therelativeuncertaintiesonthebest-fitSer-
4.2. Phase-space properties of MOND collapses
sicindicesareinallcasessmallerthan5percentandthe
averageresiduals between the data and the fits are typi- To explore the phase-space evolution of the systems
cally0.05< ∆SB <0.2,whereSB 2.5log[I(R)/Ie]. during the collapse and the following relaxation we con-
The fitting∼rhadialira∼nge 0.1 < R/Re≡<−10 is comparable sider time snapshots of the particles radial velocity (vr)
with or larger than the typ∼ical rang∼es spanned by ob- vs. radius as in Londrillo et al. (1991). In Fig. 4 we plot
servations (e.g., see Bertin, Ciotti & Del Principe 2002). fiveofthesediagramsforsimulationsD,M1andN:each
Inagreementwithpreviousinvestigations,wefoundthat plotshowsthephase-spacecoordinatesof32000particles
the Newtoniancollapseproduceda systemwellfitted by randomly extracted from the corresponding simulation,
the de Vaucouleurs (1948) law. MOND collapses led to and, as in Fig. 1, times are normalizedto the dynamical
systems with Sersic index m < 4, down to m 2 in the time t (see Section3.1). At time t=0.5t allparti-
∼ dyn dyn
case of the dMOND collapse. Figure 3 (bottom) shows cles are still collapsing in simulation N, while in MOND
the circularized (major-axis) projected density profiles simulations a minority of particles have already crossed
fortheend-productsofsimulationsD,M1andNtogether thecenterofmass,asrevealedbytheverticaldistribution
with their best-fit Sersic laws (m =2.87, m = 3.20, and ofpointsatr 0inpanelsDandM1. Att=t (time
dyn
∼
m=4.21,respectively),andthecorrespondingresiduals. ofthe peak of2K/W inthe three models), sharpshells
| |
Curiously,NLC06 foundthat low-msystems canbe also in phase space are present, indicating that particles are
obtained in Newtonian dissipationless collapses in the movinginandoutcollectively andphasemixing has not
presence of a pre-existing dark-matter halo, with Sersic taken place yet. At t = 4t is already apparent that
dyn
index value decreasing for increasing dark-to-luminous phase mixing is operating more efficiently in simulation
mass ratio. N than in simulation M1, while there is very little phase
mixing in the dMOND collapse. At significantly late
4.1.2. Kinematics times (t = 44tdyn), when the three systems are almost
virialized(2K/W 1;seeFig.1),phasemixingiscom-
We quantify the internal kinematics of the collapse | |∼
pleteinsimulationN,butphase-spaceshellsstillsurvive
end-productsbymeasuringtheangle-averagedradialand
in models M1 and D. Finally, the bottom panels show
tangential components (σ and σ ) of their velocity-
r t the phase-space diagrams at equilibrium (t = 150t ),
dyn
dispersion tensor, and the anisotropy parameter β(r)
≡ when phase mixing is completed also in the MOND and
1 0.5σ2/σ2. These quantities are shown in Fig. 2 for
− t r dMOND galaxies: note that the populated region in the
simulations D, M1, and N. We note that the σ pro-
r (r,v ) space is significantly different in MOND and in
r
file decreases more steeply in the Newtonian than in the
Newtoniangravity,consistently with the sharper decline
MONDend-products,whileitpresentsaholeintheinner
of radial velocity dispersion in the Newtonian system.
regionsofthedMONDsystem. Inaddition,thedMOND
Thus, our results indicate that phase mixing is more
galaxy is radially anisotropic (β 0.4) even in the cen-
∼ effectiveinNewtoniangravitythaninMOND3.Itisthen
tral regions, where models N and M1 are approximately
isotropic (β 0.1). All systems are strongly radially 3 Ciotti, Nipoti & Londrillo(2007) found similarresults in“ad
∼
Dissipationless collapses in MOND 9
Fig. 4.—Phase-space(radial-velocityvs. radius)diagramsforsimulationsD,M1,andNatvarioustimes. vr isinunitsofv∗d(simulation
D),andv∗n (simulationsM1andN;seeTable1).
interestingtoestimateinphysicalunitsthephase-mixing tional systems is the energy distribution N(E) (i.e. the
timescales of MOND systems. From Table 1 it follows number of particles with energy per unit mass between
that t∗n 4.7(r∗/kpc)3/2(M∗/1010M⊙)−1/2Myr = E and E+dE; e.g., Binney & Tremaine 1987; Trenti &
≃
29.8κ−3/4(M∗/1010M⊙)1/4Myr for a0 = 1.2 Bertin2005). Independently of the force law,the energy
×
10−10ms−2. For example, in the case of model M1, perunitmassofaparticleorbitingatxwithspeedvina
adoptingM∗ =1012M⊙ (andr∗ = GM∗/a0 34kpc), gravitationalpotentialφ(x)isE =v2/2+φ(x),andE is
≃
shells in phase space are still apparent after 8.3Gyr constant if φ is time-independent. In Newtonian gravity
p ∼
( 44t ). Simulation M1 might also be interpreted φisusuallysettozeroatinfinityforfinite-masssystems,
dyn
≃
as representing a dwarf elliptical galaxy of, say, so E <0 for bound particles; in MOND all particles are
M∗ = 109M⊙ (and r∗ = GM∗/a0 1.1kpc). In this bound, independently of their velocity,becauseφ is con-
≃
case44t 1.5Gyr. WeconcludethatinsomeMOND fining, and all energies are admissible. This difference
dyn ∼ p
systems substructures in phase space can survive for is reflected in Fig. 5, which plots the initial (top) and
significantly long times. final (bottom) differential energy distributions for simu-
In addition to the (r,v ) diagram, another useful di- lationsD,M1,andN.Giventhattheparticlesareatrest
r
agnostictoinvestigatephase-spacepropertiesofgravita- at t=0, the initial N(E) depends only on the structure
of the gravitational potential, and is significantly differ-
hoc”numericalsimulationsinwhichtheangularforcecomponents
ent in the Newtonian and MOND cases. We also note
werefrozentozero,sothattheevolutionwasdrivenbyradialforces
only. In fact, while phase mixing is less effective both in MOND that N(E) is basically the same in models D and M1 at
andinNewtoniansimulationswithrespecttothesimulationshere t =0, because model M1 is initially in dMOND regime.
reported,thephasemixingtimescaleinMONDisstillconsiderably
In accordance with previous studies, in the Newtonian
longerthaninNewtoniangravity.
case the final differential N(E) is well representedby an
10 Nipoti, Londrillo & Ciotti
Fig. 5.—Initial(top)andfinal(bottom)differentialenergydistributions. TheenergyperunitmassE isinunitsofE∗d/M∗ (modelD),
and E∗n/M∗ (models M1 and N). The energy zero points inmodels D and M1 are such that the most bound particles of the M1 and N
end-products havethesameenergy,andthehighestenergyparticlesofmodelsDandM1att=0havethesameenergy.
exponential function over most of the populated energy Sanders 2000). For example, when interpreting the FP
range (Binney 1982; van Albada 1982;Ciotti 1991;Lon- tilt in Newtonian gravity one can invoke a systematic
drillo et al. 1991; NLC06). In contrast, in model D the and fine-tuned increase of the galaxy dark-to-luminous
final N(E) decreases for increasing energy, qualitatively mass ratio with luminosity (e.g., Bender, Burstein &
preservingits initialshape. Inthe caseofsimulationM1 Faber 1992; Renzini & Ciotti 1993; Ciotti, Lanzoni &
it is apparent a dichotomy between a Newtonian part Renzini1996),whileinMONDthe tiltshouldbe related
at lower energies (more bound particles), where N(E) tothecharacteristicaccelerationa . Note,however,that
0
is exponential, and a dMOND part at higher energies, in MOND as well as in Newtonian gravity other impor-
wherethefinalN(E)resemblestheinitialone. Weinter- tantphysicalproperties may help to explainthe FP tilt,
pretthisresultasanothermanifestationofalesseffective suchasa systematicincreaseofradialorbitalanisotropy
phase-spacereorganizationinMONDthaninNewtonian with mass or a systematic structural weak homology
collapses. (Bertin et al. 2002). Due to the relevance of the sub-
ject, we attempt here to derive some preliminary hints.
4.3. Comparison with the observed scaling relations of In particular, for the first time, we can compare with
the scaling relations of elliptical galaxies MOND models
elliptical galaxies
producedbyaformationmechanism,yetassimpleasthe
It is not surprising that galaxy scaling relations rep-
dissipationless collapse.
resent an even stronger test for MOND than for New-
In this Section we consider the end-products of sim-
tonian gravity, due to the absence of dark matter and
ulations M1, M2, and M4. As already discussed in
the existence of the critical acceleration a with a uni-
0 Section 3, each of the three systems corresponds to a
versalvalueintheformertheory(e.g.,seeMilgrom1984;