softening
   HOME

TheInfoList



OR:

In
physics Physics is the natural science that studies matter, its fundamental constituents, its motion and behavior through space and time, and the related entities of energy and force. "Physical science is that department of knowledge which ...
and
astronomy Astronomy () is a natural science that studies celestial objects and phenomena. It uses mathematics, physics, and chemistry in order to explain their origin and evolution. Objects of interest include planets, moons, stars, nebulae, g ...
, an ''N''-body simulation is a simulation of a
dynamical system In mathematics, a dynamical system is a system in which a function describes the time dependence of a point in an ambient space. Examples include the mathematical models that describe the swinging of a clock pendulum, the flow of water i ...
of particles, usually under the influence of physical forces, such as
gravity In physics, gravity () is a fundamental interaction which causes mutual attraction between all things with mass or energy. Gravity is, by far, the weakest of the four fundamental interactions, approximately 1038 times weaker than the stro ...
(see ''n''-body problem for other applications). ''N''-body simulations are widely used tools in
astrophysics Astrophysics is a science that employs the methods and principles of physics and chemistry in the study of astronomical objects and phenomena. As one of the founders of the discipline said, Astrophysics "seeks to ascertain the nature of the h ...
, from investigating the dynamics of few-body systems like the
Earth Earth is the third planet from the Sun and the only astronomical object known to harbor life. While large volumes of water can be found throughout the Solar System, only Earth sustains liquid surface water. About 71% of Earth's sur ...
-
Moon The Moon is Earth's only natural satellite. It is the fifth largest satellite in the Solar System and the largest and most massive relative to its parent planet, with a diameter about one-quarter that of Earth (comparable to the width of ...
- Sun system to understanding the evolution of the
large-scale structure of the universe The observable universe is a ball-shaped region of the universe comprising all matter that can be observed from Earth or its space-based telescopes and exploratory probes at the present time, because the electromagnetic radiation from these obj ...
. In
physical cosmology Physical cosmology is a branch of cosmology concerned with the study of cosmological models. A cosmological model, or simply cosmology, provides a description of the largest-scale structures and dynamics of the universe and allows study of f ...
, ''N''-body simulations are used to study processes of non-linear structure formation such as galaxy filaments and galaxy halos from the influence of
dark matter Dark matter is a hypothetical form of matter thought to account for approximately 85% of the matter in the universe. Dark matter is called "dark" because it does not appear to interact with the electromagnetic field, which means it does not ...
. Direct ''N''-body simulations are used to study the dynamical evolution of
star cluster Star clusters are large groups of stars. Two main types of star clusters can be distinguished: globular clusters are tight groups of ten thousand to millions of old stars which are gravitationally bound, while open clusters are more loosely cl ...
s.


Nature of the particles

The 'particles' treated by the simulation may or may not correspond to physical objects which are particulate in nature. For example, an N-body simulation of a star cluster might have a particle per star, so each particle has some physical significance. On the other hand, a simulation of a gas cloud cannot afford to have a particle for each atom or molecule of gas as this would require on the order of particles for each mole of material (see
Avogadro constant The Avogadro constant, commonly denoted or , is the proportionality factor that relates the number of constituent particles (usually molecules, atoms or ions) in a sample with the amount of substance in that sample. It is an SI defining ...
), so a single 'particle' would represent some much larger quantity of gas (often implemented using
Smoothed Particle Hydrodynamics Smoothed-particle hydrodynamics (SPH) is a computational method used for simulating the mechanics of continuum media, such as solid mechanics and fluid flows. It was developed by Gingold and Monaghan and Lucy in 1977, initially for astrophysi ...
). This quantity need not have any physical significance, but must be chosen as a compromise between accuracy and manageable computer requirements.


Direct gravitational ''N''-body simulations

