© G. Chardin et al. 2021

1. Introduction

The dark matter enigma found its first expression in the 1930s after the observation by Fritz Zwicky (Zwicky 1933) that peculiar velocities in the Coma cluster were far too large, by more than two orders of magnitude according to Zwicky’s analysis, to account for the bound behavior of the cluster components if only the visible mass was taken into account. Following this initial and remarkably prescient observation and analysis, a long dormant period followed where the dark matter question was mostly forgotten (however, see Kahn & Woltjer 1959). In the 1970s, Vera Rubin and collaborators (Rubin et al. 1980) and Albert Bosma (Bosma 1981), measuring galactic rotation curves, noticed that they generally had a rather flat behavior at large distances from the core, and their analysis contributed very significantly to the revival of the dark matter enigma. The accumulation of galactic rotation curves then led to the gradual realization (Bertone & Hooper 2018) that, quite generally, observed galaxy rotation curves are flat at large distances from the galaxy core, which is at odds with theoretical predictions based on the assumption of dominating mass related to luminous matter.

Two main lines of hypotheses were proposed as tentative solutions to this enigma. The first was to conjecture that there really exists a large part of the matter component of our Universe that is dark and interacts very weakly, apart from its gravitational interactions. This conjecture was the dominant hypothesis for a long time, under the implementation of weakly interacting massive particles (WIMPs), and later cold dark matter (CDM), after the findings that both massive neutrinos (Drukier & Stodolsky 1984) and supersymmetric particles (Goodman & Witten 1985) provided candidates that naturally suited the existing constraints, in terms of mass and interaction cross-sections, to solve the dark matter problem – the so-called WIMP miracle. This hypothesis seemed indeed almost necessary since the then mainstream Einstein-de-Sitter (EdS) universe featured a critical density: dark matter suggested an elegant way to fill the gap between the small baryonic component, with a density derived and constrained by nucleosynthesis to less than 5% of the critical density, and the critical density of the EdS universe.

But despite extensive experimental searches, both in direct detection experiments and at the Large Hadron Collider at CERN, no such dark matter candidates have been found. Meanwhile, the tensions between the age of an Einstein-de Sitter and the age of the oldest structures in the Universe (Bolte & Hogan 1995) grew to a point in the mid-1990s where a cosmological constant or some other repulsive component was considered compulsory, a few years before the discovery in 1998 of what is now called dark energy, through the SNIa flux measurements (Riess et al. 1998; Perlmutter et al. 1999).

The standard cosmological model, ΛCDM, although an impressive fit to diverse sets of data, does in fact present some tensions: probably the highest tension originates from the local measurements of the Hubble constant H₀ on the one hand, and the determinations of this same parameter deduced from the cosmic microwave background (CMB) and baryonic acoustic oscillations (BAO) on the other (Riess 2020; Di Valentino 2021). But even on the dark matter side, additional tensions exist: there is for example no evidence for the cusped density profiles predicted by the dark matter simulations (Flores & Primack 1994), and many galaxies apparently have no need for dark matter in their central regions. For some galaxies, there is even no apparent need for dark matter altogether, requiring to explain how the dominant dark matter component may have been almost entirely ejected from these galaxies (see e.g., Guo et al. 2020, and references therein). This scenario also suffers from the fact that its two main ingredients, dark energy and dark matter, supposed to represent 95% of its content, have not been identified despite extensive investigations. The extremely rare coplanarity of the satellites around Andromeda (Ibata et al. 2013), Centaurus A (Müller et al. 2018), and the Milky Way (Pawlowski & Kroupa 2020), at present the only three galaxies with well-known satellite distributions, is an indication of additional tension in the ΛCDM scenario. Finally, we note the too-big-to-fail problem (Boylan-Kolchin et al. 2011), possibly addressed by baryonic back-reaction (Governato et al. 2012; Teyssier et al. 2013; Wetzel et al. 2016).

The second line of explanation rests on a modification of the laws of gravitation, MOdified Newtonian Dynamics less (MOND; Milgrom 1983) being its most popular expression. The MOND hypothesis was initially phenomenological, and based on the observation that the dark matter problem seemed mostly confined to regions with values of the gravitational field ≲10⁻¹⁰ m s⁻². Indeed, MOND proposes that the law of gravitation deviates from its Newtonian expression in the following way: $a=\sqrt{a_0{a}_{\mathrm{N}}}$ where a is the acceleration, a₀ is the crossover acceleration of MOND, and a_N is the acceleration predicted by Newton’s law. With this single hypothesis, MOND predicted the Tully-Fisher relation (Tully & Fisher 1977) linking the mass to the rotation velocity for structures over a wide range of mass (at least four orders of magnitude). And although the MOND crossover acceleration, a₀, seems to differ when the analysis is done at galactic or at cluster scales, a reappraisal of the systematic errors in the cluster mass profiles (Ettori et al. 2019), and hybrid models with neutrinos as hot dark matter (Angus 2009; Haslbauer et al. 2020) could alleviate this tension. Looking for a field-theoretical expression, initial versions proposed for MOND, such as TeVeS (Bekenstein 2004), were recently ruled out by the multi-messenger observations of GW170817, but a slightly modified form passes this test (Skordis & Złośnik 2019). Possible links between MOND and cosmology include the entropic expression of gravity (Pazy 2013), emergent gravity (Verlinde 2017), mimetic gravity (Vagnozzi 2017) or some expressions of quantum gravity (Smolin 2017). The fascinating proximity between the values of a₀ and cH₀ is a suggestion of yet another link between cosmology and this crossover acceleration. On the other hand, we show in the following that Dirac-Milne may provide an explanation for this coincidence. A fundamental remark to which we come back is the fact that MOND’s law looks akin to the effect of gravitational polarization (Blanchet 2007; Blanchet & Le Tiec 2009).

In the present paper, we study the gravitational polarization predicted by the Dirac-Milne (D–M) cosmology (Benoit-Lévy & Chardin 2012; Chardin & Manfredi 2018), providing an explanation for this apparent modification of the Newtonian law of gravitation. For this purpose, in Sect. 2, we recall the main characteristics of the D–M universe. In Sect. 3, we study a simple idealized analytical model, showing that in the D–M cosmology rotation curves are indeed expected to be generically flat after a characteristic distance, for which we provide an approximate relation. We discuss the law obtained for the rotation velocity, notably in relation to the shell and Birkhoff theorems (Newton 1760; Birkhoff & Langer 1923), with a more detailed analysis in Appendix A. In Sect. 4, we present preliminary results of velocity distributions as a function of mass using a modified version of the RAMSES 3D simulation code incorporating the gravitational behavior of the D–M universe, as described in Manfredi et al. (2018) and Manfredi et al. (2020). In Sect. 5, we show that D–M follows a Faber-Jackson relation with a very small scatter and an exponent ≈3 over a range of more than three orders of magnitude in mass. In Sect. 6, we discuss in more detail the MOND-like behavior of D–M and in particular the value of the predicted acceleration parameter a₀. In the final section, we summarize our findings, and provide some perspectives and possible lines of development for future work.

2. The Dirac-Milne universe

The Dirac-Milne (D–M) universe, proposed recently (Benoit-Lévy & Chardin 2012; Chardin & Manfredi 2018), suggests a radically different paradigm for the cosmology of our universe. It features a symmetric matter-antimatter universe, where matter and antimatter effectively repel each other, but where antimatter also repels itself. As is well-known, when the usual expression of the weak equivalence principle is respected for matter and antimatter, such symmetric matter-antimatter cosmologies are excluded by the non observation of a diffuse gamma-ray flux (Omnès 1972; Cohen et al. 1998). On the other hand, the D–M cosmology is a gravitational implementation of the Dirac particle-hole sea, analogous to the electron-hole system in a semiconductor. Given the fact that a repulsive and enigmatic dark energy component represents ≈70% of the energy density in the ΛCDM model, it seems interesting to test more extensively this hypothesis of repulsion between matter and antimatter. The D–M cosmology is further motivated by the remark by Price (1993) that the usual expression of the weak equivalence principle, stating that all particles must follow the same trajectories given the same initial conditions in a gravitational field (universality of free fall), must necessarily be modified in the case where particles with negative mass are considered along with particles of positive mass. Indeed, as shown by Price, two new elements appear in such bound systems: (i) Gravitational polarization appears between particles of positive and negative mass whenever they are bound by non gravitational forces, and (ii) levitation is predicted for a symmetric ( + m, −m) system, a gross violation of the usual formulation of the weak equivalence principle.

A fundamental feature of the D–M universe is that its expansion factor varies linearly with time:

$$\begin{array}{c}\hfill a(t)\propto t\end{array}$$(1)

while in D–M there is neither dark matter, nor dark energy beyond its matter and antimatter components. Being a universe that appears gravitationally empty at large scales, the initial phases of the D–M universe have timescales differing rather dramatically from the ΛCDM universe: for example, the Quark-Gluon-Plasma transition lasts for about one day, instead of about 10 microsec in the Standard Model, while nucleosynthesis lasts about 35 years, compared to three minutes in the Standard Model, and recombination occurs at an age of about 14 million years, compared to the 380 000 years of the ΛCDM model (Benoit-Lévy & Chardin 2012; Chardin & Manfredi 2018).

Despite these tremendous differences in the initial timescales, the D–M universe, with only one adjustable parameter H₀, presents several elements of concordance: its age, equal to 1/H₀, is almost equal to the age of the ΛCDM universe for H₀ ≈ 70 km s⁻¹ Mpc⁻¹, while the H(z) dependence of cosmic chronometers is nicely reproduced in coasting universes (Melia & Maier 2013), along with primordial nucleosynthesis¹ and the SNIa luminosity distance (Sethi et al. 1999; Chodorowski 2005; Benoit-Lévy & Chardin 2012). Also, its non linear structure formation mechanism appears to reproduce the main features of the matter power spectrum starting from a single scale of matter-antimatter domains at decoupling (Manfredi et al. 2018, 2020). In addition, the D–M cosmology does not suffer from the horizon problem (Benoit-Lévy & Chardin 2012) and therefore does not require inflation.

