# How can I add dark matter to my $N$-body simulation?

I've written a simple non-scientific N-body simulation for fun: http://magwo.github.io/fullofstars/

I expected to create something looking like spinning galaxies (there are two invisible very heavy objects in that simulation).

However, as I learned, for it to resemble galaxies there needs to be something more than just $N$-body simulation. http://en.wikipedia.org/wiki/Galaxy_rotation_curve

So my question is this: For my simulation to look like spinning galaxies, how should the dark matter be positioned, moving, and how should it interact with the visible particles?

Should the dark matter be at all affected by the gravity of the visible particles? Or should it be a one-way interaction? What makes the dark matter stay "in place" to preserve its "compressing" effect on the particle cloud? Is the gravity interaction from the dark matter an attraction force just like normal gravity, or is it a repelling force?

• Whoa!!! How many questions are you really asking? You have three followed by another five. No offence intended, but would you consider splitting them up? Better chance of getting answers that way. Have you read this? en.m.wikipedia.org/wiki/Dark_matter Regards – user81619 Jul 28 '15 at 10:38
• Hehe well! The followup questions are more like specific clarifying questions of the main question. Yes I've skimmed the wikipedia article but wasn't able to find any specific proposals for how dark matter interacts with regular matter, and especially, if it is affected by regular matter in the same way as regular matter is affected by other regular matter. – Magnus Wolffelt Jul 28 '15 at 10:49
• Hi Magnus, I don't think it works like that. Some of your questions are complicated, for them to be answered properly, so you can model things accurately, will take time. Search this site, there may be answers here already. Best of luck with it though. – user81619 Jul 28 '15 at 10:57
• What do your normal matter particles do other than gravitate? If the answer is "nothing," then they are not terribly different from dark matter. On large scales where everything is like a fluid, the only difference between dark and normal matter is that the dark matter is pressureless. – user10851 Jul 28 '15 at 14:33
• @ChrisWhite ok nice. But I fail to understand how normal gravity can achieve the "compression" effect needed for high-velocity high-altitude stars to not escape the galaxy fairly quickly. Surely, adding more mass is not the answer, as that would also increase the gravity on the inner stars, requiring them to move much faster. And a realistic galaxy rotation would then not be achieved. Edit: and yes, the particles do not do anything except attract each other. – Magnus Wolffelt Jul 28 '15 at 14:39

As I understand your question, and after watching your simulation, it seems that you actually only have dark matter (DM) in your simulation. An N-body code simulates DM, i.e. particles that interact only by gravitation. If you want to take it a step further, you want to include gas, i.e. particles that interact both by gravitation and by hydrodynamics.

This is not trivial, but the basic idea is that when gas particles get close enough, they build up pressure, are slowed down by viscosity forces, experience turbulence, heat up and lose energy through radiation, etc. If a region becomes sufficiently dense, it may form stars, which in turn heat up the gas again.

When a code is particle-based (as opposed to grid-based), this is called smoothed particle hydrodynamics (SPH).

However, I think you can start by doing something more simple. I may be wrong, but from looking at your simulation, it seems to me that you are calculating gravitational forces between particle $i$ and particle $j$ through $$F_{ij} = \frac{m_i m_j}{r_{ij}^2}, \qquad\qquad(\mathrm{correct \, in \, real \, life})$$ which of course in principle is correct. However, in numerical simulations, due to finite step sizes and floating point precision, particles easily come to close, resulting in diverging forces, which tend to sling out particles with unnaturally high velocities. For this reason, people use a so-called gravitational softening length $\varepsilon$, such that the force between the particles become $$F_{ij} = \frac{m_i m_j}{r_{ij}^2 + \varepsilon^2} \qquad\qquad(\mathrm{correct \, in \, simulations}).$$ The optimal value of $\varepsilon$ is a matter of choice. You may want to try different values, but it should be of the order of (or a little smaller than) the mean distance between particles in the regions you want to resolve. Your time steps should also be adapted to this, so that particles on average don't travel more than $\varepsilon/2$ per timestep (e.g. Rodionov & Sotnikova 2005). For instance, try and start by setting $\varepsilon = L/100$, where $L$ is the length of your box.

Including this parameter should make your DM particles clump more. However, if you want beautiful spiral patterns, you have to include gas.