In direct gravitational ''N''-body simulations, the equations of motion of a system of ''N'' particles under the influence of their mutual gravitational forces are integrated numerically without any simplifying approximations. These calculations are used in situations where interactions between individual objects, such as stars or planets, are important to the evolution of the system. The first direct gravitational ''N''-body simulations were carried out by Erik Holmberg at the
Lund Observatory Lund Observatory is the official English name for the astronomy department at Lund University. Between 1867-2001 "Lund Observatory" was also the name of the Observatory building, which is now referred to as the "Lund Old Observatory". As of Janu ...
in 1941, determining the forces between stars in encountering galaxies via the mathematical equivalence between light propagation and gravitational interaction: putting light bulbs at the positions of the stars and measuring the directional light fluxes at the positions of the stars by a photo cell, the equations of motion can be integrated with effort. The first purely calculational simulations were then done by
Sebastian von Hoerner Sebastian Rudolf Karl von Hoerner (15 April 1919 – 7 January 2003) was a German astrophysicist and radio astronomer. He was born in Görlitz, Lower Silesia. After the end of World War II he studied physics at University of Göttingen. He obta ...
at the Astronomisches Rechen-Institut in
Heidelberg Heidelberg (; Palatine German: ') is a city in the German state of Baden-Württemberg, situated on the river Neckar in south-west Germany. As of the 2016 census, its population was 159,914, of which roughly a quarter consisted of students ...
, Germany. Sverre Aarseth at the
University of Cambridge , mottoeng = Literal: From here, light and sacred draughts. Non literal: From this place, we gain enlightenment and precious knowledge. , established = , other_name = The Chancellor, Masters and Schola ...
(UK) has dedicated his entire scientific life to the development of a series of highly efficient ''N''-body codes for astrophysical applications which use adaptive (hierarchical) time steps, an Ahmad-Cohen neighbour scheme and regularization of close encounters. Regularization is a mathematical trick to remove the singularity in the Newtonian law of gravitation for two particles which approach each other arbitrarily close. Sverre Aarseth's codes are used to study the dynamics of star clusters, planetary systems and galactic nuclei.


General relativity simulations

Many simulations are large enough that the effects of
general relativity General relativity, also known as the general theory of relativity and Einstein's theory of gravity, is the geometric theory of gravitation published by Albert Einstein in 1915 and is the current description of gravitation in modern physics ...
in establishing a Friedmann-Lemaitre-Robertson-Walker cosmology are significant. This is incorporated in the simulation as an evolving measure of distance (or scale factor) in a comoving coordinate system, which causes the particles to slow in comoving coordinates (as well as due to the
redshift In physics, a redshift is an increase in the wavelength, and corresponding decrease in the frequency and photon energy, of electromagnetic radiation (such as light). The opposite change, a decrease in wavelength and simultaneous increase in fr ...
ing of their physical energy). However, the contributions of general relativity and the finite
speed of gravity In classical theories of gravitation, the changes in a gravitational field propagate. A change in the distribution of energy and momentum of matter results in subsequent alteration, at a distance, of the gravitational field which it produces. In ...
can otherwise be ignored, as typical dynamical timescales are long compared to the light crossing time for the simulation, and the space-time curvature induced by the particles and the particle velocities are small. The boundary conditions of these cosmological simulations are usually periodic (or toroidal), so that one edge of the simulation volume matches up with the opposite edge.


Calculation optimizations

''N''-body simulations are simple in principle, because they involve merely integrating the 6''N''
ordinary differential equation In mathematics, an ordinary differential equation (ODE) is a differential equation whose unknown(s) consists of one (or more) function(s) of one variable and involves the derivatives of those functions. The term ''ordinary'' is used in contrast ...
s defining the particle motions in
Newtonian gravity Newton's law of universal gravitation is usually stated as that every particle attracts every other particle in the universe with a force that is proportional to the product of their masses and inversely proportional to the square of the distanc ...
. In practice, the number ''N'' of particles involved is usually very large (typical simulations include many millions, the Millennium simulation included ten billion) and the number of particle-particle interactions needing to be computed increases on the order of ''N''2, and so direct integration of the differential equations can be prohibitively computationally expensive. Therefore, a number of refinements are commonly used. Numerical integration is usually performed over small timesteps using a method such as
leapfrog integration In numerical analysis, leapfrog integration is a method for numerically integrating differential equations of the form \ddot x = \frac = A(x), or equivalently of the form \dot v = \frac = A(x), \;\dot x = \frac = v, particularly in the case of a ...
. However all numerical integration leads to errors. Smaller steps give lower errors but run more slowly. Leapfrog integration is roughly 2nd order on the timestep, other integrators such as
Runge–Kutta methods In numerical analysis, the Runge–Kutta methods ( ) are a family of implicit and explicit iterative methods, which include the Euler method, used in temporal discretization for the approximate solutions of simultaneous nonlinear equations. T ...
can have 4th order accuracy or much higher. One of the simplest refinements is that each particle carries with it its own timestep variable, so that particles with widely different dynamical times don't all have to be evolved forward at the rate of that with the shortest time. There are two basic approximation schemes to decrease the computational time for such simulations. These can reduce the
computational complexity In computer science, the computational complexity or simply complexity of an algorithm is the amount of resources required to run it. Particular focus is given to computation time (generally measured by the number of needed elementary operations) ...
to O(N log N) or better, at the loss of accuracy.


Tree methods

In tree methods, such as a
Barnes–Hut simulation The Barnes–Hut simulation (named after Josh Barnes and Piet Hut) is an approximation algorithm for performing an ''n''-body simulation. It is notable for having order O(''n'' log ''n'') compared to a direct-sum algorithm which would b ...
, an octree is usually used to divide the volume into cubic cells and only interactions between particles from nearby cells need to be treated individually; particles in distant cells can be treated collectively as a single large particle centered at the distant cell's center of mass (or as a low-order
multipole A multipole expansion is a mathematical series representing a function that depends on angles—usually the two angles used in the spherical coordinate system (the polar and azimuthal angles) for three-dimensional Euclidean space, \R^3. Similarly ...
expansion). This can dramatically reduce the number of particle pair interactions that must be computed. To prevent the simulation from becoming swamped by computing particle-particle interactions, the cells must be refined to smaller cells in denser parts of the simulation which contain many particles per cell. For simulations where particles are not evenly distributed, the
well-separated pair decomposition In computational geometry, a well-separated pair decomposition (WSPD) of a set of points S \subset \mathbb^d, is a sequence of pairs of sets (A_i, B_i), such that each pair is well-separated, and for each two distinct points p, q \in S, there exists ...
methods of Callahan and
Kosaraju Kosaraju Raghavayya (23 June 1905 – 27 October 1987), known mononumously by his surname Kosaraju, was an Indian lyricist and poet known for his works in Telugu cinema. He wrote about 3,000 songs in 350 films. His lyrics are steeped in Tel ...
yield optimal O(''n'' log ''n'') time per iteration with fixed dimension.