3. Gravitational setup in the D–M universe

It is fundamental to note that although existing in nature, the Dirac particle-hole system has no Newtonian expression, even when the three Newtonian mass parameters (inertial, gravitational active, and gravitational passive) are used (Manfredi et al. 2018). On the other hand, as studied in Manfredi et al. (2018, 2020), the gravitational sector of matter and antimatter in the D–M universe can be expressed with two separate gravitational potentials, using the following coupled Poisson equations:

$$\begin{array}{cc}\hfill {\mathrm{\nabla }}^2{\varphi }_{+}& =4\pi G({\rho }_{+}-{\rho }_{-}),\hfill \end{array}$$(2)

$$\begin{array}{cc}\hfill {\mathrm{\nabla }}^2{\varphi }_{-}& =4\pi G(-{\rho }_{+}-{\rho }_{-}).\hfill \end{array}$$(2)

It should be noted that in these equations, although two potentials ϕ₊ and ϕ₋ are invoked, there is a single gravitational constant G, and not two independent constants for matter and antimatter. In particular, following Price (1993), we can derive the gravitational field for antimatter once the gravitational field for matter is known (and vice-versa). For particles with equal but opposite mass, the gravitational fields:

$$\begin{array}{c}\hfill {\mathit{g}}_{+}\equiv -\mathbf{\nabla }{\varphi }_{+}\end{array}$$(4)

and

$$\begin{array}{c}\hfill {\mathit{g}}_{-}\equiv -\mathbf{\nabla }{\varphi }_{-}\end{array}$$(5)

exerted on a particle and its antiparticle are opposite in the “Newtonian” regime (i.e., when the gravitational field created by matter is much larger than the contribution of antimatter, which can then be neglected), as the total force on the bound ( + m, −m) system is zero (Manfredi et al. 2018, 2020).

For further use, we also introduce the gravitational field acting on a particle sourced only by matter as:

$$\begin{array}{c}\hfill {\mathit{g}}_m\equiv -\mathbf{\nabla }({\varphi }_{+}-{\varphi }_{-})/2\end{array}$$(6)

and the corresponding field acting on a particle sourced only by antimatter:

$$\begin{array}{c}\hfill {\mathit{g}}_{\mathit{am}}\equiv -\mathbf{\nabla }({\varphi }_{+}+{\varphi }_{-})/2\end{array}$$(7)

3.1. General properties of matter, antimatter and depletion zone

In order to introduce the properties of rotation curves and virial velocity distributions in the D–M cosmology, let us first consider the following example of cluster configuration, represented on Fig. 1. This configuration was obtained in a 3D simulation using a modified version of the RAMSES code (Teyssier 2002), that we discuss more precisely in the next section. This figure represents a small cluster configuration in the D–M universe at redshift z ≈ 20 for self-gravitating particles, at this stage without dissipation. The full tomography of the matter, antimatter, and matter+antimatter configurations can be found as supplementary material at the following link².

Fig. 1. Slice of a matter-antimatter simulation of a small volume in the D–M universe using a modified version of the RAMSES simulation code. The thickness of the slice is 20% of the simulation box, of comoving size 1 h−1 cMpc. The condensed structures (represented in light blue) are all made of matter, while antimatter is spread in nearly homogeneous extended halos (represented in red) over half of the volume. The other half of the volume is occupied by depletion zones, surrounding the matter galaxies and clusters, isolating them from antimatter. These depletion zones are a consequence of the gravitational polarization expected between positive mass and “negative mass” objects. We note that antimatter clouds percolate, i.e., it is possible to roam at infinity without leaving the antimatter cloud, and the same property is respected for the depletion zones surrounding matter.

An animation of the formation of structures centered on a massive cluster of the simulation can be found at the following link³.

In Fig. 1, it can be seen that while matter has the usual clustering properties in planes, filaments and nodes (where globular clusters, galaxies and clusters of galaxies accumulate), and is concentrated in relatively small regions, antimatter has a completely different distribution, occupying nearly exactly half of the total volume, with a much more homogeneous density than matter. In addition, a new element is apparent: empty or very low-density depletion zones surround matter, isolating it from antimatter and occupying also about half of the total volume. We note that while it will be in most occasions possible to define rather well-defined isolated matter galaxies and clusters, this will not be the case for antimatter clouds, which are percolating. Strictly speaking, it is therefore incorrect to speak of “domains" for antimatter, since there is path continuity at all distances within the antimatter clouds. The same property of percolation is respected by the depletion zones, which are therefore organized more in tubes surrounding matter structures than in spheres.

Let us now show that both these antimatter clouds and the depletion zones between matter and antimatter lead to an additional surface gravity compared to the Newtonian expectation. For this, we first use a simple analytical model, represented on Fig. 2, where the D–M universe is represented as the periodic repetition of elementary cubic cells, in which half the volume is occupied by a galaxy (located at the center of the box) and the surrounding almost spherical depletion zone, while the antimatter cloud occupies the other half of the volume in the outer part of the box. Although this is clearly an oversimplification, since in this representation all structures are supposed to have the same mass and size, it will enable us to evidence the salient features of galactic rotation curves in the D–M universe.

Fig. 2. Schematic representation of a periodic structure of identical galaxies and their surroundings in the D–M universe. A galaxy core is represented by a condensed object at the center of a cubic cell, while a depletion zone, occupying about 50% of the cell volume, surrounds the galaxy. Antimatter is spread out beyond this depletion zone, on the outskirts of the cell, occupying also nearly 50% of the cell volume, although its volume appears smaller in the 2D projection.

More precisely, we postulate that each galactic cell is constituted by a cubic box of volume equal to the average volume per galaxy. We represent the content of this box as the sum of three contributions:

• (a)

  A galaxy concentrated in a very small volume at the center of the box, the density of the galaxy being much higher than the average matter density (typically a factor 200 within the virial radius);

• (b)

  A depletion zone, mostly empty of matter and antimatter, around the galaxy. In the galactic box, the volume of the depletion zone is about half of the total volume of the box when antimatter and matter in the box have equal and opposite “masses”, as supposed in D–M;

• (c)

  Surrounding this depletion zone, a region filled with antimatter with nearly constant density and a volume equal to half of the total volume of the box.

The whole space is then represented as a collection of such cubes in a periodic geometry.

We note that although the antimatter clouds occupy half of the volume of the box, the radius of the depletion zone extends almost to the confines of the box. In order to show this, let us derive, in our approximation of spherical symmetry, the extent of the depletion zone delimiting the zones of matter (mostly condensed) and antimatter (spread out almost uniformly). The constraints are the following:

• The total mass of the “galactic cell” is zero, i.e., m₊ = m₋ ≡ m

• The depletion zone is spherical and of radius r_d

• The linear size of the individual galactic cell is L, and its total volume is therefore L³.

• The antimatter, due to its internal repulsion, has a constant density ρ₋.

The edge of the depletion zone is defined by the condition that the total gravitational force on an antimatter particle becomes zero. In order to express the repulsive force exerted by the antimatter cloud on an antiparticle at the edge of the depletion zone, we use the fact that the force created by the antimatter cloud of uniform density would be zero were it not for the spherical depletion zone. On the other hand, creating the depletion zone by removing an homogeneous sphere of antimatter creates a force directed toward the center of the box, and its intensity can be calculated by Newton’s shell theorem stating that the situation is equivalent to that where the antimatter mass of the depletion volume is concentrated at the center. The spherical depletion zone will therefore create an attractive force as it if were concentrated at the origin and with a mass $\frac{4}{3}\pi {\rho }_{-}{r}_d^3$. Since we want the total gravitational force on an antimatter particle at the edge of the depletion zone to be zero, it means that this mass must be equal to m₊.

On the other hand, the remaining mass of antimatter, outside the depletion zone, is:

$$\begin{array}{c}\hfill m_{-}={\rho }_{-}(L^3-\frac{4}{3}\pi r_d^3)\end{array}$$(8)

By our first condition, it must be equal to the positive mass m₊ of the galaxy.

Overall, we have therefore the conditions relating ρ₋ and r_d:

$$\begin{array}{c}\hfill {\rho }_{-}(L^3-\frac{4}{3}\pi r_d^3)=\frac{4}{3}\pi {\rho }_{-}{r}_d^3=m_{+}=m_{-}\end{array}$$(9)

Therefore, in our approximation of spherical symmetry, the antimatter cloud occupies half the volume of a simulation cube, and the depletion zone extends until:

$$\begin{array}{c}\hfill r_{\mathrm{d}}=\sqrt[3]{\frac{3}{8\pi }}L\approx 0.492L\end{array}$$(10)

This means that, counterintuitively, the depletion zone almost reaches the confines of the average individual box around a galaxy. It is important to note that the situation that we have represented on Fig. 2 has not a spherical symmetry, notably concerning the antimatter cloud. Also, in the actual 3D equilibrium configuration, the depletion zone of a cell connects with that of the adjacent cells on its three axis.

3.2. Analytical approximation of the rotation velocity

We now estimate the rotation velocity for a matter test particle as a function of the distance to the galaxy in the preceding configuration. The galaxy of mass m obviously contributes, for a test body located at the distance r from the galaxy center, to an acceleration equal to Gm/r², as in the Newtonian case.

For a matter particle, the gravitational force due to the configuration of the antimatter cloud and the depletion zone surrounding the galactic core can be estimated as follows. We first consider the configuration where antimatter of “negative mass” (see Manfredi et al. 2018 for a precise definition), is uniformly distributed in the box. The gravitational force exerted by this uniform antimatter cloud is then zero everywhere, by symmetry. Then, we represent the depletion zone as the sum of a homogeneous sphere of positive density (equal and opposite to that of antimatter) covering the volume of the depletion zone, which will compensate the uniform background of the negative mass fluid.

For each individual cubic cell, the situation shown in Fig. 2 can therefore be expressed as the superposition of three cubes (Fig. 3):

