Back to the main page Back to the documentation 1. Overview: Random fields 2. Gaussian random field 3. Local and non-local non-gaussianity 4. Gaussian white noise

We generate a gaussian random field in Fourier space by drawing the norm from a Raylight distribution and the phase from a uniform distribution. In practice this is done by drawing $A,\theta\in U(0,1)$ and setting $$\delta(k) = \sqrt{-\log(A)\frac{P(k)}{V}}e^{2\pi i \theta}$$ where the dimensionless ratio $\frac{P(k)}{V}$ is provided by the user. If the option amplitude_fixed is set then $A$ is fixed to $e^{-1}$ so $|\delta|^2 \equiv \frac{P(k)}{V}$, i.e. only the phases are random. This process is straight forward, however we must be careful to ensure that $\delta(k) = \delta(-k)^*$ (since $\delta(x)$ is real) and the Fourier conjugate to $k$, $-k$, might reside on a different task (this is solved by seeding the RNG the same way on each task based on $k$).

// Generate a gaussian random field in Fourier space
template<int N>
  void generate_gaussian_random_field_fourier(
      FFTWGrid<N> & grid,
      RandomGenerator *rng,
      std::function<double(double)> & Pofk_over_V_of_kB,
      bool fix_amplitude);

// This is just the method above + a FFT to realspace
template<int N>
  void generate_gaussian_random_field_real(
      FFTWGrid<N> & grid,
      RandomGenerator *rng,
      std::function<double(double)> & Pofk_over_V_of_kB,
      bool fix_amplitude);

The function Pofk_over_V_of_kB takes the dimensionless argument $k\cdot B$ where $B$ is the boxsize and returns the dimensionless $\frac{P(k)}{V}$. For example if $P(k) = k^3$ in 3D then

double Boxsize = 100.0; // Boxsize in whatever units you want to use
std::function<double(double)> & Pofk_over_V_of_kB = [&](double kBox){
  double k = kBox / Boxsize;
  double pofk = k * k * k;
  double volume = Boxsize * Boxsize * Boxsize;
  return pofk / volume;
};

Local non-gaussianity in its simplest form is given by $$\Phi = \phi + f_{\rm NL}(\phi^2 - \left<\phi^2\right>)$$ where $\phi$ is a gaussian random field. More generalized we can define (in Fourier space) $$\Phi = \phi + f_{\rm NL}(K(\phi,\phi) - \left<K(\phi,\phi)\right>)$$ for some non-local kernel $K$ which (in fourier space) can be written $$F[K(\phi,\phi)] \equiv \int K(k_1,k_2)\phi(k_1)\phi(k_2)\delta(k - k_1-k_2)d^3k_1d^3k_2$$ for some Fourier space kernel $K(k_1,k_2)$. We provide some methods for computing such fields. Defining $$[\mathcal{P}_mA](x) = \int e^{-ikx}\left(\frac{P(k)}{V}\right)^m A(k)d^3k$$ where $P(k)$ is the power-spectrum and $V$ is the volume of the box (just to have correct units, but this factors out below) then we can evaluate $$K[\phi,\phi] = A\phi^2 + B \mathcal{P}_{2/3}( [\mathcal{P}_{-1/3}\phi]^2 ) + C \mathcal{P}_{1/3}( \phi[\mathcal{P}_{-1/3}\phi] ) + D \mathcal{P}_{2/3}( \phi[\mathcal{P}_{-2/3}\phi] ) + E \mathcal{P}_{3/3}( \phi[\mathcal{P}_{-3/3}\phi] ) + F \mathcal{P}_{3/3}( [\mathcal{P}_{-2/3}\phi][\mathcal{P}_{-1/3}\phi] )$$ for a given set of coefficients $A,B,C,D,E,F$. Specific choices for these coefficients gives us local, equilateral and orthogonal configurations (the three options of type_of_fnl below). If you want to specify your own values than use type_of_fnl=generic and provide them in kernel_values (in the order $A,C,D,E$ and as currently written we assume $B=-D$ and $F=-E$).

template<int N>
  void generate_nonlocal_gaussian_random_field_fourier(
      FFTWGrid<N> & phi_fourier,
      RandomGenerator *rng,
      std::function<double(double)> & DeltaPofk,
      bool fix_amplitude,
      double fNL,
      std::string type_of_fnl,
      double u = 0.0,
      std::vector<double> kernel_values = {});

For more info about this see Scoccimarro et al. 2012 and Schmidt and Kamionkowski 2013

We can generate a white noise field in real space. It is normalized such that the power-spectrum of this field is unity.

template<int N>
  void generate_white_noise_field_real(
      FFTWGrid<N> & grid,
      RandomGenerator * rng);