Particle mesh method

Another possibility is the particle mesh method in which space is discretised on a mesh and, for the purposes of computing the
gravitational potential In classical mechanics, the gravitational potential at a location is equal to the work (energy transferred) per unit mass that would be needed to move an object to that location from a fixed reference location. It is analogous to the electric ...
, particles are assumed to be divided between the surrounding 2x2 vertices of the mesh. Finding the potential energy Φ is easy, because the
Poisson equation Poisson's equation is an elliptic partial differential equation of broad utility in theoretical physics. For example, the solution to Poisson's equation is the potential field caused by a given electric charge or mass density distribution; with ...
:\nabla^2\Phi=4\pi G,\, where ''G'' is
Newton's constant The gravitational constant (also known as the universal gravitational constant, the Newtonian constant of gravitation, or the Cavendish gravitational constant), denoted by the capital letter , is an empirical physical constant involved in t ...
and is the density (number of particles at the mesh points), is trivial to solve by using the
fast Fourier transform A fast Fourier transform (FFT) is an algorithm that computes the discrete Fourier transform (DFT) of a sequence, or its inverse (IDFT). Fourier analysis converts a signal from its original domain (often time or space) to a representation in ...
to go to the
frequency domain In physics, electronics, control systems engineering, and statistics, the frequency domain refers to the analysis of mathematical functions or signals with respect to frequency, rather than time. Put simply, a time-domain graph shows how a s ...
where the Poisson equation has the simple form :\hat= -4\pi G\frac,\, where \vec is the comoving wavenumber and the hats denote Fourier transforms. Since \vec = -\vec\Phi, the gravitational field can now be found by multiplying by -i\vec and computing the inverse Fourier transform (or computing the inverse transform and then using some other method). Since this method is limited by the mesh size, in practice a smaller mesh or some other technique (such as combining with a tree or simple particle-particle algorithm) is used to compute the small-scale forces. Sometimes an adaptive mesh is used, in which the mesh cells are much smaller in the denser regions of the simulation.