Fig. 3. Our approximation of a galaxy represented on Fig. 2 can be represented as the sum of three cubes: a cube (a) with uniform “negative mass” density of total mass −2m, a cube (b) with only the galaxy of mass +m at the center, a cube (c) with a sphere occupying half the volume of the cube and of homogeneous positive density of matter, equal and opposite to that of the first cube, that together with (a) will create the empty depletion zone.

(a) A cube with “negative mass” equal to −2m, i.e., twice the mass of the galaxy, and uniform density (zero gravitational fields g₊ and g₋ everywhere);

(b) A cube with the mass m of the galaxy at the center, i.e., half of the total mass of the antimatter in the first box. This galactic point-like mass will produce the classical Newton force field g₊ = −Gmr/r³ (and g₋ = −g₊);

(c) A cube with the same positive mass +m as the galaxy, but uniformly distributed over the spherical depletion zone. This positive mass will compensate the negative mass −m of the first cube contained in the same sphere, thus creating the empty depletion zone, with a volume 50% of the box.

As is well known, the homogeneous matter sphere in the third cube will generate a harmonic restoring force, which is superimposed on the Newtonian force of the galaxy, supposed to be a point galaxy. We note that the second and third cubes have a total mass of 2m, compensating the “negative mass” −2m of antimatter in the first cube. Space is supposed to be covered by adjacent such cubes (see Fig. 2).

At a distance r from the center of the box (and the center of the galaxy), the gravitational field created by this homogeneous sphere of positive density will then be:

$$\begin{array}{c}\hfill {\mathit{g}}_{+}=-\frac{Gm\mathit{r}}{r_{\mathrm{d}}^3}\end{array}$$(11)

where m is the positive mass added to create the depletion zone and r is the radial distance vector to the galactic center. The total gravitational field g₊ exerted on a matter particle orbiting at a distance r from the galaxy center is therefore:

$$\begin{array}{c}\hfill {\mathit{g}}_{+}=-\frac{Gm\mathit{r}}{r^3}-\frac{Gm\mathit{r}}{r_{\mathrm{d}}^3}\end{array}$$(12)

while the relation between orbital velocity and this force is:

$$\begin{array}{c}\hfill \frac{v^2}{r}=\frac{\mathit{Gm}}{r^2}+\frac{\mathit{Gmr}}{r_{\mathrm{d}}^3}\end{array}$$(13)

The orbital velocity as a function of the distance r to the center of the galaxy can be simply deduced from the above expression, and reads as:

$$\begin{array}{c}\hfill v(r)=\sqrt{\frac{\mathit{Gm}}{r}+\frac{Gm{r}^2}{r_{\mathrm{d}}^3}}\end{array}$$(14)

Remarkably, this function has a derivative that vanishes at r = 2^−1/3r_d ≈ 0.79 r_d, and is almost constant for all r ≳ r_d/2. This function is represented in Fig. 4, where we also plot the expected orbital velocity assuming there were no depletion zones nor antimatter clouds on the outskirts. This figure shows that for such idealized galaxies, the D–M orbital velocity reaches a nearly constant plateau, even rising very slightly near the end of the depletion zone in our approximation of spherical symmetry. We note, however, that for r > r_d, the spherical approximation breaks down as we have seen in Sect. 3.3. Additionally, beyond the depletion zone, the antimatter cloud is present.

Fig. 4. Rotation velocity predicted in the D–M universe for a point mass located at r = 0 and surrounded by a depletion zone and the same mass of antimatter (with opposite “sign”). Instead of the expected Newtonian velocity (see Eq. (11)), decreasing to zero, D–M predicts that the rotation velocity for a condensed object of positive mass will be nearly constant (see Eq. (14)) above the characteristic distance approximately equal to half the radius rd of the depletion zone, or ≈3 virial radii.

Since the depletion zones occupy about 50% of the total volume and since the baryonic density in a D–M universe is defined by the baryon/photon ratio η ≈ 8 × 10⁻⁹ (Sethi et al. 1999; Benoit-Lévy & Chardin 2012), an observer assuming that the ΛCDM cosmology is valid will interpret the depletion zone as a halo with a density equal to the matter density of the D–M universe, that is:

$$\begin{array}{c}\hfill \frac{{\rho }_m}{{\rho }_{\mathrm{c}}}=\frac{{\eta }_{\mathrm{D}-\mathrm{M}}}{{\eta }_{\mathrm{\Lambda }\mathrm{CDM}}}\times {\mathrm{\Omega }}_{\mathrm{baryon}}^{\mathrm{\Lambda }\mathrm{CDM}}\approx \frac{8\times {10}^{-9}}{6\times {10}^{-10}}\times 0.049\approx 0.65\end{array}$$(15)

where η = n_B/n_γ is the number of baryons per photon in the corresponding cosmology, and ρ_c is the critical density for the ΛCDM universe at the present epoch.

Now, the precise characterization of the virial radius depends on the cosmological model used and notably on the inflows between the redshift of initial collapse of a structure and the present epoch. The overdensity Δ_c at the virial radius for ΛCDM is usually approximated as Δ_c ≈ 100ρ_c (Shull 2014), so that a nearly constant rotation velocity will be observed for distances:

$$\begin{array}{c}\hfill r\gtrsim \frac{r_{\mathrm{d}}}{2}\approx \frac{1}{2}\sqrt[3]{\frac{{\mathrm{\Delta }}_{\mathrm{c}}}{{\rho }_{\mathrm{m}}}}{r}_{\mathrm{v}}\approx \frac{1}{2}\sqrt[3]{150}{r}_{\mathrm{v}}\approx 2.6r_{\mathrm{v}}\end{array}$$(16)

where r_v is the virial radius. Therefore, the flat rotation curves, instead of being due to an invisible halo of slowly decreasing density, are in fact due in D–M to the asymmetric configuration of matter and antimatter that we have described above.

Indeed, the D–M universe tells us that in order to derive the galactic rotation curves and the virial velocities in a cluster of galaxies, it is necessary to take into account not only the matter present between us and the center of the galaxy or the cluster, i.e., matter “below our feet”, but also (anti)matter “above our heads”, which may be interpreted as an external field effect (EFE; Chae et al. 2020). This is indeed very counterintuitive, and may even be considered in contradiction with the Newtonian shell theorem (Newton 1760), stating that inside a spherical shell, spacetime must be Minkowskian, i.e., the gravitational field must be zero. We further discuss this apparent contradiction in Appendix A.

The approximation of spherical symmetry allowed us to show that flat rotation curves are expected within a large fraction of the volume of the depletion zone. Using a numerical simulation, we now explore in more precision the edge effects created by the antimatter clouds, which are in fact asymmetric.

3.3. Numerical study of the depletion zone

In this subsection, we use our modified version of RAMSES (see Sect. 4 for details) to study numerically the spherical approximation of the previous subsection. We consider a galaxy shell configuration following Sect. 3.1, as described in Fig. 5.

Fig. 5. Configuration studied in this simulation. A matter galaxy, located at the center of the box, is surrounded by a depletion zone, extending over half the volume of the box, and by an antimatter cloud (represented in magenta), extending to the outskirts of the box over the other half of the volume.

More precisely, we considered an initial configuration where a central galaxy core, confined in a radius equal to ≈5% of the simulation box length, presents an isothermal density and velocity profile. The same mass of antimatter is uniformly distributed at r > r_d ≈ 0.5L, (see Eq. (10)), but its initial velocity is zero. The boundary conditions for the simulation box are periodic, which means that the actual configuration simulated is that of Fig. 2, repeated out to infinity.

Using RAMSES on a 128³ grid, we checked that this initial configuration does not evolve noticeably. The gravitational field g_am created by antimatter is represented in Fig. 6. Due to the non spherical distribution, the antimatter field is seen to be stronger on the diagonals of the cube for a plane with normal parallel to one of the axes of the cube, for example with a normal vector (1, 0, 0). In the bottom panel of the same figure, we have represented the antimatter field for a plane along one of the diagonals of the cube, for example with a normal vector $\frac{(1,1,1)}{\sqrt{3}}$. For such planes, the confining field created by antimatter is more symmetrical, in a ring configuration, and stronger than the field present in the plane parallel to one of the faces of the cube.

Fig. 6. Gravitational field gam created by antimatter on two planes defined by their normal vector (1,0,0) for the top panel and $\frac{(1,1,1)}{\sqrt{3}}$ for the bottom panel. The color scaling indicates the strength of the field. The confining gravitational field created by the antimatter cloud is more symmetric and has a higher intensity in the bottom panel, which represents the field in one of the diagonal planes of the cube.

In order to quantify this statement, we represent in Fig. 7 the total gravitational field experienced by a matter particle, created by both matter and antimatter, as a function of the gravitational field created by matter alone. These gravitational fields are calculated on a 128³ Cartesian grid (about 2 million cells), and the quantities estimated on the matter region, as defined by the distance to the galactic center, which is the only cut applied for the selection. A Newtonian regime can be observed in the high field region of the figure, with a gradual transition from the Newtonian regime to a MOND-like regime beginning at an acceleration of ≈1.0 × 10⁻⁹ m s⁻². Fitting an a₀ parameter according to MOND interpolating formulas listed in Famaey & McGaugh (2012) results in values of the a₀ parameter ranging between ≈1.5 and ≈5 × 10⁻¹⁰ m s⁻², depending on the parametrization chosen. Using as an example the parametrization of Lelli et al. (2017), we present in the bottom panel of Fig. 7 the ratio of the gravitational field acting on a matter particle, g₊, to the gravitational field created by matter alone, i.e., the Newtonian expectation, g_m. Superimposed on this 2D distribution, we have plotted the analytical expression of Lelli et al. (2017) for values of the a₀ parameter between 1.5 and 3.0 × 10⁻¹⁰ m s⁻². We have also indicated the upper and lower limits of 0.11 dex (1-σ error) reported by Lelli et al. (2017) in their analysis.

