Stellar Bar evolution in the absence of dark matter halo
Mahmood Roshan
astro-ph.GA gr-qc
  

 





apj
packed_enum

  
  
  











ad hoc
cf.
e.g.,
et al.
etc.
i.e.
vs.
BT08


#10=#1
  -0.25em0-0
  .05em0-0
  -0.025em.0433em0
 a
 x
 v

#1
#1#1









Department of Physics, Ferdowsi University of Mashhad, P.O. Box 1436, Mashhad, Iran; mroshan@um.ac.ir

We study the stellar bar growth in high resolution numerical galaxy models with and without dark matter halos. In all models the galactic disk is exponential and the halos are rigid or live Plummer spheres. More specifically, when there is no dark matter halo, we modify the gravitational force between point particles. To do so we use the weak field limit of an alternative theory of dark matter known as MOG in the literature. The galaxy model in MOG has the same initial conditions as in galaxy models with dark matter halo. On the other hand, the initial random velocities and the Toomre's local stability parameter are the same for all the models. We show that the evolution and growth of the bar in MOG is substantially different from the standard cases including dark matter halo. More importantly, we find that the bar growth rate and its final magnitude is smaller in MOG. On the other hand, the maximum value of the bar in MOG is smaller than the Newtonian models. It is shown that although the live dark matter halo may support the bar instability, MOG has stabilizing effects. Furthermore, we show that MOG supports fast pattern speeds, and unlike in the dark matter halo models pattern speed does not decrease with time. Theses differences, combined with the relevant observations, may help to distinguish between dark matter an modified gravity in galactic scales.  





§ INTRODUCTION


It is well-known that bulge-less numerical galactic models which initially are in an equilibrium rotational state, are violently unstable to a pressure dominated system, for example for primary simulations see <cit.>; <cit.>. This fact is known as the bar instability in the literature. More specifically, the disk changes its shape to a rotating bar in a small time scale compared to the rotational period of the outermost particles. However, it was understood that existence of a rigid spherical component can suppress the bar instability <cit.>. After Ostriker & Peebles influential work, the role of dark matter halo on the global stability of the galaxies has been extensively investigated, for example see <cit.>; <cit.>;
<cit.>; <cit.> and <cit.> for a recent review of the subject. 


Although, historically, the dark matter halo introduced to prevent the bar formation, observations show that one-third of disk galaxies are strongly barred. On the other hand, one-third of them are weakly barred, and the remainder do not contain bar. In fact one of the main unsolved problem with these galaxies is the frequency of bars. Several attempts have been done to explain the presence or absence of a bar, for very short review on the subject see <cit.>. In this paper we do not attempt to present a new explanation for the observed bar frequency. We are interested to bar evolution in models when there is no dark matter halo and, instead, the gravitational force is different from the Newtonian case. 

It is important to mention that the bar instability is directly linked to the dark matter problem. On the other hand modified gravity is another path to address this problem. Therefore, the global stability of disk galaxies is an important issue which might help to distinguish between dark matter and modified gravity. In fact, the dark matter particles have not yet been observed while there are several laboratories throughout the world searching for these particles using completely different techniques, see <cit.>. This means that modified gravity is still an important approach to address the problem. Theses theories can be divided to two main categories: dark energy models and alternative theories which dispense with the dark matter paradigm. For a review on dark energy models we refer the reader to <cit.>. On the other hand Modified Newtonian Dynamics (MOND) is one of the well studied modified theories for dark matter (; ). Scalar-tensor-vector theory of gravity known as MOG in the literature is another theory which presented to address the dark matter problem <cit.>. 

Before explaining our aim in this paper, let us briefly mention some well established results on the stability of disk galaxies in both viewpoints, i.e. in particle dark matter and modified gravity. It is well understood that disks in live dark matter halos form bars more readily than in rigid halos <cit.>. In fact angular momentum transfer between the bar and the halo, in principle, can increase the growth rate. This means that the so-called disk stability criteria that derived using rigid halos <cit.> can  not be used for real galaxies <cit.>. The rotation of the halo can also influence the growth rate. For example, <cit.> reported that if the halo rotates in the same direction as the disk, the bar growth will increase. On the other hand, the counter-rotating halos may suppress the instability. 

Another important result of stellar N-body simulations is that the bar pattern speed decreases with time <cit.>. However observations do not confirm this point <cit.>. Therefore it is still a challenge for standard dark matter halo models.

On the other hand, less investigation has been done on the stability of disk galaxies in modified gravity theories. The local stability of disk galaxies has been investigated in <cit.>. The global stability of disk galaxies in MOND has been investigated using N-body simulations in <cit.>; <cit.> and <cit.>. Using a low resolution simulation, it has been shown in <cit.> that disk galaxies are more stable in MOND than equivalent models in Newtonian gravity. A similar result has been reported in <cit.>. However, advanced simulations with larger number of particles show that the bar growth rate in MOND is larger than in Newtonian models <cit.>. This means that galactic disks without dark matter halo in MOND develop a bar sooner than in Newtonian models with dark matter halo. However, in MOND the bar starts to weaken while in the dark matter model it does not stop growing. In this case, the final magnitude of the bar at the end of the simulation in MOND is smaller than the dark matter halo models. 