Special-case optimizations

Several different gravitational perturbation algorithms are used to get fairly accurate estimates of the path of objects in the solar system. People often decide to put a satellite in a
frozen orbit In orbital mechanics, a frozen orbit is an orbit for an artificial satellite in which natural drifting due to the central body's shape has been minimized by careful selection of the orbital parameters. Typically, this is an orbit in which, over ...
. The path of a satellite closely orbiting the Earth can be accurately modeled starting from the 2-body elliptical orbit around the center of the Earth, and adding small corrections due to the oblateness of the Earth, gravitational attraction of the Sun and Moon, atmospheric drag, etc. It is possible to find a frozen orbit without calculating the actual path of the satellite. The path of a small planet, comet, or long-range spacecraft can often be accurately modeled starting from the 2-body elliptical orbit around the sun, and adding small corrections from the gravitational attraction of the larger planets in their known orbits. Some characteristics of the long-term paths of a system of particles can be calculated directly. The actual path of any particular particle does not need to be calculated as an intermediate step. Such characteristics include Lyapunov stability,
Lyapunov time In mathematics, the Lyapunov time is the characteristic timescale on which a dynamical system is chaotic. It is named after the Russian mathematician Aleksandr Lyapunov. It is defined as the inverse of a system's largest Lyapunov exponent. Use T ...
, various measurements from
ergodic theory Ergodic theory ( Greek: ' "work", ' "way") is a branch of mathematics that studies statistical properties of deterministic dynamical systems; it is the study of ergodicity. In this context, statistical properties means properties which are expr ...
, etc.


Two-particle systems

Although there are millions or billions of particles in typical simulations, they typically correspond to a real particle with a very large mass, typically 109
solar mass The solar mass () is a standard unit of mass in astronomy, equal to approximately . It is often used to indicate the masses of other stars, as well as stellar clusters, nebulae, galaxies and black holes. It is approximately equal to the mass ...
es. This can introduce problems with short-range interactions between the particles such as the formation of two-particle binary systems. As the particles are meant to represent large numbers of dark matter particles or groups of stars, these binaries are unphysical. To prevent this, a softened Newtonian force law is used, which does not diverge as the inverse-square radius at short distances. Most simulations implement this quite naturally by running the simulations on cells of finite size. It is important to implement the discretization procedure in such a way that particles always exert a vanishing force on themselves.


Softening

Softening is a numerical trick used in N-body techniques to prevent numerical
divergence In vector calculus, divergence is a vector operator that operates on a vector field, producing a scalar field giving the quantity of the vector field's source at each point. More technically, the divergence represents the volume density of ...
s when a particle comes too close to another (and the
force In physics, a force is an influence that can change the motion of an object. A force can cause an object with mass to change its velocity (e.g. moving from a state of rest), i.e., to accelerate. Force can also be described intuitively as a ...
goes to infinity). This is obtained by modifying the regularized
gravitational potential In classical mechanics, the gravitational potential at a location is equal to the work (energy transferred) per unit mass that would be needed to move an object to that location from a fixed reference location. It is analogous to the electric ...
of each particle as :\Phi = - \frac, (rather than 1/r) where \epsilon is the softening parameter. The value of the softening parameter should be set small enough to keep
simulation A simulation is the imitation of the operation of a real-world process or system over time. Simulations require the use of Conceptual model, models; the model represents the key characteristics or behaviors of the selected system or proc ...
s realistic.


Incorporating baryons, leptons and photons into simulations