Fig. 7. Top panel: total gravitational field acting on a matter particle as a function of the gravitational field created by matter alone. A Newtonian regime can be observed for the high field region of the figure, with a gradual transition from the Newtonian regime to a MOND-like regime at an acceleration of ≈1.0 × 10−9 m s−2. Middle panel: total gravitational field acting on a matter particle as a function of the gravitational field created by antimatter alone. The antimatter field is rather peaked, with an average value of ≈2 × 10−10 m s−2. Bottom panel: ratio between the total gravitational field acting on a matter particle and the Newtonian field (created by matter only), as a function of the gravitational field created by matter only. Almost exactly unity in the Newtonian regime (at high values of the field), this ratio gradually increases and reaches a factor of ≈2 at the end of the depletion zone. The MOND interpolating function used in Lelli et al. (2019) for a value of the a0 parameter 1.85 × 10−10 m s−2 has been superimposed on the simulation data, while the two adjacent curves represent the 1-σ error of 0.11 dex found by these authors.

It should be noted, however, that the proximity with the value of the MOND parameter equal to ≈1.2 × 10⁻¹⁰ m s⁻² is coincidental, as other values of the mass of the central structure in our simulation would have changed the value of the fitted parameter a₀. We further discuss this crucial aspect in Sect. 6, which is dedicated to the simulations of “clusters” of mass between 2 × 10¹⁰ M_⊙ and 1 × 10¹⁶ M_⊙, and where we study the evolution and value of the a₀ parameter predicted by D–M as a function of redshift and total mass.

3.4. Self-consistent model of the depletion zone

In the present section, in order to confirm the qualitative and numerical analysis of the depletion zone presented in the two preceding subsections, we extend the analysis of Manfredi et al. (2020) to incorporate the notion of depletion zone and inhomogeneity of the antimatter cloud. Using the Poisson equations (2) and (3), we consider a spherically symmetric geometry, where all quantities depend only on the radius r and the Laplacian operator is defined as ${\mathrm{\nabla }}^2\varphi =\frac{1}{r^2}{\partial }_r(r^2{\partial }_r\varphi )$, where ∂_r is a radial derivative. We want to solve these Poisson equations for a typical situation where a positive high-density mass (“galaxy”) is located in a small spatial region located near r = 0 and is described by the distribution ρ₊(r). In contrast, negative masses are supposed to be thermalized (at low temperature T) and described by a Boltzmann distribution:

$$\begin{array}{c}\hfill {\rho }_{-}(r)={\rho }_{0}exp\left(\frac{-m{\varphi }_{-}+\mathrm{\mu }}{k_{\mathrm{B}}T}\right),\end{array}$$(17)

where k_B is Boltzmann’s constant, ρ₀ is a reference density, and μ is a chemical potential that will be chosen so that positive and negative masses are present in equal amounts in the computational box, i.e., 4π∫ρ₊r²dr = 4π∫ρ₋r²dr. We note that the mass m appearing in Eq. (17) is the passive gravitational mass, which in the D–M model is always positive (Manfredi et al. 2018). Then, Eq. (3) becomes a non linear Poisson equation for ϕ₋(r), which can be solved self-consistently, and the result substituted into Eq. (2) in order to obtain ϕ₊(r).

However, as discussed at the end of Sect. 3, the external field effects require a careful treatment of the boundary conditions. Indeed, as they are written, the Poisson Eqs. (2)–(3) describe an isolated system in otherwise empty space. But this is not the case for a cosmological setting for which the positive masses are localized in small regions, and the negative masses spread out almost uniformly across about 50% of the available space.

We therefore solve Eqs. (2)–(3) in a spherical region between r = 0 and r = R, using the following boundary conditions:

$$\begin{array}{c}\hfill \varphi {\prime }_{+}(R)=2\frac{4\pi G{\int }_0^R{\rho }_{+}{r}^{2}dr}{R^2},\varphi {\prime }_{-}(R)=0\end{array}$$(18)

where the apex stands for differentiation with respect to r.

For comparison, the Newtonian case is computed by solving: ∇²ϕ_Newt = 4πGρ₊, with boundary condition: ϕ_Newt(R) = 0.

In this subsection, we use units in which 4πG = 1, and take R = 20 and temperature T = 0.001. The depletion zone is clearly seen in the simulation results shown in Fig. 8, and extends up to r_d ≈ 15.8. Thus, the ratio of the total volume to the volume of the depletion zone is approximately $R^3/r_{\text{d}}^3=2.04\approx 2$, as expected. On the outer border of the depletion zone, the gravitational field (defined as: g = −∂_rϕ) for matter is twice the Newtonian value, while the antimatter field vanishes within the entire depletion zone (see Fig. 9), justifying our hypothesis that the antimatter cloud is both very cold and homogeneous. The rotation speeds are defined as: $\mathit{v}(r)=\sqrt{r|\varphi \prime (r)|}$. For the D–M model, the rotation curve flattens between approximately r_d/2 and r_d, also in accordance with our model (see Fig. 10).

Fig. 8. Matter density ρ+(r) (black dashed line), antimatter density ρ−(r) (blue line), and pseudo-dark-matter halo density: ${\rho }_{\mathrm{halo}}=2\overline{\rho }-{\rho }_{-}$ (red line).

Fig. 9. Gravitational fields: g− (blue line), −g+ (black dashed line), and −gNewt (red line).

Fig. 10. Rotation speeds for the D–M model (black line) and for standard Newtonian gravity (red line). The decrease of the rotation velocity beyond r/rd = 1 is non physical, as the test particle would enter the antimatter cloud.

Finally, we point out that the Poisson equations (2)–(3), with the boundary condition (18), are equivalent to adding a constant density $2\overline{\rho }$ on the right-hand side of both equations (where $\overline{\rho }\equiv ⟨{\rho }_{+}⟩=⟨{\rho }_{-}⟩$ is the average matter or antimatter density) and using Dirichlet boundary conditions ϕ_±(R) = 0, i.e., fixing the value of the potentials on the sphere (r = R). The conditions on the gradients (the fields), i.e., Eq. (18), will be automatically satisfied because of Gauss’s theorem. Hence, the potential ϕ₊ acting on matter results from two sources: the central “galaxy" ρ₊ and the additional distribution ${\rho }_{\mathrm{halo}}\equiv 2\overline{\rho }-{\rho }_{-}$. The latter is analogous to the dark matter halo postulated in the standard CDM theory. This halo distribution is also represented in Fig. 8.

4. Simulation results

In order to validate the analytical approximations studied in Sect. 2, we modified the Adaptive Mesh Refinement code RAMSES (Teyssier 2002) in order to implement the gravitational behavior of matter and antimatter present in D–M, following Eqs. (2)–(3). We note that in most cosmological simulations, the average density is first calculated and subtracted from the local density in order to calculate the evolution of the scale factor a(t), so that the Poisson equation is usually written ∇²ϕ = 4πGa²δρ. In D–M, this is directly the case, since $\overline{\rho }=0$ in this matter/antimatter symmetric universe, so that $\delta \rho \equiv \rho (t,\mathit{x})-\overline{\rho }(t)=\rho (t,\mathit{x})$. The new element introduced in this D–M cosmological version of RAMSES is an extra set of particles with negative mass, which in turn introduces the new aspect of a depletion zone.

The acceleration of these negative mass particles is given by the gravitational potential of Eq. (3) (respectively, Eq. (2) for positive mass particles). The mass density of both particle species are projected separately onto the mesh with a cloud-in-cell interpolation, so that the mesh contains both the positive and negative mass density of the corresponding positive and negative mass particles. Once these two separate mass densities on the mesh are obtained, the positive and negative gravitational potentials are derived using the standard conjugate gradient algorithm of RAMSES. The plus and minus accelerations are computed with a simple finite difference and attributed to the corresponding set of positive and negative mass particles with the cloud-in-cell interpolation of the accelerations computed on the mesh. All particles share a common level-based timestep obtained from the smallest Courant condition from particle velocities and free fall time. The cosmological time is a linear function of the scale factor a(t), as in the Milne geometry (Milne 1933).

Numerical Setup. The first simulation that we present involves a small simulation volume, of dimensions 1 h⁻¹ cMpc, discretized on a grid of 256³ cells. We allowed for up to two refinement levels, therefore leading to an effective resolution of Δx ≈ 1 h⁻¹ ckpc. According to Benoit-Lévy & Chardin (2012), Manfredi et al. (2020), the initial size of matter and antimatter domains is on the order of 100 pc at z ≈ 1080 (i.e., 100 ckpc at the present epoch). Thus the resolution will capture these domains and their contrast at the initial times. We ran the simulation from z = 1080 to z = 0. At a redshift of z = 20, 7769 halos were identified using the AdaptaHOP algorithm (Aubert et al. 2004), the development of structures in D–M occurring mostly before z = 3. This allows to sample structure characteristics such as size and velocity dispersion over a range of about three orders of magnitude in mass (10⁵ − 2 × 10⁸ M_⊙). The total mass of matter present in this simulation is ≈3.3 × 10¹⁰ M_⊙ (and the same mass of antimatter), which gives an individual equal mass for each particle of ≈3.9 × 10³ M_⊙. For the initial conditions, we used a Gaussian velocity distribution with dispersion: ${\mathit{v}}_I=\sqrt{\frac{\mathit{Gm}}{r}}\approx 2\mathrm{km}{\mathrm{s}}^{-1}=\mathcal{O}(1)$ km s⁻¹.

We checked that, even when the particles are initially distributed on a uniform grid with v = 0, the qualitative behavior of our results is not significantly modified. This is due to the fact that the initial contrast in density in D–M is already of order unity immediately after the CMB transition, leading to a very efficient virialization of the first structures within typically one Hubble time at that epoch, a few million years.

In the simulation, we use the value H₀ = 70 km s⁻¹ Mpc⁻¹ for the Hubble parameter at the present time. We note that, unlike standard large-scale structure simulations, we did not impose a given initial power spectrum but rather a specific matter/antimatter pattern. In this first simulation, the matter-antimatter pattern for the initial condition was generated using an Ising code (with two states per spin, used here to represent matter and antimatter), on a 256³ Cartesian grid evolved for a few time steps starting from a random distribution, and using a temperature well below the second order transition. This procedure was used in order to create a fine-grained distribution where matter and antimatter both percolate in analogy with an emulsion. The characteristic size of the emulsion (≈100 pc at z = 1080) results from the annihilation and polarization of the matter and antimatter emulsion at much higher temperatures (Benoit-Lévy & Chardin 2012), and we do not describe this phase of the evolution here. This aspect will be treated elsewhere.