<cit.> investigated the local stability of the galactic disks in MOG. In other words, the generalized version of the Toomre criterion has been found and applied to a sample of spiral galaxies. The result shows that there is no significant difference between MOG and the standard picture in the local stability issue. More specifically, the effects of MOG appear at large distances and, naturally, does not affect the dynamics of local perturbations. <cit.> investigated the global stability of the Maclaurin disk in MOG using semi analytic method. On the other hand the global stability of exponential and Mestel disk galaxies in MOG has been investigated in <cit.> using a low resolution N-body simulations. Although there is no meaningful difference in the local stability of disk galaxies in MOG and the dark matter model, the global stability is substantially different. More specifically, it is found in <cit.> that disks are more stable in MOG and the growth rate is smaller. 

However the maximum number of particles in the simulations of <cit.> is N=2× 10^4. Therefore the results may suffer from shot noise and numerical artifacts. Although some tests have been done to ensure the reliability of the results, it still seems necessary to perform high resolution simulations. Furthermore galactic models in <cit.> do not include dark matter halos. On the other hand, dark matter halo is one of the main components of a spiral galaxy in the standard picture. Therefore, without adding dark matter halo it is not possible to compare appropriately galactic dynamics in MOG with the standard case. Therefore we generalize simulations of <cit.> by using another N-body code which allows large number of particles. In fact in this paper we will use 2× 10^6-10^7 particles for each live component. On the other hand we add the dark matter component and make more realistic galactic models. 



The outline of this paper is as follows: In section <ref> we briefly introduce MOG and discuss its weak filed limit. In section <ref> we describe the numerical method for setting the initial conditions and evolving the point particles. In fact we make three models which have identical initial conditions, one model in MOG without dark matter halo and two models in the standard case including dark matter halo. In section <ref> we compare the bar evolution in these three models and discuss the significant differences and their relevance to current observations. Finally we discuss the future developments.



§ WEAK FIELD LIMIT OF MOG



MOG is a relativistic generalization of General Relativity (GR) which beside the metric tensor has two scalar fields, i.e. G and μ, and one massive Proca vector field ϕ^α. In other words MOG is a scalar-tensor-vector theory of gravity. In particle physics language, MOG possesses a spin two graviton, i.e. the metric tensor g_μν, a spin one massive graviton, namely the vector field ϕ^α, and a spin zero massless graviton, i.e. the scalar field G.

Existence of these fields enables the theory to address the dark matter problem without invoking dark matter particles. More specifically, in the weak field limit, MOG explains the flat rotation curve of spiral galaxies <cit.>, as well as the mass discrepancy in galaxy clusters <cit.>. This theory has been applied also to Bullet Cluster and Train Wreck Cluster in order to explain the unusual lensing data, see <cit.>. On the other hand the cosmological aspects of the theory has been widely investigated, for example see <cit.>; <cit.>; <cit.>; <cit.>; <cit.> and <cit.>. It is important mentioning that MOG has a viable sequence of cosmological epochs. On the other hand, in its original form, it can not be considered as a dark energy model. In other words, this theory is designed to address the dark matter problem and uses a cosmological constant Λ to be consistent with the dark energy observations. 



In order to study the dynamics of a disk galaxy, the weak field limit of the theory is required. The weak field limit of MOG has been studied in
<cit.> and <cit.>. In the following we briefly review the weak field limit. The field equations of the theory, which are obtained from varying the generic action of the theory with respect to the fields, i.e. the metric tensor g_μν, G, ϕ_α, μ, are given by 
 
    G_μν+g_μν(2/G^2∇_αG∇^αG-□ G/G+Gωμ^2 ϕ_αϕ^α
          -4π∇_αμ∇^αμ/μ^2)-2∇_μG∇_νG/G^2 + ∇_μ∇_νG/G
          -8π(1/4πGωμ^2ϕ_αϕ^α-∇_μG∇_νG/G^2-
    ∇_μμ∇_νμ/μ^2)
          
    +Gω/4π(B^α_ μB_να+1/4B^αβB_αβg_μν)=-8π G T_μν

 
    □ G+G/2(∇^αμ∇_αμ/μ^2+R/8π)
        -3/2G∇^αG∇_αG=0

 
    ω∇_μB^μν+ωμ^2ϕ^ν=4π J^α

 
    □μ-∇_αμ∇^αμ/μ-∇_αμ∇^αG/G+1/4πGωμ^3ϕ^αϕ_α=0

where □=∇_μ∇^μ, B_μν=∇_μϕ_ν-∇_νϕ_μ, G_μν is the Einstein tensor, R is the Ricci scalar, T_μν is the energy-momentum tensor, ω and κ are coupling constants, J^α is a "matter current" obtained from the variation of the matter action with respect to the vector field, and its time component is written as J^0=κρ. Only in this section, we use units in which the speed of light is unity. In order to find the weak field limit, let us perturb the Minkowski metric η_μν as follows
 
    g_μν=η_μν+h_μν=(
                  [       1+2 Φ           0;           0 (-1+2Ψ)δ_ij;             ])