Many simulations simulate only
cold dark matter In cosmology and physics, cold dark matter (CDM) is a hypothetical type of dark matter. According to the current standard model of cosmology, Lambda-CDM model, approximately 27% of the universe is dark matter and 68% is dark energy, with only a sm ...
, and thus include only the gravitational force. Incorporating
baryon In particle physics, a baryon is a type of composite subatomic particle which contains an odd number of valence quarks (at least 3). Baryons belong to the hadron family of particles; hadrons are composed of quarks. Baryons are also classifie ...
s,
lepton In particle physics, a lepton is an elementary particle of half-integer spin (spin ) that does not undergo strong interactions. Two main classes of leptons exist: charged leptons (also known as the electron-like leptons or muons), and neutr ...
s and
photon A photon () is an elementary particle that is a quantum of the electromagnetic field, including electromagnetic radiation such as light and radio waves, and the force carrier for the electromagnetic force. Photons are massless, so they alwa ...
s into the simulations dramatically increases their complexity and often radical simplifications of the underlying physics must be made. However, this is an extremely important area and many modern simulations are now trying to understand processes that occur during galaxy formation which could account for galaxy bias.


Computational complexity

Reif and Tate prove that if the ''n''-body reachability problem is defined as follows – given ''n'' bodies satisfying a fixed electrostatic potential law, determining if a body reaches a destination ball in a given time bound where we require a poly(''n'') bits of accuracy and the target time is poly(''n'') is in PSPACE. On the other hand, if the question is whether the body ''eventually'' reaches the destination ball, the problem is PSPACE-hard. These bounds are based on similar complexity bounds obtained for ray tracing.


Example Simulations


Common Boilerplate Code

The simplest implementation of N-body simulations where n\geq 3 is a naive propagation of orbiting bodies; naive implying that the only forces acting on the orbiting bodies is the gravitational force which they exert on each other. In
object-oriented programming Object-oriented programming (OOP) is a programming paradigm based on the concept of "objects", which can contain data and code. The data is in the form of fields (often known as attributes or ''properties''), and the code is in the form of ...
languages, such as C++, some
boilerplate code In computer programming, boilerplate code, or simply boilerplate, are sections of code that are repeated in multiple places with little to no variation. When using languages that are considered ''verbose'', the programmer must write a lot of boile ...
is useful for establishing the fundamental mathematical structures as well as data containers required for propagation; namely state vectors, and thus vectors, and some fundamental object containing this data, as well as the mass of an orbiting body. This method is applicable to other types of N-body simulations as well; a simulation of point masses with charges would use a similar method, however the force would be due to attraction or repulsion by interaction of electric fields. Regardless, acceleration of particle is a result of summed force vectors, divided by the mass of the particle: \vec=\frac\sum\vec An example of a programmatically stable and
scalable Scalability is the property of a system to handle a growing amount of work by adding resources to the system. In an economic context, a scalable business model implies that a company can increase sales given increased resources. For example, a ...
method for containing kinematic data for a particle is the use of fixed length arrays, which in optimised code allows for easy memory allocation and prediction of consumed resources; as seen in the following C++ code: struct Vector3 ; struct OrbitalEntity ; Note that contains enough room for a state vector, where: e_ = x , the projection of the objects position vector in Cartesian space along \left \;0\;0\right e_ = y , the projection of the objects position vector in Cartesian space along \left \;1\;0\right e_ = z , the projection of the objects position vector in Cartesian space along \left \;0\;1\right e_ = \dot , the projection of the objects velocity vector in Cartesian space along \left \;0\;0\right e_ = \dot , the projection of the objects velocity vector in Cartesian space along \left \;1\;0\right e_ = \dot , the projection of the objects velocity vector in Cartesian space along \left \;0\;1\right Additionally, contains enough room for a mass value.


Initialisation of Simulation Parameters