Results. Analogous to observations of electrons and holes in semiconductors (Tsidil’kovskii 1975), this small-scale simulation exhibits a nearly empty depletion zone between the condensed clumps of matter and the extended clouds of antimatter. Figure 1 illustrates this geometric distribution, asymmetric between matter and antimatter. We also refer the reader to the supplementary movie available at the following weblink⁴ showing that unlike in ΛCDM, structure formation in D–M starts very shortly after the CMB transition, as was already predicted by our earlier 1D simulations (Manfredi et al. 2018, 2020). The configuration represented in Fig. 1 is due to the fact that matter and antimatter repel each other, leading to a zone where almost no particles or antiparticles are present. As antimatter particles repel each other, antimatter spreads as much as it can without going close to matter, by which it is also repelled. Matter clusters and forms halos as in the conventional gravitational scenario. Indeed, by looking only at the matter zones, it is difficult at first glance to distinguish, with its planes, filaments and nodes, the matter configuration in D–M compared to ΛCDM. Obviously, this qualitative statement deserves a more detailed study.

On the other hand, we can look at the common features of matter halos present in the 3D simulation. In Fig. 11, we show the logarithmic density profile of matter and antimatter in 25 halos, with individual masses on the order of 10⁷ M_⊙ in our simulation. These density profiles share common features such as a bulge region where the matter density is steeply decreasing with the distance from the center of the halo, while the density of antimatter in the same region is zero. After a few virial radii, the antimatter density becomes non zero and the enclosed mass of antimatter rapidly increases and becomes of the same order as the mass of matter present in the halo, with nearly equal average densities. These extended antimatter clouds, present beyond the depletion zone, create an approximately harmonic restoring gravitational field, which adds its contribution to that of matter present in the halo, and lead to nearly flat rotation curves. This behavior is generic to most halos present in our simulation, although some isolated halos can present a somewhat different behavior when their vicinity includes more massive structures. Depending on their antimatter environment, they might then appear as being completely devoid of dark matter when their close environment does not contain antimatter, with accelerations described by the expected Newtonian behavior or, on the contrary, as being submitted to strongly confining fields, created by the configuration of surrounding matter and antimatter clouds, as in the bottom panel of Fig. 6, and appear as almost entirely composed of dark matter.

Fig. 11. Logarithmic spherical density as a function of radius for both matter and antimatter for a set of 25 halos with mass ≈107 M⊙ in our RAMSES simulation. While the density of matter, represented in blue, is steeply decreasing from the halo center, the density of antimatter, represented in orange, is zero in this region. The logarithmic density is scaled to the average matter and antimatter density, equal in D–M at large scales, and converge (equal average densities for matter and antimatter) for almost all halos. It can be seen that the antimatter density becomes non zero after a few virial radii, and that the mass of antimatter becomes on the order of the halo matter mass after typically twice the radius at which antimatter density becomes non zero. The distribution for 25 consecutive halos has been represented in the figure in order to show that this behavior is quite generic in D–M, although some variations can be seen in the distribution of antimatter around the matter halos.

As discussed in Manfredi et al. (2020), the highly non linear structure formation in D–M appears to give at the present epoch the same order of magnitude for the matter power spectrum at the peak scale as its ΛCDM counterpart. However, at higher redshifts, the situation is quite different and bound structures appear much earlier in D–M than in ΛCDM. This may indeed alleviate the constraints met by ΛCDM in this regard (Kovács et al. 2019; Asencio et al. 2021; Lelli et al. 2021). It should be noted that in the evolution for redshifts between z = 100 and z = 20, the formation of structures in the upper range in mass may be underestimated as the size of the simulation box is very limited. Note also that in the present simulations, particles are dissipationless (behaving as dark matter particles), which may lead to a modification of the power spectrum at small scales.

A few additional remarks are in order: the average velocity of antimatter clouds is much smaller (typically by one order of magnitude, depending on redshift) than that of matter. In addition, these peculiar velocities for antimatter correspond to nearly global flows of cold antihydrogen (and antihelium) clouds, and the actual temperature of the antimatter clouds is much colder than that of matter. Also, as mentioned previously, antimatter clouds percolate with one another, which means that it is possible to travel continuously at infinity without leaving the antimatter region. Similarly, the depletion zones also exhibit percolation, meaning that matter structures are confined in tubes rather than in depletion spheres. In order to better apprehend the gravitational influence of these extended antimatter clouds, we now turn to the study of the Tully-Fisher and the Faber-Jackson relations.

5. Tully-Fisher and Faber-Jackson relations in the D–M universe

The Tully-Fisher relation (TFR; Tully & Fisher 1977) and Faber-Jackson relation (FJR; Faber & Jackson 1976) are phenomenological scaling relations between the velocities inside galaxies and clusters and the baryonic mass of the structure considered. The scatter around a relation of the type:

$$\begin{array}{c}\hfill m\propto v^{\alpha }\end{array}$$(19)

is surprisingly low, and seems to indicate that either there is no dark matter, or that there is a strict correlation between the dark matter and the baryonic matter component, which seems difficult to justify in the ΛCDM model. MOND, on the other hand, predicts such a strict correlation and an exponent α = 4 (Famaey & McGaugh 2012).

Let us now show that such power-law relations are a natural consequence of the D–M cosmology. We have seen that, due to the presence of the antimatter component on the outskirts of galaxies and clusters of galaxies, rotation curves are almost flat beyond a characteristic distance, which, to a good approximation, is about half the size of the depletion zone (see Fig. 2), and about 2.5r_v, where r_v is the virial radius as usually defined in the ΛCDM cosmology. At distances r > r_d/2 ≈ 2.5r_v from a condensed object, the D–M rotational velocity is predicted to be almost constant and equal to ${\mathit{v}}_{\mathrm{d}}=\sqrt{2Gm/r_{\mathrm{d}}}$.

The effective radius of the depletion zone⁵ is in turn related to the mass m of the structure considered, at least on average, and will be approximated here by the following simple equation between m and r_d: $m\approx (4/3)\pi r_{\text{d}}^3$ $\overline{\rho }$ where $\overline{\rho }$ is the average density of matter. Replacing this in the above expression for v_d, one finds that the asymptotic and nearly constant velocity is approximately proportional to r_d, and the following approximate relation holds: $m\propto v_{\text{d}}^3$.

In this very crude approximation, where we have assumed spherical symmetry, the exponent of the TFR preferred in D–M is ≈3, significantly smaller than the index 4 predicted by MOND (see e.g., McGaugh et al. 2000; Lelli et al. 2019). On the other hand, we note that while the analysis of Lelli et al. (2019) finds a range of exponents between 3.5 and 4, with a preferred value of 3.85 ± 0.09, the analysis of Ponomareva et al. (2017) finds values of the exponent close to 3.0, with slope errors between 0.11 and 0.31 (see e.g., Table 4 and Fig. 12 of Ponomareva et al. 2017). It can also be noted that the estimate of the exponent α in the TFR relation varies rather widely depending on the photometric band used to infer the mass, although the near-infrared band using Spitzer photometry (Lelli et al. 2019) may help to reduce this uncertainty.

In order to take into account more realistic situations than the idealized “spherical” galaxy studied above, we have first analyzed the FJR (Faber & Jackson 1976), using the RAMSES simulation described in the preceding section. Using the set of 7769 halos identified by the AdaptaHOP algorithm (Aubert et al. 2004), we calculated the velocity dispersion within the virial radius. We plot this quantity as a function of the halo mass within the same radius. The resulting scatter plot is shown in Fig. 12, where the scatter beyond the power-law relation appears to be very small.

Fig. 12. FJR for halos with mass between ≈105 and ≈2 × 108 M⊙ in the D–M universe, using simulations of a modified version of the RAMSES code. The mass versus virial velocity relation exhibits a power-law behavior with a very small scatter. The exponent in the power-law relation m ∝ vα is ≈2.6 when interpolating across the whole mass range, with α ≈ 2.45 and ≈2.85 in the lower and the upper part of the plot, respectively. This behavior is similar to the power-law exponent of 3 favored by ΛCDM, but somewhat smaller than the exponent 4 predicted by MOND.

We consider this first analysis to be very encouraging, but we note that the FJR can be obtained also in the context of ΛCDM, with a similar exponent. An important difference lies in the fact that in D–M, there is no dark matter, and that therefore the scaling relation between baryonic mass and virial velocity is tighter. On the other hand, we clearly need more realistic simulations, in particular regarding hydrodynamics and feedback, which will allow us to test the TFR for spiral galaxies, in addition to the FJR.

6. MOND-like behavior in the D–M universe

As mentioned previously, the MOND phenomenology can be justified if there exists gravitational polarization (Blanchet 2007; Blanchet & Le Tiec 2009), while D–M predicts such a polarization. We have indeed seen in Sect. 4 that the D–M universe exhibits a MOND-like behavior for idealized, “spherical” galaxies. On the other hand, the question remains open whether the value of the additional field created by antimatter and predicted by D–M can justify the characteristic acceleration a₀ ≈ 1.2 × 10⁻¹⁰ m s⁻² that seems to best fit the MOND behavior. In the present section, we study the value of this additional gravitational acceleration and show that this value cannot be reduced in D–M to a constant and that, in particular, it depends on the redshift at which it is measured.

In order to study the intensity of the additional gravitational field created by the antimatter clouds on the outskirts of galaxies and clusters, we compare the modulus |g_m| of the gravitational field produced by matter alone given in Eq. (6) to |g₊|, the total gravitational field acting on matter, and produced both by matter and antimatter as in Eq. (4). Similarly, we use the gravitational field produced by antimatter alone |g_am| to quantify the additional gravitational field present in D–M, due to the gravitational polarization between matter and antimatter, mimicking a MOND acceleration transition.