where |h_μν|≪ 1. On the other hand, using the field equations it is easy to show that the background values of J^α and ϕ^α are zero. Therefore, the linear perturbation of other fields and the matter current can be written as
  
    ϕ^μ=ϕ_1^μ(𝐫,t),      G=G_0+G_1 (𝐫,t)
       μ=μ_0+μ_1(𝐫,t)                   J^μ=J_1^μ(𝐫,t)

 where subscripts "0" and "1" stand for background and perturbed quantities respectively. Now, assuming a perfect fluid energy-momentum tensor and substituting equations (<ref>) and (<ref>) into the field equations (<ref>)-(<ref>), we find the linearised field equations

    ∇^2 Ψ=4π G_0 (16π-2/16π-3)ρ


    ∇^2 Φ=4π G_0 (16π/16π-3)ρ


    ∇^2ϕ_1^μ-μ_0 ϕ_1^μ=-4π/ω J_1^μ


    ∇^2 G_1=8π G_0^2/16π-3ρ


    ∇^2 μ_1=0

Now in order to find the modified version of the Poisson equation, we need to linearise the equations of motion of a test particle, or equivalently the geodesic equation, in MOG. the geodesic equation in MOG differs from GR and is given by 

    d^2x^i/dτ^2+Γ^i_αβdx^α/dτdx^β/dτ=κ B^i_ αdx^α/dτ

in which Γs are the Christoffel symbols. It is easy to show that in the weak field limit this equation takes the following form  

    d^2𝐫/dt^2≃-∇Φ_eff

where the effective potential Φ_eff is defined as

    Φ_eff=Φ+κϕ_1^0