Commonly, N-body simulations will be systems based of some type of
equations of motion In physics, equations of motion are equations that describe the behavior of a physical system in terms of its motion as a function of time.''Encyclopaedia of Physics'' (second Edition), R.G. Lerner, G.L. Trigg, VHC Publishers, 1991, ISBN (V ...
; of these, most will be dependent on some initial configuration to "seed" the simulation. In systems such as those dependent on some gravitational or electric potential, the force on a simulation entity is independent on its velocity. Hence, to seed the ''forces'' of the simulation, merely initial positions are needed, but this will not allow propagation- initial velocities are required. Consider a planet orbiting a star- it has no motion, but is subject to gravitational attraction to its host star. As a time progresses, and time ''steps'' are added, it will gather velocity according to its acceleration. For a given instant in time, t_, the resultant acceleration of a body due to its neighbouring masses is independent of its velocity, however, for the time step t_, the resulting change in position is significantly different due the propagation's inherent dependency on velocity. In basic propagation mechanisms, such as the symplectic euler method to be used below, the position of an object at t_ is only dependent on its velocity at t_, as the shift in position is calculated via \vec_=\vec_+\vec_ Without acceleration, \vec_is static, however, from the perspective of an observer seeing only position, it will take two time steps to see a change in velocity. A solar-system-like simulation can be accomplished by taking average distances of planet equivalent point masses from a central star. To keep code simple, a non-rigorous approach based on semi-major axes and mean velocities will is used. Memory space for these bodies must be reserved before the bodies are configured; to allow for scalability, a
malloc C dynamic memory allocation refers to performing manual memory management for dynamic memory allocation in the C programming language via a group of functions in the C standard library, namely , , , and . The C++ programming language includes t ...
command may be used: OrbitalEntity* orbital_entities; orbital_entities = (OrbitalEntity*)malloc(sizeof(OrbitalEntity) * (9 + N_ASTEROIDS)); orbital_entities = ; // a star similar to the sun orbital_entities = ; // a planet similar to mercury orbital_entities = ; // a planet similar to venus orbital_entities = ; // a planet similar to earth orbital_entities = ; // a planet similar to mars orbital_entities = ; // a planet similar to jupiter orbital_entities = ; // a planet similar to saturn orbital_entities = ; // a planet similar to uranus orbital_entities = ; // a planet similar to neptune where is a variable which will remain at 0 temporarily, but allows for future inclusion of significant numbers of asteroids, at the users discretion. A critical step for the configuration of simulations is to establish the time ranges of the simulation, t_ to t_ , as well as the incremental time step dt which will progress the simulation forward: double t_0 = 0; double t = t_0; double dt = 86400; double t_end = 86400 * 365 * 10; // approximately a decade in seconds double BIG_G = 6.67e-11; // gravitational constant The positions and velocities established above are interpreted to be correct for t=t_ . The extent of a simulation would logically be for the period where t_\leq t .


Propagation

An entire simulation can consist of hundreds, thousands, millions, billions, or sometimes trillions of time steps. At the elementary level, each time step (for simulations with particles moving due to forces exerted on them) involves * calculating the forces on each body * calculating the accelerations of each body (\vec) * calculating the velocities of each body (\vec_=\vec_+\vec_) * calculating the new position of each body (\vec_=\vec_+\vec_) The above can be implemented quite simply with a while loop which continues while t exists in the aforementioned range: while (t < t_end) Focusing on the inner four rocky planets in the simulation, the trajectories resulting from the above propagation is shown below:


See also

* Millennium Run *
Large-scale structure of the cosmos The observable universe is a ball-shaped region of the universe comprising all matter that can be observed from Earth or its space-based telescopes and exploratory probes at the present time, because the electromagnetic radiation from these obj ...
*
GADGET A gadget is a mechanical device or any ingenious article. Gadgets are sometimes referred to as ''gizmos''. History The etymology of the word is disputed. The word first appears as reference to an 18th-century tool in glassmaking that was develop ...
*
Galaxy formation and evolution The study of galaxy formation and evolution is concerned with the processes that formed a heterogeneous universe from a homogeneous beginning, the formation of the first galaxies, the way galaxies change over time, and the processes that have gen ...
*
Natural units In physics, natural units are physical units of measurement in which only universal physical constants are used as defining constants, such that each of these constants acts as a coherent unit of a quantity. For example, the elementary charge ma ...
* Virgo Consortium *
Barnes–Hut simulation The Barnes–Hut simulation (named after Josh Barnes and Piet Hut) is an approximation algorithm for performing an ''n''-body simulation. It is notable for having order O(''n'' log ''n'') compared to a direct-sum algorithm which would b ...
* Bolshoi Cosmological Simulation


References


Further reading

* * * * * *. {{Numerical integrators Physical cosmology Gravity Simulation Cosmological simulation Articles containing video clips Computational physics