In the top panel of Fig. 13, the total gravitational field |g₊| acting on a matter particle, created by both matter and antimatter, is plotted as a function of the gravitational field |g_m| that would be created by matter if it were alone. These gravitational fields are calculated on a Cartesian 256³ grid (about 16 million cells), and the quantities estimated on the matter region (occupying ≈50% of the grid volume), which is the only cut applied for the selection. The figure clearly shows two regimes: the high field part on the right part of the figure shows the Newtonian regime, where the gravitational field |g₊| matches almost exactly the expected Newtonian gravitational field |g_m| created by matter. It can be seen that for accelerations smaller than ≈3 × 10⁻¹¹ m s², a non Newtonian behavior appears, with a flattening of the observed acceleration, similar to the MOND behavior (but, on this figure, for a different value of the acceleration parameter a₀). Although the dispersion increases at low accelerations, it can be seen that, in this situation where the halo mass range exceeds three orders of magnitude, the vast majority of the points lie above the diagonal, indicating that the antimatter field |g_am| reinforces quite generally the matter field |g_m|, mimicking a MOND behavior.

Fig. 13. Top panel: scatter diagram showing the relation between the gravitational field |gm| that would be created by matter if it were alone, on the x-axis, as a function of the total gravitational field |g+| acting on a matter particle, created by both matter and antimatter, on the y-axis. While the high field part on the right part of the figure shows the expected Newtonian behavior, it can be seen that for accelerations smaller than ≈3 × 10−11 m s−2, a non Newtonian behavior appears, with a flattening of the acceleration observed, in a behavior analogous to the MOND behavior. The MOND interpolating function used in Lelli et al. (2019) with a value of the a0 parameter equal to 0.4 × 10−11 m s−2 has been superimposed on the simulation data, and the two adjacent curves represent the 1-σ error of 0.11 dex found by these authors. Bottom panel: scatter diagram showing the relation between the gravitational field |gam| that would be created by antimatter if it were alone, on the x-axis, and the total gravitational field |g+| acting on a matter particle, created by both matter and antimatter, on the y-axis. It can be seen on this figure that the gravitational field created by antimatter is much more uniform, with values smaller by typically one order of magnitude, than the Newtonian regime of the previous figure, created for the most part by matter.

In the bottom panel of Fig. 13, we represent the modulus of the total gravitational field |g₊| acting on a matter particle, created by both matter and antimatter, as a function of the gravitational field |g_am| that would be created by antimatter if it were alone. We can see that the additional acceleration created by antimatter in the matter and depletion zones is significantly smaller on average than the matter gravitational field observed in the Newtonian regime, and has a relatively peaked distribution, at ≈5 × 10⁻¹² m s⁻² in the cluster simulation presented.

Clearly, the MOND-like transition observed in Fig. 13 at an acceleration of ≈3 × 10⁻¹¹ m s⁻² is here significantly smaller than the value of ≈1.2 × 10⁻¹⁰ m s⁻² of the fundamental acceleration a₀ favored by MOND. Therefore, one may ask whether the D–M matter-antimatter scenario has any relevance as a possible explanation of MOND. In order to answer this question, we note that the average modulus of the gravitational field created by antimatter depends on the redshift and on the largest scale of structure at this redshift. Fundamentally, the additional confining gravitational field created by antimatter, acting most of the time coherently with the matter field (see Fig. 13), is almost independent of the mass of the individual galaxy considered, and equal to the antimatter field created at the scale of the largest structures.

In order to emphasize this fundamental aspect, we have plotted in Fig. 14 the evolution as a function of the scale factor a(t) of the distribution of the gravitational field |g_am| created by antimatter for two simulations, the first one of “galactic” size and total matter mass of about 2.2 × 10¹⁰ M_⊙, the second one of “cluster” size and total matter mass of about 10¹⁶ M_⊙. These two figures present the following features: the distribution of the antimatter additional field is rather peaked at all redshifts. But although its overall shape is almost independent of redshift, reflecting the hierarchical build-up of gravitational structures, it is clear that the value of the peak depends not only on the redshift, but also, surprisingly, on the size of the simulation box, as long as it is smaller than a few homogeneity scales.

Fig. 14. Evolution as a function of the scale factor a(t) of the distribution of the modulus of the gravitational field |gam| created by antimatter for two simulations. In the first simulation, on the top panel, the total mass of the simulation is of galactic size, ≈2.2 × 1010 M⊙. The second simulation, on the bottom panel, is of cluster size and has a total mass of ≈1016 M⊙. Both distributions are rather peaked at all redshifts, and with an overall shape and width largely independent of redshift. On the other hand, the numerical value of the peak of the distribution clearly depends on the redshift. For both simulations, of limited size and mass, the peak value of the antimatter field differs from the fundamental constant postulated by MOND, a0 ≈ 1.2 × 10−10 m s−2.

The latter implies that the value of the extra gravitational field due to antimatter in our Universe, much higher than in these two simulations, is still not numerically converged. We are currently designing a larger simulation with a 4096³ resolution and linear comoving dimensions of ≈1 Gpc, which is expected to sample the homogeneity scale in D–M, and therefore to estimate more precisely the antimatter confining field predicted by D–M in our Universe. Ideally, the simulation volume would extend over a distance of a few gigaparsecs, which is inaccessible at present.

Using the non linear structure formation of D–M as a function of redshift shown in Fig. 4 of Manfredi et al. (2020), we provide at this stage an estimate of the evolution of the average modulus |g_am| of the gravitational field created by antimatter as a function of redshift for a simulation of a region of a D–M universe of cosmological dimensions, extending much beyond the homogeneity scale of ≈200 Mpc.

The quantity:

$$\begin{array}{c}\hfill |{\mathit{g}}_{\mathrm{am}}|\approx \frac{4\pi G2\rho (z)}{3}\frac{\pi }{k_{\mathrm{peak}}(z)}\end{array}$$(20)

where ρ(z) is the average density at redshift z, and π/k_peak(z) is the scale of the largest structures at redshift z, provides an estimate of the gravitational field created by antimatter at this redshift z. It should be noted that the density contrast between matter and antimatter regions is on the order of unity, and for this reason ρ enters directly in the above formula.

The evolution of this quantity is represented in Fig. 15. This figure shows two striking features: the modulus of the additional gravitational field created by antimatter is on the order of 2 × 10⁻¹¹ m s⁻² at the present epoch (z ≈ 0). The conversion of this antimatter gravitational field to the a₀ parameter of MOND (using the parametrizations of Famaey & McGaugh 2012) and the present uncertainty on the size of the initial domains of matter and antimatter results in a range of values of a₀ between ≈4 × 10⁻¹¹ m s⁻² and ≈2 × 10⁻¹⁰ m s⁻², a rather striking similarity with MOND. It should be noted that this quantity clearly depends on the redshift, with a decreasing trend in the future, reflecting the fact that the hierarchical formation of large-scale structures has almost completely stopped at the present epoch in the D–M universe, while the adiabatic expansion continues. The prediction on the variation of this transition acceleration, which can be extracted from Fig. 15, might be used in observations of high-redshift galaxies to discriminate between the MOND and D–M scenarios. More precisely, once averaged over the range a = [0.1, 1.0], the scale dependence of the antimatter gravitational field seems well approximated by a power-law g_am ∝ a^−1.8. This might be tested by rotation curve data at relatively high redshifts, for example using the data from Genzel et al. (2017).

Fig. 15. Evolution as a function of the scale factor a(t) of the average modulus of the gravitational field created by antimatter in the D–M universe. The scale factor is normalized such that at the present epoch the scale factor a = 1. At our epoch, the value of the average modulus of the antimatter gravitational field is of a few 10−11 m s−2, leading to a value of the a0 acceleration parameter similar to that postulated by MOND.

7. Conclusions and perspectives

The present study has evidenced new elements of concordance between our Universe and the D–M universe. In particular, it proposes an explanation for the observation of flat rotation curves in galaxies, which is usually attributed to the presence of dark matter or to a modification of the laws of gravitation akin to MOND. We have seen in particular that, due to the combined influence of the depletion zone and the antimatter clouds, flat rotation curves are generic in the D–M universe, leading to a systematic overestimate of the mass present in galaxies and clusters beyond a few (⪆2.5) virial radii.

The additional force experienced in D–M differs both from the interpretation of MOND, with its modified expression for the gravitational field using a fundamental acceleration constant a₀, or from the local expression of the Λ cosmological term conjectured by Gurzadyan (1985) and Gurzadyan & Stepanian (2019). On the other hand, in the same spirit as Gurzadyan, we propose a common explanation to the tentative dark energy and dark matter components of ΛCDM, using in D–M a single constant G, instead of two constants G and Λ. This provides an explanation for the otherwise rather extraordinary and fine-tuned coincidence of the dark energy and dark matter densities by linking them to the dynamical evolution of the matter and antimatter components, leading in D–M to a coasting universe, with a(t)∝t (cf. Eq. (1)).

We have also noted that the gravitational polarization between the positive and negative mass components (Price 1993) is at the origin of a MOND-like behavior. Blanchet (2007) and Blanchet & Le Tiec (2009) have indeed shown that gravitational polarization could explain the MOND phenomenology, although these authors did not have in mind that this gravitational polarization could be due to antimatter. On the other hand, our analysis differs significantly from that of Hajdukovic (2011, 2014, 2020), which seems to lead to unobserved effects in the Solar system (Banik & Kroupa 2020) and from the analysis by Penner (2016), who conjectured that MOND could be justified by the gravitational polarization of the vacuum, without taking into account the gravitational polarization of matter and antimatter structures.

Additionally, using both a simple analytical model and RAMSES simulations, we showed that the D–M cosmology predicts a rather well-defined power law between mass on the one hand, and rotation or virialized velocity, on the other, in a large mass range of gravitational structures. This may provide an explanation for the impressive correlation evidenced in the TFR and FJR in galaxies and clusters, although the exponent α in the relation m ∝ v^α was shown to be closer to 3, rather smaller than the preferred value of 3.85 ± 0.09 found by Lelli et al. (2019), and the value of 4 predicted by MOND (see also McGaugh (2012) for an earlier work focusing on gas-rich galaxies).