• I think that the softening of the newtonian force is a zero-order kind of regularization. It is pretty easy to implement and simple to grasp, but I don't see how one can estimate how good this approximation is. I think the "matter of choice" argument is good enough for non-physical "nice-to-see" application, but is of course of no help if one is concerned with real physics. – sintetico Jul 29 '15 at 10:22
• @sintetico: There are ways of estimating an optimal $\varepsilon$. For instance you can argue that it should be high enough to prevent divergence, but low enough that it affects test benchmark results only to within some chosen tolerance. Plenty of authors that do real physics has delved into this, e.g. these, these, and these. – pela Jul 29 '15 at 12:03
• So, if I am right, corrections are evaluated by a confrontation between known "benchmark". This is a kind of euristic, but how can one be sure that the approximation will work also in the general case. I mean, is there a way to estimate the correction? Also consider that in most of the cases gravitational systems are chaotic, which means that a small computational error or approximation can have huge effects on the outcome. Please consider that I am not an expert on this problem, and that I am just curious on this problem. – sintetico Jul 29 '15 at 12:21
• I mean, your "softening" method is much better than regularization, because is much simpler. So can one use the "softening" in any case? Are there caveats to consider to apply this method? – sintetico Jul 29 '15 at 12:23
• @sintetico: Okay, I had a closer look, and the Plummer model is not the same as "mine" (in lack of a better name). The Plummer model, having $F\propto r/(r^2+\epsilon^2)^{3/2}$, increases for $r\rightarrow0$, but then decreases for very small $r$, ending at zero for $r=0$. In contrast, "my" model, having $F\propto1/(r^2+\epsilon^2)$, increases for $r\rightarrow0$ but ends at a finite value. – pela Aug 3 '15 at 9:19

The answer of Acid Jazz is correct in every point, in my opinion. However, I think that the problem in your simulations is not dark matter, but something more subtle. Perhaps I oversimplify here, but dark matter will behave in your model exactly like normal matter, with the exception of being dark, not visible. It is in fact exactly the case of the two invisible objects of your simulation. In fact, you are already considering two dark matter halos in the simulation.

For what I see, in your simulation the "stars" are accelerated and expelled at huge speeds after that they experience a close encounter with another star (or with one of the invisible objects). I presume that you are solving directly and separately the dynamics of each of these star, taking into account the mutual attraction. However, the newtonian potential has a singularity when two masses approach at short distance $r\approx0$ (the attraction is proportional to $r^{-2}$). This constitutes a computational problem in gravitational simulations, and is probably why you see these anomalies in your simulation, which prevent the emergence of large scale structures, like galaxies. A method to deal with this problem is known as Regularization (although I think it is not the only method).

I am not an expert in the field, I just found this link that seems to be a useful introduction to non-experts: Basics of regularization theory

My best guess for these questions: (but I am no expert so it's a guess, right?)

Should the dark matter be at all affected by the gravity of the visible particles?

Yes

Or should it be a one-way interaction?

No, gravity should act equally on visible and dark matter.

What makes the dark matter stay "in place" to preserve its "compressing" effect on the particle cloud?

Normal gravity probably (open to correction though), I am not exactly sure what else you mean, it is distributed in a halo around galaxies but extends much further out than visible matter. Also, it may be distributed in filaments and clumps without any visible matter nearby.

Is the gravity interaction from the dark matter an attraction force just like normal gravity, or is it a repelling force?

Attractive only

• Thanks for your reply! This clarifies some things. However, if the dark matter interacts with regular matter just like regular matter (if more weakly), then what makes it stay in this halo distribution around the galaxy? It seems like some sort of pressure/buoyancy effect is in place, no? – Magnus Wolffelt Jul 28 '15 at 11:10
• angular momentum – Bort Jul 28 '15 at 11:11
• And also, if the dark matter is in a halo around the galaxy, with attractive force, why is the speed of high-altitude stars so unexpectedly high? If anything, there needs to be a force pushing these high-velocity, high-altitude stars into the galaxy rather than attracting them out of it. – Magnus Wolffelt Jul 28 '15 at 11:11
• See what I mean, it's complicated :) so set them up as individual questions if you need to, but you are getting there. – user81619 Jul 28 '15 at 11:13

I've written a few small gravity simulations, and I would bet that the problem is that the particles aren't colliding. I know this doesn't answer your original question, but I think that's the answer to fixing your simulation.

## protected by Qmechanic♦Jul 29 '15 at 6:57

Thank you for your interest in this question. Because it has attracted low-quality or spam answers that had to be removed, posting an answer now requires 10 reputation on this site (the association bonus does not count).