On the other hand, equations (<ref>) and (<ref>) can be straightforwardly solved as

    Φ=-16π G_0  /16π-3∫ρ(𝐫')/|𝐫-𝐫'|d^3x'


    ϕ_1^0=κ/ω∫e^-μ_0|𝐫-𝐫'|/|𝐫-𝐫'|ρ(𝐫')d^3x'

It is clear from Yukawa equation (<ref>) that the scalar field μ appear as a mass for the vector field. Using these solutions, we rewrite equation (<ref>) as 

    Φ_eff(𝐫)=-G_N∫ρ(𝐫')/|𝐫-𝐫'|(1+α-α e^-μ_0|𝐫-𝐫'|)d^3x'

where α is a dimensionless parameter defined as α=κ/ω G_N,   and G_N is the Newton's gravitational constant. Also, in order to recover the Newtonian gravity in small distances, it is necessary to set 16π G_0/(16π-3)-κ^2/ω=G_N. Now using equations (<ref>), (<ref>), (<ref>) and (<ref>), one can readily verify that the Poisson equation in Newtonian gravity is replaced with 

    ∇^2 Φ_eff(𝐫)=4π G_N ρ +αμ_0^2 G_N ∫e^-μ_0 |𝐫-𝐫'|/|𝐫-𝐫'|ρ(𝐫')d^3x'


Hereafter, for simplicity, we drop the subscript "N" in the Newton's constant G_N. As it is clear we have two free parameters α and μ which should be fixed using relevant dark matter observations. From rotation curve data of spiral galaxies we have α=8.89 ± 0.34 and μ_0=0.042 ± 0.004  kpc^-1 <cit.>.  

Now let us find the gravitational potential of a test particle M located at 𝐫'=0. In this case the density is given by ρ(𝐫')=M δ(𝐫'), where δ is the Dirac delta function. By substituting this density into equation (<ref>) we arrive at

    Φ_eff(r)=-G M/r[1+α-α e^-μ_0 r].

 Consequently, the gravitational force experienced by another unit mass located at distance r from the point mass M is given by

    𝐅=-G M/r^2[1+α-α(1+μ_0 r)e^-μ_0 r]𝐫/r

At small radii the Newtonian force is recovered. On the other hand at larger scales the Yukawa term is suppressed and the the gravitational force in MOG is enhanced by factor (1+α) compared to the Newtonian force. 


As in other alternative theories of  gravity for dark matter, the gravitational force is stronger than the Newtonian force. It is interesting to mention that in strong field limit, MOG predicts new effects which are absent in GR. For example, MOG leads to gravito-electrical and gravito-magnetic effects which are directly related to the existence of the vector field, for more details see <cit.> and <cit.>. However, we do not expect this effect in the galactic scale where almost all the dynamics can be explained in the Newtonian limit.
ccccccccccccccc
15
Model information

    Disk     Disk     Disk     Initial      Halo     Halo     Halo   Halo   Time   Cylindrical      Spherical   Softening 

Run     Type     Edge     Mass   Q   Type   Mass    Scale     r_ max   Step   Polar grid   δ z   Grid   Length   



EN200     Exp   4     1     1.5     none     ...     ...     ...     0.01     390× 407× 135   0.08   ...   0.16   

EM200     Exp   4     1     1.5     none     ...     ...     ...     0.01     193× 224× 125   0.08   ...   0.16   

EPR200     Exp   4     1     1.5     Plum     7.5     5.5     ...     0.01     99× 128× 125   0.08   ...   0.16   

EPL200     Exp   4     1     1.5     Plum     8     10     30     0.01     100× 128× 125   0.08   1001   0.16   

Column 1: Letter identification for the simulation. Column 2: Type
   of initial disk; “Exp” is exponential disk. Column 3:
   Initial outer radius of the disk in units of R_d.  Column 4:
   Initial disk mass M_d.  Column 5: Initial Toomre's Q parameter.  Column 6: Type of halo, “Plum” is the Plummer isotropic sphere.  Column 7: Mass of the halo
   component M_b.  Column 8: Radial scale of the halo
   component b in units of R_d. Column 9: Outer edge of halo, r_
     max in units of R_d. Column 10: the basic time step in units of τ_0. Column 11: number of rings, spokes, and planes in the cylindrical polar grid. Column 12: vertical distance between grid planes. Column 13: the number of shells in the spherical grid. This grid extends to R=180 R_d. Column 14: the gravity softening length in units of R_d..
 





§ NUMERICAL METHOD


We use the N-body code developed by <cit.>. For galaxy model with live halo, we use the hybrid N-body code described by <cit.> in which two grids are used to compute the gravitational field: a 3D cylindrical polar grid for the disk component, and a spherical grid for the live halo component. One may find the full description of the code in the online manual <cit.>. 

For galaxy models in MOG, we have modified the code in order to take into account the modified gravity corrections to the gravitational force. In our galaxy model in MOG, as expected, there is no dark matter halo. The energy conservation has been tested, and the error in total energy is always smaller than 1%. 


 §.§ Initial conditions



All the models we present in this paper start with an exponential disk given by the following surface density

    Σ(R)=M_d/2π R_d^2 e^-R/R_d

where M_d is the disk mass, and R_d is its length scale. The extend of the disk is limited by tapering the surface density with a function that varies as a cubic polynomial from unity at R=3.2  R_d to zero at R=4  R_d. Using the epicyclic approximation, in which particles move on nearly circular orbits, we set the radial velocity dispersion σ_R such that to ensure about the local stability of the disk. In fact the local stability criterion is written as Q>1 where the Toomre's stability parameter Q <cit.> is given by 

    Q(R)=σ_R/σ_R,crit      where    σ_R,crit=3.36 G Σ/κ

and κ(R) is the epicyclic frequency. Using the epicycle relation, we also find the initial azimuthal velocity dispersion as σ_ϕ=κσ_R/2Ω, where Ω(R)=v_c/R and v_c is the circular velocity. On the other hand the disk has a Gaussian density profile in the vertical direction with the scale height z_0=0.05 R_d. Furthermore, the vertical velocities are determined by integrating the vertical one dimensional Jeans equation as

    σ_z^2(R,z)=1/ρ(R,z)∫_z^∞ρ(R,z')∂Φ/∂ z'dz'


For the halo component, we have chosen the isotropic Plummer sphere, which has the following density profile

    ρ(r)=3 M_h/4π b^3[1+(r/b)^2]^-5/2

where M_h is the halo mass and b is the characteristic radius. The values of these parameters are selected so that the Newtonian model has the same rotation curve as the disk in modified gravity (which has no spherical component). The halo is truncated at r=3 b. It is important mentioning that this mass distribution is consistent with the rising rotation curves in late type galaxies of lower luminosity <cit.>. 

It is necessary to mention that we have two main reasons for choosing this halo. In fcat with this choice, it is easy to match the rotation curves of the dark matter and modified gravity models. On the other hand  this profile has been used also in <cit.> in order to compare bulge formation in MOND and standard Newtonian dynamics. Also it has been used in <cit.> to study the stellar bar evolution in MOND. Therefore this choice of halo enables us to make a direct comparison with  <cit.>, and conseuently a direct comarison between MOND and MOG. 



 <cit.> gives expressions for the distribution function (DF) of the Plummer sphere, characterized by a parameter q. We set q=0 which correspond to a isotropic sphere. Therefore we start with a sphere with known DF when it is in isolation. However, when a disk component is added, one needs to find the equilibrium DF of the spherical component in the composite halo and disk model. To do so we use a halo compression algorithm described by <cit.> and <cit.>. This algorithm uses the fact that both radial action and the total angular momentum are conserved during compression. Finally the particles will be selected from the new DF function as described in <cit.>. For more details see <cit.> in which the iterative solution for a single spheroid with an added disk were studied using the above mentioned procedure. 
 
We emphasize that no bulge is included in our models. More specifically, we study two single component models "EN" and "EM" with an exponential disk. EN is the exponential disk model in Newtonian gravity and EM is the corresponding model in MOG. On the other hand "EPR" and "EPL" are two component models including an exponential disk and a Plummer halo. In EPR the halo is rigid, and in EPL is live and responsive. The properties of these models have been summarized in Table <ref>. For each responsive components we employ two million particles. For example, in EPL we employed two million particle for the disk and two million particle for the live halo. Other models have two million particles. In order to check that the results do not depend on the number of particles, we also increase the number of particles by a factor of 10 in all the models.



 §.§ Units

We use M_d as the unit of mass and R_d as the unit of length throughout of this paper. In other words we use units such that R_d=M_d=1. Also we scale the Newtonian gravitational constant as G=1. In this case, the unit of velocity is V_0=(G M_d/R_d)^1/2, and time unit is τ_0=R_d/V_0=(R_d^3/G M)^1/2. For example a suitable choice for our models can be R_d=2.6 kpc and τ_0=10 Myr, which yields M_d≃ 4 × 10^10 M_⊙ and V_0≃ 254 km/s. For all the models we compute the evolution for 500 τ_0. We use a time step Δ t= 0.01 τ_0 and the time step is increased by successive factors of 2 in three radial zones. The cubic spline softening length <cit.> ϵ is set to 0.16 R_d in all the models. Also, except in the model EPR which has a rigid component, the grid is recentered every 16 time step. This dynamic recentering helps to avoid numeric artifacts.


 §.§ Rotation curves

Using the above mentioned units, the rotation curves of our four models have been illustrated in Fig. <ref>. It is necessary to mention that for an exponential disk the rotation curve in MOG is computed for a given α and μ_0 first, i.e. model EM, and then equivalent Newtonian systems are built with a dark matter halo added, in order to obtain the same rotation curve. More specifically by changing the halo parameters, i.e. mass and the characteristic length, in both the rigid and live cases, we construct models EPR and EPL which their rotation curve is similar to the model EM. It should be emphasized that in order to compare the stellar bar evolution in MOG and Newtonian gravity, it is necessary to start with the same initial rotation curves. In this case, the dispersion velocities are also almost the same. One should note that in model EM we have set α=8.9 and μ_0=0.015 R_d^-1, which α is close to the observed value. However, we are constructing toy galactic models and naturally we do not have to use the observational values. On the other hand, we know that these parameters change from galaxy to galaxy. 



 §.§ Bar amplitude


As we mentioned before, the main aim of this paper is to study and compare the bar growth in Newtonian gravity and MOG. Therefore, the main quantity to be measured in our simulations is the bar amplitude. It is important to mention that the code determines the gravitational forces on the cylindrical polar grid using sectoral harmonics 0≤ m≤ 8, and surface harmonics 0≤ l≤ 4 on the spherical grid. For more details see the online manual available at <http://www.physics.rutgers.edu/ sellwood/manual.pdf>.  

The amplitude of the non-axisymmetric features can be measured by computing

    A_m(t)=∑_j μ_j e^imϕ_j(t)

where μ_j is the mass of the  jth particle particle and ϕ_j is its cylindrical polar angle at time t. One should note that this summation is over the disk particles only. So the bar amplitude, for which m=2, would be obtained by the ratio A_2/A_0. 







 






§ RESULTS




 §.§ Stabilizing effects of MOG

In the recent paper <cit.>, we have compared two galaxy models which are similar to EM and EN. In fact without adding a spherical halo, it has been shown that MOG has stabilizing effects. More specifically, increasing the MOG free parameters, supports the global stability of the exponential disk. Albeit this point has been checked also for Mestel disk. In fact, the procedure in <cit.> is similar to what presented in <cit.>, in which it has been proved that a rigid halo has stabilizing effects. The maximum number of particles in <cit.> is 2× 10^4. Therefore the results may suffer from numeric artifacts caused by employing small number of particles. Here we have increased the particle number to 2× 10^7 and used a code which employs a completely different technique to compute the gravitational field. In model EN, the disk rapidly extends to larger radii. Therefore we have to enlarge the grid in order to avoid particles to escape from the grid. In the following let us briefly show that the main results of <cit.> do not change by increasing the number of particles. Then we will discuss more realistic models by adding dark matter halos to the Newtonian models.

In the model EM we have chosen α=8.9 and μ_0=0.015 R_d^-1. The corresponding rotation curve is shown in Fig. <ref>. As it is expected the circular velocity at large distances from the center of the disk is higher than the Newtonian model EN. The evolution of bar amplitude has been shown in Fig. <ref>. In both top and bottom panels in this figure, the black dashed curve belongs to model EN and the red solid curve belongs to model EM. We have used a logarithmic scale so that the growth rate can be estimated from the slop of the curve in the period of linear growth. The small boxes inside each panel show the evolution of the bar amplitude at early times τ< 100. It is clear from the top panel that the exponential instability growth starts around τ≃ 25 in both models EN and EM. However, in the model EM the growth rate is smaller than the Newtonian case. This point is more clear from the bottom panel where we have employed larger number of particles. In both cases, after reaching a maximum the bar decays slowly at large times. The slop of this decrease is smaller in the bottom panel. However the qualitative behavior is almost the same. It is clear from the top panel that the final magnitude of the bar in EM is smaller than EN. In the bottom panel this will happen at larger times. This fact has also been reported in <cit.>. 

The bar pattern speed has been shown in Fig. <ref>. Similar to the results presented in <cit.>, the pattern speed in EM is larger than EN model. 



 



 §.§ Bar growth


As mentioned before it is necessary to compare MOG with Newtonian models which include dark matter halo. In this section we compare model EM with two component models EPR and EPL. EPR has a rigid Plummer halo and EPL includes a live Plummer halo. In fact we have chosen the halos physical properties such that the rotation curve coincides with that of the model EM. As it is clear from the Table <ref>, the mass and the characteristic length scale in the rigid halo is different from the live halo. In the case of live halo, if we use the same parameters used in the rigid halo, the rotation curves will not match. In fact in model EPL the halo is compressed by the growth of the initial disk and leads to a different rotation curve than that of model EPR. For example see models "GR" and "GL" in <cit.>.
 

It is clear from Fig. <ref> that the bar grows rapidly in model ERL. On the other hand the growth rate in model EPR is substantially lower than both EM and EPL models. Also there are rapid oscillations in the bra amplitude of the EPR model. However it does not reveal a physical process. In fact for models with rigid components we do not use moving grids. Therefore the center of mass of the system rotates around the grid center and gives rise to these unreal oscillations in the bar growth.

This is a well known fact that disks in live halos form bars readily than in the rigid halos <cit.>. It has been shown in <cit.> that the increased growth rate of the bar instability in systems with live halo results from angular momentum exchange between the halo and the disk. We see from the top panel of Fig. <ref> that the exponential growth in both models EPL and EM starts around τ≃ 25 and reaches its maximum at τ≃ 100. However the instability growth rate is smaller in the model EM. More specifically, we found that the exponential growth rate, i.e. e^ωτ, in EPL and EM models are ω≃ 0.067 and ω≃ 0.054 respectively. 

It is important to mention that there is a significant difference between the bar growth in MOG and MOND. It is shown in <cit.> that the growth rate in MOND is substantially higher than the dark matter halo case. It is not the case in MOG, and MOG postpones the occurrence of the bar instability. Albeit although the halo density in <cit.> is as ours, the disk density in <cit.> is different and given by a Miyamoto-Nagai disk. Therefore, for a more careful comparison between MOG and MOND it is necessary to construct equivalent galactic models with the same surface densities. We leave this  work as a future study.

On the other hand, it is interesting to mention that the maximum value of the bar magnitude in model EM is smaller than in model EPL. This means that MOG leads to weaker bars during the bar instability. This behavior is completely in contrary with MOND which leads to stronger bars during the instability. It is clear in Fig. <ref> that the final bar magnitude in MOG is smaller than EPL. This is also the case in MOND. In other words, although the maximum magnitude of the bar is larger in MOND compared with the dark matter halo model, its final value is smaller <cit.>. 

The projected positions of the particles on the x-y, x-z and y-z planes for EM and EPL models are shown in Fig. <ref>. The left two columns belong to model EM and two right columns belong to model EPR. In both models, a strong two-arm spiral develops at τ≃ 100 and rapidly disappear and give place to a bar. The bar in model EPL starts to grow slowly. However, in EM model the bar magnitude decreases and at the same time oscillates.

This oscillation can be related to the existence of some beating between different modes. This oscillatory behavior in model EM can also be seen in the pattern speed evolution. In fact one should note A_2 is not only a bar strength, but also a spiral strength. In principle, these waves pattern speeds can be slightly different, and this can give rise to the above mentioned oscillations.

In Fig. <ref> we have plotted the relative overdensities associated to m=2 mode at different times. The left four columns belong to model EPL, and the right four columns to model EM. In the model EPL, after τ≃ 250 the bar mode retains its dominance. The twofold symmetry is clear until the end of simulations. On the other hand, in the EM model the evolution is completely different. As it is seen in the right panel of Fig. <ref>, the twofold symmetry fades and appears frequently. This is completely in agreement with the period of oscillations observed in Fig. <ref>.


 §.§ Power spectra

The clearest way in support of existence of beating between different modes in model EM, is to find the power spectrum. As we already mentioned, we compute the azimuthal Fourier coefficients of the mass distribution at each grid radius and save it at regular time intervals. In order to find the power in terms of frequency and radius for each sectoral harmonic m, we compute the Fourier transformation in time of the coefficients.

In Fig. <ref> we have shown the power spectra for three sectoral harmonics m=2, 3 and 4 for model EM with 2× 10^7 particles. The top row is for the first part of the evolution 22<τ<260, and the bottom row is for the second part of the evolution. In each panel, the solid (red) curve marks m Ω_c and the dashed (blue) curves m Ω_c ±κ. Where Ω_c(R) is the angular frequency for the circular motion and κ(R) is the Linblad epicycle frequency. Since each horizontal ridge indicates a coherent density wave, it is clear from the bottom left panel that there are two m=2 waves in the second half of the simulation with pattern speeds mΩ_1≃ 0.52 and mΩ_2≃ 0.3. In this case the beat angular frequency is Ω_b=m(Ω_1-Ω_2)= 0.22. Consequently the beat period is τ_b≃ 28.6, which is completely in agreement with the period of oscillations in the bar magnitude shown in the bottom panel of Fig. <ref>. Therefore the power spectra confirms the existence of beating between modes.

By plotting the magnitude of m=3 and m=4 modes, it turns out that they also have oscillatory behvior. Silimilar to the two-fold symmetric mode m=2, it is clear form the bottom row of Fig. <ref> that beating occurs also for m=3 and m=4 modes. For these modes there are at least four different beating waves with different frequencies. 
 



 §.§ Vertical structure

Particles can resonate both in the plane of disk and perpendicular to it. In this case the vertical oscillations of the particles can be amplified and a peanut shape forms <cit.>. One may infer from the top panel of Fig. <ref> that the bar strength drops substantially in both models EPL and EM at τ≃ 150. On the other hand, at this time the vertical thickness of the disk grows significantly. To see more clearly, we have plotted the root mean square (rms) thickness, i.e. √(⟨ z^2⟩), in Fig. <ref> for all the models for the same instants presented in Fig. <ref>. Therefore the drop in the bar amplitude coincides with the peanut resonance in both models.  The vertical resonances thicken the disk and consequently stabilize the disk against bar formation. In this case, it is natural to expect that the bar weakens. 
 On the other hand, in the EPR model there is no drop in the evolution of the bar strength, and consequently the thickness does not grow and the peanut shape does not appear. This confirms that the drop in the bar magnitude is due to peanut occurrence. 
 
It is seen in Fig. <ref> that the position of the peanut lobe in both models EM and EPL is in the interval 1.5  R_d<R<2  R_d and does not shift significantly during the disk evolution. Therefore let us measure the rms thickness at R=1.5  R_d with respect to time. The result has been shown in Fig. <ref>. The bar instability heats the disk and thickens it continuously. The steps correspond to the peanut occurrence when the particles leave suddenly the plane of the disk. The buckling occurs around τ≃ 80 in model EPL, and τ≃ 100 in the EM model. This is consistent with the above mentioned claim that weakening of the bar coincides with the appearance of the peanut.

 
In both models buckling occurs in a short time interval compared to the simulation time. In model EPL, during the buckling, the rms thickness increases from 0.05  R_d to 0.3  R_d. On the other hand in the model EM, the increase rate in slower and the rms thickness increases from 0.05  R_d to 0.2  R_d. Therefore this simulation shows that the dark matter stellar model produce stronger peanuts compared with the model in MOG. This is also the case in stellar models in MOND <cit.>.  


















 §.§ Pattern speed Ω_p

It is instructive to calculate the pattern speed Ω_p(t). The result has been illustrated in Fig. <ref>. There are 2× 10^6 particles in the each live components of the top panel. On the other hand, for comparison, there are 3× 10^6 particles in the bottom panel. As reported in <cit.> the pattern speed in model EM is larger than the EN model. On the other hand, Ω_p(t) oscillates in the model EM. This is not the case in other models. 

There are some important features in Fig. <ref>. The pattern speed starts with an almost same value for all models. However, although all models begin with same initial conditions, the pattern speed's magnitude and evolution is different when the disks evolve. More specifically Ω_p(t) oscillates in EM but its mean magnitude is almost constant. However in the model EPL, when the bar is formed at τ≃ 100, the pattern speed  starts to decrease at τ≃ 125. In fact the particles of the halo slow down the bar by dynamical friction. It is clear that there is no such a decrease in the model EPR where the halo is rigid. So, one may ensure that the substantial decrease in the model EPL is due to the dynamical friction induced by the live halo particles. As expected dynamical friction does not appear in the MOG model, i.e. in model EM. This is the case also in MOND, see <cit.> for more details. Absence of dynamical friction in modified gravity can be a crucial point to distinguish between dark matter halo and modified gravity from observational point of view. For example inefficiency of the dynamical friction can substantially reduce the merger frequency of the galaxies. As a consequence, bulge formation and its morphology and kinematics, in principle, can be different from the standard view <cit.>.


Let us put emphasize on an important observational fact which is relevant to the pattern speed of disk galaxies. Observations show that independent of the Hubble type of disk galaxies, bars have been formed and then evolve with time as fast rotators, see <cit.> and references therein. On the other hand these observations are not in complete agreement with the numerical simulations. For example, a recent ΛCDM cosmological hydrodynamical simulation by <cit.> verifies that bars slow down as they grow. Therefore existence of fast bars can be considered as a challenge to ΛCDM models <cit.>. 

In order to compare the pattern speeds in models EM and EPL, let us use the parameter ℛ defined as

    ℛ=D_L/a_B

 where D_L is the corotation radius and a_B is the bar semi-major axis. Bars are fast if ℛ≃ 1, and is slow if ℛ≫ 1. Finding the corotation radius in the simulation is not hard in the sense that we have the bar pattern speed and the angular frequency Ω(R) with respect to time. Therefore one may simply estimate the corotation radius. On the other hand in order to estimate the bar semi-major axis at a given time τ, we assume a rectangle with width Δ L and length 5 around the line y = tanϕ(τ) x, where ϕ(τ) is the angular displacement of the bar. Then, we divide this rectangle to small
elements, and choose the bar length to be the length at which the density of the particles is less than b per cent of the central element. We vary Δ L between 0.2 and 1, and choose two values for b, i.e. 10 and 20. Naturally, these choices lead to an error bar in the final value of ℛ at each time. 
 
 
 The result has been shown in Fig. <ref> for both models EM and EPL. In this figure we have computed ℛ for models with 2× 10^7 particles in each live component, at four different times. In fact we realize that the corotation radius in model EPL increases with time. More specifically, at τ=500 the corotation radius is D_L≃ 3.6, while at τ=800 it reaches to D_L≃ 5.3. On the other hand D_L stays constant in model EM around 2.7. Furthermore the bar length for both models is almost constant for τ> 500. Therefore, we expect that ℛ increases with time in model EPL, and stays constant in model EM. This behavior can be seen in Fig. <ref>. It is important mentioning that the bar is faster in model EM compared to the dark matter halo model. In fact, in reality most observed bars appear to lie in the range 0.9≲ℛ≲ 1.4, and consequently observed bars are fast. From this perspective, it is clear from Fig. <ref> that the model EM is closer to the observed interval, while both models lie outside the interval.  
  
 
 
However more careful studies are needed to compare MOG with the relevant observations. More specifically, cosmological hydrodynamical simulations would be extremely helpful to perform a statistical study for comparing MOG and observations. To the best of our knowledge, cosmological simulations in the context of alternative theories. which dispense with the need for dark matter particles, have not been developed.  



§ DISCUSSION AND CONCLUSION


In this paper, we have presented N-body simulations of isolated galaxies in both approaches, of Newtonian gravity with dark matter halo and MOG. The main aim is to compare the bar growth in the presence and absence of dark matter halo. More specifically, in the absence of dark matter halo, we used the modified gravitational force introduced in MOG. We construct three galactic models with almost the same initial conditions, i.e. EM, EPR and EPL. The model EM belongs to MOG and has no spherical component. On the other hand, model EPR is a galactic model with an exponential disk and a rigid Plummer spherical component. The model EPL is more realistic and its halo is a live Plummer sphere. In fact we have a fourth toy model named as EN. This model is an exponential disk in Newtonian gravity without dark matter halo. Comparing this model with EM, we recovered the main results already reported in <cit.> using low resolution simulations.



Our results show that the bar growth in model EM is slower than in model EPL, and faster than in EPR model. In fact, MOG behaves neither as live dark matter halo nor rigid dark matter halo. However when compared with EN model, it is clear that MOG has stabilizing effects and suppresses the bar instability for a while. Also we found that the final strength of the bar in MOG is smaller than both EPR and EPL models. Moreover the maximum magnitude of the bar in MOG is smaller than the Newtonian models. As we have already mentioned, this is in contrary with MOND where stellar simulations lead to stronger bars. However, although the bar keeps its maximum strength unchanged for a relatively long time in MOND, it keeps its maximum value for a short time in MOG.

We also found that the evolution of the pattern speed in these models is different. More specifically, the pattern speed in model EM, oscillates around an almost constant value and does not change substantially. However the pattern speed in model EPL, decreases significantly by dynamical friction induced by the dark matter particles. In other words, not only the morphology of the bar in our modified gravity theory is different from the standard cases, but also its dynamics is different. We explicitly showed that model EM leads to faster bars compared to the dark matter model.

Furthermore, the rms thickness of the disk in model EM at the end of the simulation is substantially smaller than model EPL. In other words the peanut formation is more effective in model EPL than in model EM. 


More investigations are still required for galactic dynamics in MOG. For example, it would be interesting to study the bar growth in the presence of a massive spherical bulge and galactic gas component. In fact the gas component have crucial effects on the system. For example, bars exert torque on the gas component, which inflows towards the center, and consequently the gas angular momentum is given to the bar. In this case, the bar is weakened, and finally the bar will be faster and shorter <cit.>. 

Also, the stabilizing effects of MOG should be compared with other well-known dark matter halo profiles. Furthermore, a careful study is required to compare bar evolution in MOG and MOND. In this case, the initial conditions, i.e. the disk density profile and the velocities, should be the same in both theories. Our results here show that, MOG has some features which are completely different from MOND. We leave these issues as subjects of study for future works.



To summarize the results, let us make a quick comparison between bar evolution in MOG and the dark matter halo models. The main differences are: 1) The bar growth rate is smaller in MOG. 2) The final magnitude of the bar is smaller in MOG. 3) The pattern speed is oscillatory in MOG. Also it does not decrease with time. This is the main difference between MOG and the dark matter models studied in this paper. 4) There is no dynamical friction experienced by the bar in MOG. 5) The maximum value of the bar strength is smaller in MOG. We reiterate that the bar magnitude in MOG stays in its maximum value for a relatively short time interval. In other words, MOG predicts that the strongly barred galaxies are not quiet frequent today. This fact is in agreement with the relevant observations which show that less than 20 % of barred galaxies are strongly barred <cit.>. 


As our final remark, we stress that in order to make a precise comparison between observation and MOG predictions for the abundance, size and
pattern speeds of bars, as well as for their evolution with redshift, it is necessary to perform high resolution cosmological simulations. For instance for a recent work in the standard cosmological model see <cit.>. However, there are no similar works in modified theories of gravity designed to dispense with the need for dark matter particles, even in MOND and its relativistic versions. Therefore, a large numerical effort is still required to be able to compare modified gravity predictions for galactic properties with relevant observations.
 


I would like to appreciate Francoise Combes for insightful and constructive comments. I am particularly grateful to Jerry Sellwood for guidances to run the GALAXY code. The main part of the calculations have been done on ARGO cluster in International Center for Theoretical Physics (ICTP).  I would like to thank ICTP for a very kind hospitality, during which the main part of this work has been done. This work is supported by Ferdowsi University of Mashhad under Grant No. 2/44972 (18/07/1396).


apj