In future studies, we intend to realize a more realistic treatment of structure formation by including hydrodynamics and feedback in our RAMSES simulations, that will allow us to study the TFR, in addition to the FJR. We intend also to study larger simulation volumes, extending beyond the homogeneity scale (≳200 Mpc) predicted by D–M, and observed in our Universe (see, however, Keenan et al. 2013 and Haslbauer et al. 2020 for an example of inhomogeneity extending somewhat beyond this scale). Such large-scale simulations are a challenge for the D–M cosmology, as the “domains” of matter and antimatter assembled at decoupling, which initiate the formation of larger structures, are of limited geometrical extension (on the order of 100 parsec at z = 1080, i.e., ≈100 ckpc). It is important to note that the domains are in the non linear regime (density contrast of order unity) almost immediately after the CMB transition, which we will study in a forthcoming publication. This requires a very high resolution compared to the usual cosmological simulations. These simulations at scales larger than a few hundred Mpc will allow us to test the deviations from local to large-scale measurements of the Hubble parameter H₀, a question of paramount interest in the present context, where the tensions on this parameter seem to reach or even exceed the 5σ level between “local” and cosmological measurements (Riess 2020; Di Valentino 2021).

Concerning the average value of the additional gravitational field created by the antimatter clouds, we have shown that the estimate of the transition acceleration a₀ at our epoch (z = 0) is on the order of 10⁻¹⁰ m s⁻², a striking similarity with the MOND formalism. Fundamentally, we note that the value of this additional “MOND-like” gravitational field is determined by the antimatter field created at the largest structure scale, on the order of 200 Mpc. Also, we have shown that the modulus of this parameter depends on the redshift and is therefore not a fundamental constant, differing fundamentally from the MOND formalism.

Finally, we note that both nucleosynthesis and the almost purely non linear structure formation in the D–M universe (Benoit-Lévy & Chardin 2012; Manfredi et al. 2018, 2020) set strong constraints on the size of initial “domains” of matter and on the later hierarchical (largely bottom-up) development of structures. A complementary way to test our hypothesis will be to predict the mass distribution of stars and black holes resulting from the very early collapse of such matter domains. Compared to the black hole mass distribution derived by the LIGO and Virgo collaborations (The LIGO Scientific Collaboration & The Virgo Collaboration 2021), now in possession of about 70 candidate events (mostly binary black holes), this could represent an important additional test of the D–M scenario.

Acknowledgments

We are indebted to the anonymous referee and to Benoit Famaey, James Rich and Yves Sacquin for their thorough reading of the manuscript and their insightful comments. Needless to say, they are not responsible for the errors and approximations remaining in this paper. This work has made use of the Horizon Cluster hosted by the Institut d’Astrophysique de Paris. The work of the YT (Turk et al. 2010), IPython (Perez & Granger 2007), Matplotlib (Hunter 2007), NumPy (van der Walt et al. 2011) and SciPy (Virtanen et al. 2020) development teams is also gratefully acknowledged. The halo catalogs have been computed using the AdaptaHOP algorithm (Aubert et al. 2004).

References

Appendix A: Are there external field effects (EFE) in the D–M cosmology

The previous discussion on the rotation curves and the virial velocity distributions in the D–M universe has evidenced the surprising property that in this cosmology, to a good approximation, around a central massive object, a harmonic restoring force is felt in addition to the usual Gm/r² Newtonian force, leading to nearly flat rotation curves. At distances comparable to the radius of the depletion zone, this additional harmonic restoring force is far from negligible since it is on average equal to the force of the central galaxy (see Fig. 7). It may seem that this violates blatantly the shell theorem of Newton (1760) and, in General Relativity, the Birkhoff theorem (Birkhoff & Langer 1923). In the present appendix, we discuss this question, which presents interesting features.

The first element of answer comes from noticing that, in order to describe the gravitational behavior of matter and antimatter in the D–M universe, two coupled Poisson equations are required:

$$\begin{array}{cc}\hfill {\mathrm{\nabla }}^2{\varphi }_{+}& =4\pi G({\rho }_{+}-{\rho }_{-}),\hfill \end{array}$$(A.1)

$$\begin{array}{cc}\hfill {\mathrm{\nabla }}^2{\varphi }_{-}& =4\pi G(-{\rho }_{+}-{\rho }_{-}).\hfill \end{array}$$(A.2)

The possibility to express the D–M gravitational behavior with coupled equations using Laplacian operators implies, through the Gauss theorem, that in a situation of spherical symmetry, the gravitational field in an empty region must be zero, since the mass enclosed is zero and spacetime in this region should be Minkowskian. However, there is a caveat to this assertion, related to the boundary conditions of the mass configuration. The previous statement on the shell theorem being respected may indeed seem at odds with our decomposition in terms of three spheres (see Fig. 3) of the average galaxy environment in the D–M universe, where a harmonic restoring force is observed in addition to the usual Newtonian force within the depletion zone. The short answer to this apparent contradiction is that the situation lacks spherical symmetry and that the harmonic restoring force that we have derived just approximates this asymmetric configuration of the antimatter cloud surrounding a spherical depletion zone (see in particular Fig. 2). Note in particular that in actual configurations, the depletion zones are not spherical, but percolate with surrounding depletion zones, with a similar percolation property for the antimatter clouds. Also, the approximation used is valid only at distances smaller than r_d, where r_d is the approximate size of the depletion zone surrounding the massive structure, and becomes increasingly inaccurate when we exceed radial distances larger than ≈r_d/2.

Although this asymmetry correctly answers the question on the violation of the shell theorem in galactic configurations in the D–M universe, it may hide some interesting elements of discussion, which we summarize in the following.

The first element comes from considering the situation resulting from the superposition of only the two cubes (a) and (c) of Fig. 3, that is the cube with uniform repulsive background and the cube containing a sphere with uniform positive mass density compensating, within the volume of the sphere, the negative mass background of the first cube. As soon as we have accepted the property that the gravitational field created by cube (a) with uniform density is necessarily zero everywhere, which seems unavoidable by symmetry, it is also clear that the contribution of the second cube will create a harmonic restoring force, although the inner sphere, in the superposition of the two cubes, is now empty.

The situation is even stranger when we consider the configuration with complete spherical symmetry, where the whole space is filled with a uniform negative background, to which we superimpose a sphere centered on the origin with a positive uniform density, compensating the negative mass fluid inside the sphere (and only there). This time we cannot invoke the asymmetry of the situation and it seems that we have a gross violation of the shell and Birkhoff theorems (Newton 1760; Birkhoff & Langer 1923) since the gravitational field appears to be nonzero in an empty region with exact spherically symmetry.

However, the expression of Birkhoff’s theorem (Birkhoff & Langer 1923) only states that any spherically symmetric solution of the vacuum field equations must be static and asymptotically flat, and represented by a Schwarzschild metric. We must then note that we have in fact filled out the entire space with a negative mass fluid of constant density, and therefore with infinite negative global mass. This configuration is necessarily not static, so the “cosmological” aspects must now be taken into account in the dynamical situation. Our decomposition into three cubes with respective masses −2m, +m and +m, on the other hand, restores a total mass zero, and has not a diverging mass and potential at infinity (Seeliger 1895), but involves a configuration without spherical symmetry.

Two important additional comments can be made. The first comment is based on the remarkable analysis by Gurzadyan as early as 1985 (Gurzadyan 1985), in the early days of dark matter searches, which in several respects reaches conclusions similar to those of the present analysis⁶.

In his first publication on this topic (Gurzadyan 1985), Gurzadyan notes that the Newtonian 1/r² force law could be extended by requiring the property, realized in the Newtonian case, that for a configuration with spherical symmetry, the gravitational action of a mass can be reduced to the situation where all the mass is concentrated at the origin. Remarkably, this requirement leads not only to the Newtonian potential, but also to an additional harmonic force, attractive or repulsive depending on the sign of Λ, which appears as the fluid analog of the cosmological constant, and we reproduce here, with a slight change of notation, Eq. (5) of Gurzadyan (1985):

$$\begin{array}{c}\hfill F(r)=A{r}^{-2}+\mathrm{\Lambda }r.\end{array}$$(A.3)

In this sense, the introduction of a cosmological constant in Newtonian cosmology finds a natural justification in the Milne cosmology (Milne 1933). Gurzadyan notes that this second component has the property that the shell theorem (zero gravitational field inside an empty spherical shell) is not respected. More importantly, in recent publications, (see e.g., Gurzadyan & Stepanian 2019, and references therein), Gurzadyan conjectures that the flat rotation curves observed in galaxies can be explained by the same (local) expression as the cosmological repulsive term. However, with his expression of a generalized gravitational potential using two gravitational constants G and Λ, Gurzadyan does not consider the possibility that these two terms could be expressed with a single constant G but with a sign reversal. For example, the potential term in the GR metric in Eq. (1) of Gurzadyan & Stepanian (2019) leads directly in Eq. (3) to the two components, Newtonian and harmonic restoring force, of the gravitational field that we have derived previously for D–M.

The second comment is related to the fact that, as mentioned previously, cosmological simulations of self-gravitating structures routinely use effective negative mass without explicitly stating it. Indeed, as mentioned previously, cosmological simulations first calculate the average density in order to derive the average cosmological expansion of a(t), the cosmological scale parameter. After subtraction of this average density, a symmetric (in terms of sign of mass) mass distribution is then obtained, usually with symmetrical overdensity and underdensity distributions. Piran indeed remarked (Piran 1997), following previous numerical simulations of gravitational structures with Dubinski and collaborators (Dubinski et al. 1993), that while we are used to representing gravitational structures in terms of (collapsing) positive mass, it is also possible and useful to consider the problem in terms of (expanding) voids of negative mass. In their expansion, these voids will effectively develop into the largest structures in the universe. As studied later by Sheth & Van De Weygaert (2004), now followed by the work of several authors (for a review, see for example Pisani et al. 2019 and references therein), the study of voids constitutes today a powerful test of the cosmology at play in our Universe. Indeed, the characteristics of the KBC void (Keenan et al. 2013) appear to rule out the ΛCDM cosmology at the 6σ level (Haslbauer et al. 2020).

