FML
|
This namespace deals with N-body simulations. More...
Typedefs | |
template<int N> | |
using | FFTWGrid = FML::GRID::FFTWGrid< N > |
template<class T > | |
using | MPIParticles = FML::PARTICLE::MPIParticles< T > |
Functions | |
template<int N, class T > | |
void | DriftParticles (FML::PARTICLE::MPIParticles< T > &part, double delta_time, bool periodic_box=true) |
template<int N, class T > | |
void | DriftParticles (T *p, size_t NumPart, double delta_time, bool periodic_box) |
This moves the particles according to \( x_{\rm new} = x + v \Delta t \). More... | |
template<int N, class T > | |
void | KickParticles (std::array< FFTWGrid< N >, N > &force_grid, MPIParticles< T > &part, double delta_time, std::string interpolation_method) |
This moves the particle velocities according to \( v_{\rm new} = v + F \Delta t \). More... | |
template<int N, class T > | |
void | KickParticles (std::array< FFTWGrid< N >, N > &force_grid, T *p, size_t NumPart, double delta_time, std::string interpolation_method) |
This moves the particle velocities according to \( v_{\rm new} = v + F \Delta t \). More... | |
template<int N> | |
void | compute_force_from_density_real (const FFTWGrid< N > &density_grid_real, std::array< FFTWGrid< N >, N > &force_real, std::string density_assignment_method_used, double norm_poisson_equation) |
Take a density grid in real space and returns the force \( \nabla \phi \) where \( \nabla^2 \phi = {\rm norm} \cdot \delta \) Different choices for what kernel to use for \( \nabla / \nabla^2\) are availiable, see the function body (is set too be a compile time option). More... | |
template<int N> | |
void | compute_force_from_density_fourier (const FFTWGrid< N > &density_grid_fourier, std::array< FFTWGrid< N >, N > &force_real, std::string density_assignment_method_used, double norm_poisson_equation) |
Take a density grid in fourier space and returns the force \( \nabla \phi \) where \( \nabla^2 \phi = {\rm norm} \cdot \delta \) Different choices for what kernel to use for \( \nabla / \nabla^2\) are availiable, see the function body (is set too be a compile time option). More... | |
template<int N, class T > | |
void | KickDriftKickNBodyStep (int Nmesh, MPIParticles< T > &part, double delta_time, std::string density_assignment_method, double norm_poisson_equation) |
Take a N-body step with a simple Kick-Drift-Kick method (this method serves mainly as an example for how one can do this). More... | |
template<int N, class T > | |
void | YoshidaNBodyStep (int Nmesh, MPIParticles< T > &part, double delta_time, std::string density_assignment_method, double norm_poisson_equation) |
Take a N-body step with a 4th order symplectic Yoshida method. More... | |
template<int N, class T > | |
void | DriftParticles (MPIParticles< T > &part, double delta_time, bool periodic_box) |
This moves the particles according to \( x_{\rm new} = x + v \Delta t \). More... | |
template<int N, class T > | |
void | NBodyInitialConditions (MPIParticles< T > &part, int Npart_1D, double buffer_factor, const FFTWGrid< N > &delta_fourier, std::vector< FFTWGrid< N > > &phi_nLPT_potentials, int LPT_order, double box, double zini, std::vector< double > velocity_norms) |
Generate particles from a given initial density field using Lagrangian perturbation theory. More... | |
template<int N, class T > | |
void | NBodyInitialConditions (MPIParticles< T > &part, int Npart_1D, double buffer_factor, int Nmesh, bool fix_amplitude, FML::RANDOM::RandomGenerator *rng, std::function< double(double)> Pofk_of_kBox_over_Pofk_primordal, std::function< double(double)> Pofk_of_kBox_over_volume_primordial, int LPT_order, std::string type_of_random_field, double fNL, double box, double zini, std::vector< double > velocity_norms) |
Generate particles from a given power-spectrum using Lagrangian perturbation theory. More... | |
template<int N> | |
void | compute_delta_fifth_force (const FFTWGrid< N > &density_fourier, FFTWGrid< N > &density_mg_fourier, std::function< double(double)> coupling_factor_of_kBox) |
This method computes the fifth-force potential for modified gravity models using the linear approximation This computes \( \delta_{\rm MG}(k) \) where the total force in fourier space is \( F(k) \propto \frac{\vec{k}}{k^2}[\delta(k) + \delta_{\rm MG}(k)] \) by solving \( \nabla^2 \phi = m^2 \phi + F^{-1}[g(k) \delta(k)] \) where \( \delta_{\rm MG}(k) = -k^2\phi(k) \). More... | |
template<int N> | |
void | compute_delta_fifth_force_potential_screening (const FFTWGrid< N > &density_fourier, FFTWGrid< N > &density_mg_fourier, std::function< double(double)> coupling_factor_of_kBox, std::function< double(double)> screening_factor_of_newtonian_potential, double poisson_norm) |
This method computes the fifth-force potential for modified gravity models which has a screening mechanism using the approximate method of Winther & Ferreira 2015. More... | |
template<int N> | |
void | compute_delta_fifth_force_density_screening (const FFTWGrid< N > &density_fourier, FFTWGrid< N > &density_mg_fourier, std::function< double(double)> coupling_factor_of_kBox, std::function< double(double)> screening_factor_of_density, double smoothing_scale, std::string smoothing_method) |
This method computes the fifth-force potential for modified gravity models which has a screening mechanism using the approximate method of Winther & Ferreira 2015. More... | |
double | compute_sigma_of_R (std::function< double(double)> Pofk_of_kBox_over_volume, double R_mpch, double boxsize_mpch) |
This computes the standard deviation of linear fluctuations smoothed over a sphere of radius, \( \sigma(R) \), by integrating the power-spectrum convolved with a top-hat windowfunction of radius \( R \). More... | |
This namespace deals with N-body simulations.
Computing forces, moving particles and generating initial conditions.
using FML::NBODY::FFTWGrid = typedef FML::GRID::FFTWGrid<N> |
using FML::NBODY::MPIParticles = typedef FML::PARTICLE::MPIParticles<T> |
void FML::NBODY::compute_delta_fifth_force | ( | const FFTWGrid< N > & | density_fourier, |
FFTWGrid< N > & | density_mg_fourier, | ||
std::function< double(double)> | coupling_factor_of_kBox | ||
) |
This method computes the fifth-force potential for modified gravity models using the linear approximation This computes \( \delta_{\rm MG}(k) \) where the total force in fourier space is \( F(k) \propto \frac{\vec{k}}{k^2}[\delta(k) + \delta_{\rm MG}(k)] \) by solving \( \nabla^2 \phi = m^2 \phi + F^{-1}[g(k) \delta(k)] \) where \( \delta_{\rm MG}(k) = -k^2\phi(k) \).
For example in \( f(R) \) gravity we have \( g(k) = \frac{1}{3}\frac{k^2}{k^2 + m^2}\) and in DGP we have \( g(k) = \frac{1}{3\beta} \) (independent of scale).
N | The dimension we work in. |
[in] | density_fourier | The density contrast in fourier space. |
[out] | density_mg_fourier | The force potential. |
[in] | coupling_factor_of_kBox | The coupling factor \( g(k) \) |
Definition at line 1072 of file NBody.h.
void FML::NBODY::compute_delta_fifth_force_density_screening | ( | const FFTWGrid< N > & | density_fourier, |
FFTWGrid< N > & | density_mg_fourier, | ||
std::function< double(double)> | coupling_factor_of_kBox, | ||
std::function< double(double)> | screening_factor_of_density, | ||
double | smoothing_scale, | ||
std::string | smoothing_method | ||
) |
This method computes the fifth-force potential for modified gravity models which has a screening mechanism using the approximate method of Winther & Ferreira 2015.
This computes \( \delta_{\rm MG}(k) \) where the total force in fourier space is \( F(k) \propto \frac{\vec{k}}{k^2}[\delta(k) + \delta_{\rm MG}(k)] \) by solving \( \nabla^2 \phi = m^2 \phi + f(\rho)F^{-1}[g(k) \delta(k)] \) where \( \delta_{\rm MG}(k) = -k^2\phi(k) \) and \( \rho \) is the density in units of the mean density. For example in DGP gravity we have \( g(k) = \frac{1}{3\beta} \) (independent of scale) and the screening function is \( f(\rho) = 2\frac{\sqrt{1 + C} - 1}{C} \) where \( C = \frac{8\Omega_M(r_cH_0)^2}{9\beta^2}\rho \) If you don't want screening then simply pass the function \( f \equiv 1 \) and the equation reduces to the one in the linear regime.
N | The dimension we work in. |
[in] | density_fourier | The density contrast in fourier space. |
[out] | density_mg_fourier | The force potential. |
[in] | coupling_factor_of_kBox | The coupling factor \( g(k) \) |
[in] | screening_factor_of_density | The screening factor \( f(\rho) \) Should be in \( [0,1] \) and go to 1 for \( \rho \to 0 \) and to 0 for \( \rho \to \infty \). |
[in] | smoothing_scale | The smoothing radius in units of the boxsize. |
[in] | smoothing_method | The k-space smoothing filter (gaussian, tophat, sharpk). |
Definition at line 1224 of file NBody.h.
void FML::NBODY::compute_delta_fifth_force_potential_screening | ( | const FFTWGrid< N > & | density_fourier, |
FFTWGrid< N > & | density_mg_fourier, | ||
std::function< double(double)> | coupling_factor_of_kBox, | ||
std::function< double(double)> | screening_factor_of_newtonian_potential, | ||
double | poisson_norm | ||
) |
This method computes the fifth-force potential for modified gravity models which has a screening mechanism using the approximate method of Winther & Ferreira 2015.
This computes \( \delta_{\rm MG}(k) \) where the total force in fourier is given by \( F(k) \propto \frac{\vec{k}}{k^2}[\delta(k) + \delta_{\rm MG}(k)] \) by solving \( \nabla^2 \phi = m^2 \phi + f(\Phi)F^{-1}[g(k) \delta(k)] \) where \( \delta_{\rm MG}(k) = -k^2\phi(k) \) For example in \( f(R) \) gravity we have \( g(k) = \frac{1}{3}\frac{k^2}{k^2 + m^2} \) and the screening function is \( f(\Phi) = \min(1, \left|\frac{3f_R}{2\Phi}\right|) \) If you don't want screening then simpy pass the function \( f \equiv 1 \) and the equation reduces to the one in the linear regime
N | The dimension we work in. |
[in] | density_fourier | The density contrast in fourier space. |
[out] | density_mg_fourier | The force potential. |
[in] | coupling_factor_of_kBox | The coupling factor \(g(k)\) |
[in] | screening_factor_of_newtonian_potential | The screening factor \( f(\Phi_N) \) Should be in \( [0,1] \) and go to 1 for \( \Phi_N \to 0 \) and 0 for very large \( \Phi_N \) |
[in] | poisson_norm | The factor \( C \) in \( \nabla^2\Phi = C\delta \) to get the potential in the metric (not the code-potential) so \( C = \frac{3}{2}\Omega_M a \frac{(H_0B)^2}{a^2} \) |
Definition at line 1122 of file NBody.h.
void FML::NBODY::compute_force_from_density_fourier | ( | const FFTWGrid< N > & | density_grid_fourier, |
std::array< FFTWGrid< N >, N > & | force_real, | ||
std::string | density_assignment_method_used, | ||
double | norm_poisson_equation | ||
) |
Take a density grid in fourier space and returns the force \( \nabla \phi \) where \( \nabla^2 \phi = {\rm norm} \cdot \delta \) Different choices for what kernel to use for \( \nabla / \nabla^2\) are availiable, see the function body (is set too be a compile time option).
Fiducial choice is the continuous greens function \( 1/k^2\), but we can also choose to also devonvolve the window and discrete kernels (Hamming 1989; same as used in GADGET) and Hockney & Eastwood 1988. See e.g. 1603.00476 for a list and references.
N | The dimension of the grid |
[in] | density_grid_fourier | The density contrast in fourier space. |
[out] | force_real | The force in real space. |
[in] | density_assignment_method_used | The density assignement we used to compute the density field. Needed only in case kernel_choice (defined in the body of this function) is not CONTINUOUS_GREENS_FUNCTION. |
[in] | norm_poisson_equation | The prefactor (norm) to the Poisson equation. |
Definition at line 253 of file NBody.h.
void FML::NBODY::compute_force_from_density_real | ( | const FFTWGrid< N > & | density_grid_real, |
std::array< FFTWGrid< N >, N > & | force_real, | ||
std::string | density_assignment_method_used, | ||
double | norm_poisson_equation | ||
) |
Take a density grid in real space and returns the force \( \nabla \phi \) where \( \nabla^2 \phi = {\rm norm} \cdot \delta \) Different choices for what kernel to use for \( \nabla / \nabla^2\) are availiable, see the function body (is set too be a compile time option).
Fiducial choice is the continuous greens function \( 1/k^2\), but we can also choose to also devonvolve the window and discrete kernels (Hamming 1989; same as used in GADGET) and Hockney & Eastwood 1988. See e.g. 1603.00476 for a list.
N | The dimension of the grid |
[in] | density_grid_real | The density contrast in real space. |
[out] | force_real | The force in real space. |
[in] | density_assignment_method_used | The density assignement we used to compute the density field. Needed only in case kernel_choice (defined in the body of this function) is not CONTINUOUS_GREENS_FUNCTION. |
[in] | norm_poisson_equation | The prefactor (norm) to the Poisson equation. |
Definition at line 222 of file NBody.h.
double FML::NBODY::compute_sigma_of_R | ( | std::function< double(double)> | Pofk_of_kBox_over_volume, |
double | R_mpch, | ||
double | boxsize_mpch | ||
) |
This computes the standard deviation of linear fluctuations smoothed over a sphere of radius, \( \sigma(R) \), by integrating the power-spectrum convolved with a top-hat windowfunction of radius \( R \).
We do this by solving \( \sigma^2(R) = \int \frac{k^3P(k)}{2\pi^2} |W(kR)|^2 d\log k \)
[in] | Pofk_of_kBox_over_volume | Dimensionless power-spectrum \( P / V \) as function of the dimensionless scale \( kB \) where \(B\) is the boxsize and \(V = B^3\) is the box volume. |
[in] | R_mpch | R in units of Mpc/h |
[in] | boxsize_mpch | Boxsize in units of Mpc/h |
Definition at line 1294 of file NBody.h.
void FML::NBODY::DriftParticles | ( | FML::PARTICLE::MPIParticles< T > & | part, |
double | delta_time, | ||
bool | periodic_box = true |
||
) |
void FML::NBODY::DriftParticles | ( | MPIParticles< T > & | part, |
double | delta_time, | ||
bool | periodic_box | ||
) |
This moves the particles according to \( x_{\rm new} = x + v \Delta t \).
Note that we assume the velocities are in such units that \( v \Delta t\) is a dimensionless shift in [0,1).
N | The dimension of the grid |
T | The particle class |
[out] | part | MPIParticles containing the particles. |
[in] | delta_time | The size of the timestep. |
[in] | periodic_box | Is the box periodic? |
Definition at line 395 of file NBody.h.
void FML::NBODY::DriftParticles | ( | T * | p, |
size_t | NumPart, | ||
double | delta_time, | ||
bool | periodic_box | ||
) |
This moves the particles according to \( x_{\rm new} = x + v \Delta t \).
Note that we assume the velocities are in such units that \( v \Delta t\) is a dimensionless shift in [0,1). NB: after this methods is done the particles might have left the current task and must be communicated (this is done automatically if you use the MPIParticles version of this method).
N | The dimension of the grid |
T | The particle class |
[out] | p | Pointer to the first particle. |
[in] | NumPart | The number of local particles. |
[in] | delta_time | The size of the timestep. |
[in] | periodic_box | Is the box periodic? |
Definition at line 427 of file NBody.h.
void FML::NBODY::KickDriftKickNBodyStep | ( | int | Nmesh, |
MPIParticles< T > & | part, | ||
double | delta_time, | ||
std::string | density_assignment_method, | ||
double | norm_poisson_equation | ||
) |
Take a N-body step with a simple Kick-Drift-Kick method (this method serves mainly as an example for how one can do this).
N | The dimension of the grid. |
T | The particle class. |
[in] | Nmesh | The gridsize to use for computing the density and force. |
[out] | part | The particles |
[in] | delta_time | The time \( \Delta t \) we move forward. |
[in] | density_assignment_method | The density assignement method (NGP, CIC, TSC, PCS or PQS). |
[in] | norm_poisson_equation | A possible prefactor to the Poisson equation |
Definition at line 88 of file NBody.h.
void FML::NBODY::KickParticles | ( | std::array< FFTWGrid< N >, N > & | force_grid, |
MPIParticles< T > & | part, | ||
double | delta_time, | ||
std::string | interpolation_method | ||
) |
This moves the particle velocities according to \( v_{\rm new} = v + F \Delta t \).
This method assumes the force is normalized such that \( F \Delta t \) has the same units as your v. If the flag free_force_grids is set in the source then we free up memory of the force grids after we have used them. The defalt is false.
N | The dimension of the grid |
T | The particle class |
[in] | force_grid | Grid containing the force. |
[out] | part | MPIParticles containing the particles. |
[in] | delta_time | The size of the timestep. |
[in] | interpolation_method | The interpolation method for interpolating the force to the particle positions. |
Definition at line 485 of file NBody.h.
void FML::NBODY::KickParticles | ( | std::array< FFTWGrid< N >, N > & | force_grid, |
T * | p, | ||
size_t | NumPart, | ||
double | delta_time, | ||
std::string | interpolation_method | ||
) |
This moves the particle velocities according to \( v_{\rm new} = v + F \Delta t \).
This method assumes the force is normalized such that \( F \Delta t \) has the same units as your v. If the flag free_force_grids is set in the source then we free up memory of the force grids after we have used them. The defalt is false.
N | The dimension of the grid |
T | The particle class |
[in] | force_grid | The force \( \nabla \Phi \). |
[out] | p | Pointer to the first particle. |
[in] | NumPart | The number of local particles. |
[in] | delta_time | The size of the timestep. |
[in] | interpolation_method | The interpolation method for interpolating the force to the particle positions. |
Definition at line 513 of file NBody.h.
void FML::NBODY::NBodyInitialConditions | ( | MPIParticles< T > & | part, |
int | Npart_1D, | ||
double | buffer_factor, | ||
const FFTWGrid< N > & | delta_fourier, | ||
std::vector< FFTWGrid< N > > & | phi_nLPT_potentials, | ||
int | LPT_order, | ||
double | box, | ||
double | zini, | ||
std::vector< double > | velocity_norms | ||
) |
Generate particles from a given initial density field using Lagrangian perturbation theory.
We generate particles in [0,1) and velocities are given by \( v_{\rm code} = \frac{a^2 \frac{dx}{dt}}{H_0 B} \)
N | The dimension we are working in. |
T | The particle class. Must have methods get_pos, get_vel, get_D_1LPT and get_D_2LPT. But only get_pos data is required to exist. Return a nullptr if the data does not exist in the particle. |
[out] | part | Particle container for particles we are to create. |
[in] | Npart_1D | Number of particles per dimension (i.e. total is \( {\rm Npart}_{\rm 1D}^N \)) |
[in] | buffer_factor | How many more particles to allocate? |
[in] | delta_fourier | The initial density field \( \delta(k,z_{\rm ini})\) in fourier space |
[in] | phi_nLPT_potentials | Return the LPT potentials: 2LPT, 3LPTa, 3LPTb, ... If the vector has zero size then nothing will be returned. |
[in] | LPT_order | The LPT order (1 or 2) |
[in] | box | The boxsize (only for prining maximum displacement) |
[in] | zini | The initial redshift |
[in] | velocity_norms | A vector of the factors we need to multiply the nLPT displacement fields by to get velocities. E.g. \( 100 {\rm Box_in_Mpch} f_i(z_{\rm ini}) H(z_{\rm ini})/H_0 \cdot a_{\rm ini} \) to get peculiar velocities in km/s and \( f_i(z_{\rm ini}) H(z_{\rm ini})/H_0 \cdot a_{\rm ini}^2 \) to get the velocities we use as the fiducial choice in N-body. |
Definition at line 694 of file NBody.h.
void FML::NBODY::NBodyInitialConditions | ( | MPIParticles< T > & | part, |
int | Npart_1D, | ||
double | buffer_factor, | ||
int | Nmesh, | ||
bool | fix_amplitude, | ||
FML::RANDOM::RandomGenerator * | rng, | ||
std::function< double(double)> | Pofk_of_kBox_over_Pofk_primordal, | ||
std::function< double(double)> | Pofk_of_kBox_over_volume_primordial, | ||
int | LPT_order, | ||
std::string | type_of_random_field, | ||
double | fNL, | ||
double | box, | ||
double | zini, | ||
std::vector< double > | velocity_norms | ||
) |
Generate particles from a given power-spectrum using Lagrangian perturbation theory.
We generate particles in [0,1) and velocities are given by \( v_{\rm code} = \frac{a^2 \frac{dx}{dt}}{H_0 B} \)
N | The dimension we are working in. |
T | The particle class. Must have methods get_pos, get_vel, get_D_1LPT and get_D_2LPT. But only get_pos data is required to exist. Return a nullptr if the data does not exist in the particle. |
[out] | part | Particle container for particles we are to create. |
[in] | Npart_1D | Number of particles per dimension (i.e. total is \( {\rm Npart}_{\rm 1D}^N \)) |
[in] | buffer_factor | How many more particles to allocate? |
[in] | Nmesh | The grid to generate the IC on |
[in] | fix_amplitude | Amplitude fixed? Only random phases if true. |
[in] | rng | Random number generator |
[in] | Pofk_of_kBox_over_Pofk_primordal | The ratio of the power-spectrum (for delta) at the time you want the density field to be created at to the primordial one (the function above). |
[in] | Pofk_of_kBox_over_volume_primordial | The dimensionless function \( P/V\) where \( V = B^N\) is the box volume as function of the dimensionless wavenumber \( kB \) where \( B \) is the boxsize and \( P(k) \) is the primordial power-spectrum for \(\Phi\). |
[in] | LPT_order | The LPT order (1 or 2) |
[in] | type_of_random_field | What random field: gaussian, local, equilateral, orthogonal |
[in] | fNL | If non-gaussianity the value of fNL |
[in] | box | The boxsize (only for prining maximum displacement) |
[in] | zini | The initial redshift |
[in] | velocity_norms | A vector of the factors we need to multiply the nLPT displacement fields by to get velocities. E.g. \( 100 {\rm Box_in_Mpch} f_i(z_{\rm ini}) H(z_{\rm ini})/H_0 \cdot a_{\rm ini} \) to get peculiar velocities in km/s and \( f_i(z_{\rm ini}) H(z_{\rm ini})/H_0 \cdot a_{\rm ini}^2 \) to get the velocities we use as the fiducial choice in N-body. The order is: 1LPT, 2LPT, 3LPTa, 3LPTb |
Definition at line 609 of file NBody.h.
void FML::NBODY::YoshidaNBodyStep | ( | int | Nmesh, |
MPIParticles< T > & | part, | ||
double | delta_time, | ||
std::string | density_assignment_method, | ||
double | norm_poisson_equation | ||
) |
Take a N-body step with a 4th order symplectic Yoshida method.
This method is mainly an illustration, for using this with cosmology we should take into account that norm_poisson_equation is a function of time
N | The dimension of the grid. |
T | The particle class. |
[in] | Nmesh | The gridsize to use for computing the density and force. |
[out] | part | The particles |
[in] | delta_time | The time \( \Delta t \) we move forward. |
[in] | density_assignment_method | The density assignement method (NGP, CIC, TSC, PCS or PQS). |
[in] | norm_poisson_equation | A possible prefactor to the Poisson equation |
Definition at line 149 of file NBody.h.