We further note, following Piran, that the matter surrounding an underdense region will create in its expansion a high density ridge along the rim of the underdense region, and density even diverging at shell crossing (see e.g., Figs. 1 and 2 of Dubinski et al. 1993). A region nearly empty, or even totally empty, will therefore, in a cosmological context, effectively “repel” the surrounding matter, due to its faster expansion compared to its surroundings. Indeed, the Hubble tension might be resolved by taking into account the outflow from the KBC void (Haslbauer et al. 2020).

Piran goes as far as to suppose that there could exist a species of matter with negative gravitational mass and positive inertial mass, violating the usual expression of the weak equivalence principle. But curiously, Piran describes the interactions of this new species endowed with negative gravitational mass and positive inertial mass as attractive between themselves. We note that this is not the behavior actually observed in “voids” as the negative density (underdensity) flattens out instead of becoming more negative. Of course, under the ordinary assumption that only positive mass particles exist, the maximum “negative” mass (underdense) region is limited by the condition of zero matter density. This leads to the fundamental new element introduced by the introduction of “negative” mass particles: a depletion zone develops, and we have seen the fundamental role that it plays in terms of mimicking the behavior of extended dark matter clouds or, alternatively, providing a MOND-like behavior. As we have studied in detail in Manfredi et al. (2018), if we want to describe the behavior of the repulsive underdense regions using a second species of “negative” mass, this is not possible in a Newtonian description, even with the three parameters of inertial, passive gravitational and active gravitational mass. In order to implement the gravitational behavior of the Dirac particle-hole system, let us stress again that it is necessary to use a bimetric description (Manfredi et al. 2018).

In conclusion, the harmonic restoring force predicted by the D–M model, adding its contribution to the usual Newtonian force term, is an unavoidable consequence of the combined influence of the antimatter cloud and empty depletion zone developing in the D–M universe around a localized structure of positive mass. But unlike the additional Λ term conjectured by Gurzadyan, this harmonic restoring force does not violate the shell theorem. This unexpected restoring force is simply explained in terms of the asymmetric distribution of the positive and negative mass components in the D–M cosmology.

All Figures

Fig. 1. Slice of a matter-antimatter simulation of a small volume in the D–M universe using a modified version of the RAMSES simulation code. The thickness of the slice is 20% of the simulation box, of comoving size 1 h−1 cMpc. The condensed structures (represented in light blue) are all made of matter, while antimatter is spread in nearly homogeneous extended halos (represented in red) over half of the volume. The other half of the volume is occupied by depletion zones, surrounding the matter galaxies and clusters, isolating them from antimatter. These depletion zones are a consequence of the gravitational polarization expected between positive mass and “negative mass” objects. We note that antimatter clouds percolate, i.e., it is possible to roam at infinity without leaving the antimatter cloud, and the same property is respected for the depletion zones surrounding matter.

In the text

Fig. 2. Schematic representation of a periodic structure of identical galaxies and their surroundings in the D–M universe. A galaxy core is represented by a condensed object at the center of a cubic cell, while a depletion zone, occupying about 50% of the cell volume, surrounds the galaxy. Antimatter is spread out beyond this depletion zone, on the outskirts of the cell, occupying also nearly 50% of the cell volume, although its volume appears smaller in the 2D projection.

In the text

Fig. 3. Our approximation of a galaxy represented on Fig. 2 can be represented as the sum of three cubes: a cube (a) with uniform “negative mass” density of total mass −2m, a cube (b) with only the galaxy of mass +m at the center, a cube (c) with a sphere occupying half the volume of the cube and of homogeneous positive density of matter, equal and opposite to that of the first cube, that together with (a) will create the empty depletion zone.

In the text

Fig. 4. Rotation velocity predicted in the D–M universe for a point mass located at r = 0 and surrounded by a depletion zone and the same mass of antimatter (with opposite “sign”). Instead of the expected Newtonian velocity (see Eq. (11)), decreasing to zero, D–M predicts that the rotation velocity for a condensed object of positive mass will be nearly constant (see Eq. (14)) above the characteristic distance approximately equal to half the radius rd of the depletion zone, or ≈3 virial radii.

In the text

	Fig. 5. Configuration studied in this simulation. A matter galaxy, located at the center of the box, is surrounded by a depletion zone, extending over half the volume of the box, and by an antimatter cloud (represented in magenta), extending to the outskirts of the box over the other half of the volume.
In the text

Fig. 6. Gravitational field gam created by antimatter on two planes defined by their normal vector (1,0,0) for the top panel and $\frac{(1,1,1)}{\sqrt{3}}$ for the bottom panel. The color scaling indicates the strength of the field. The confining gravitational field created by the antimatter cloud is more symmetric and has a higher intensity in the bottom panel, which represents the field in one of the diagonal planes of the cube.

In the text

Fig. 7. Top panel: total gravitational field acting on a matter particle as a function of the gravitational field created by matter alone. A Newtonian regime can be observed for the high field region of the figure, with a gradual transition from the Newtonian regime to a MOND-like regime at an acceleration of ≈1.0 × 10−9 m s−2. Middle panel: total gravitational field acting on a matter particle as a function of the gravitational field created by antimatter alone. The antimatter field is rather peaked, with an average value of ≈2 × 10−10 m s−2. Bottom panel: ratio between the total gravitational field acting on a matter particle and the Newtonian field (created by matter only), as a function of the gravitational field created by matter only. Almost exactly unity in the Newtonian regime (at high values of the field), this ratio gradually increases and reaches a factor of ≈2 at the end of the depletion zone. The MOND interpolating function used in Lelli et al. (2019) for a value of the a0 parameter 1.85 × 10−10 m s−2 has been superimposed on the simulation data, while the two adjacent curves represent the 1-σ error of 0.11 dex found by these authors.

In the text

	Fig. 8. Matter density ρ+(r) (black dashed line), antimatter density ρ−(r) (blue line), and pseudo-dark-matter halo density: ${\rho }_{\mathrm{halo}}=2\overline{\rho }-{\rho }_{-}$ (red line).
In the text

	Fig. 9. Gravitational fields: g− (blue line), −g+ (black dashed line), and −gNewt (red line).
In the text

	Fig. 10. Rotation speeds for the D–M model (black line) and for standard Newtonian gravity (red line). The decrease of the rotation velocity beyond r/rd = 1 is non physical, as the test particle would enter the antimatter cloud.
In the text

Fig. 11. Logarithmic spherical density as a function of radius for both matter and antimatter for a set of 25 halos with mass ≈107 M⊙ in our RAMSES simulation. While the density of matter, represented in blue, is steeply decreasing from the halo center, the density of antimatter, represented in orange, is zero in this region. The logarithmic density is scaled to the average matter and antimatter density, equal in D–M at large scales, and converge (equal average densities for matter and antimatter) for almost all halos. It can be seen that the antimatter density becomes non zero after a few virial radii, and that the mass of antimatter becomes on the order of the halo matter mass after typically twice the radius at which antimatter density becomes non zero. The distribution for 25 consecutive halos has been represented in the figure in order to show that this behavior is quite generic in D–M, although some variations can be seen in the distribution of antimatter around the matter halos.

In the text

Fig. 12. FJR for halos with mass between ≈105 and ≈2 × 108 M⊙ in the D–M universe, using simulations of a modified version of the RAMSES code. The mass versus virial velocity relation exhibits a power-law behavior with a very small scatter. The exponent in the power-law relation m ∝ vα is ≈2.6 when interpolating across the whole mass range, with α ≈ 2.45 and ≈2.85 in the lower and the upper part of the plot, respectively. This behavior is similar to the power-law exponent of 3 favored by ΛCDM, but somewhat smaller than the exponent 4 predicted by MOND.

In the text

Fig. 13. Top panel: scatter diagram showing the relation between the gravitational field |gm| that would be created by matter if it were alone, on the x-axis, as a function of the total gravitational field |g+| acting on a matter particle, created by both matter and antimatter, on the y-axis. While the high field part on the right part of the figure shows the expected Newtonian behavior, it can be seen that for accelerations smaller than ≈3 × 10−11 m s−2, a non Newtonian behavior appears, with a flattening of the acceleration observed, in a behavior analogous to the MOND behavior. The MOND interpolating function used in Lelli et al. (2019) with a value of the a0 parameter equal to 0.4 × 10−11 m s−2 has been superimposed on the simulation data, and the two adjacent curves represent the 1-σ error of 0.11 dex found by these authors. Bottom panel: scatter diagram showing the relation between the gravitational field |gam| that would be created by antimatter if it were alone, on the x-axis, and the total gravitational field |g+| acting on a matter particle, created by both matter and antimatter, on the y-axis. It can be seen on this figure that the gravitational field created by antimatter is much more uniform, with values smaller by typically one order of magnitude, than the Newtonian regime of the previous figure, created for the most part by matter.

In the text

Fig. 14. Evolution as a function of the scale factor a(t) of the distribution of the modulus of the gravitational field |gam| created by antimatter for two simulations. In the first simulation, on the top panel, the total mass of the simulation is of galactic size, ≈2.2 × 1010 M⊙. The second simulation, on the bottom panel, is of cluster size and has a total mass of ≈1016 M⊙. Both distributions are rather peaked at all redshifts, and with an overall shape and width largely independent of redshift. On the other hand, the numerical value of the peak of the distribution clearly depends on the redshift. For both simulations, of limited size and mass, the peak value of the antimatter field differs from the fundamental constant postulated by MOND, a0 ≈ 1.2 × 10−10 m s−2.

In the text

Fig. 15. Evolution as a function of the scale factor a(t) of the average modulus of the gravitational field created by antimatter in the D–M universe. The scale factor is normalized such that at the present epoch the scale factor a = 1. At our epoch, the value of the average modulus of the antimatter gravitational field is of a few 10−11 m s−2, leading to a value of the a0 acceleration parameter similar to that postulated by MOND.

In the text