Open Lagrang3 opened 2 years ago
As seen in the picture. The power spectrum of T00 and k^2 Phi do not match at the small and again at the large scales. Why?
This problem is obtained at this commit 716f5d25ae1, running Gadget.
First idea: produce a minimal working gevolution main, that loads a snapshot (or produces one), computes the potentials, prints the power spectrum and exits.
Next. Let us find the exact conversion factor from T00 to Phi.
The first key functions to look at are relativistic_pm::sample
and relativistic_pm::compute_potential
Question: if $\tilde A(\vec k) = \tilde B (\vec k) / (|\vec k|)$, then what is the relation between $P(\tilde A, k)$ and $P(\tilde B,k)$?
In gevolution, the T00 component of the energy momentum tensor is substracted the background value of T00 and multiplied by a factor fourpiG dx^2 / a, where fourpiG = 4 pi G in gevolution units, which is numerically equal to 1.5 L^2/c^2, where L is the boxsize in Mpc/h and c is the speed of light in km/s. dx = L/N, where N is the number of grid points per dimension. Out of the previous operation gevolution obtains a new field S00 which becomes the source of the potential Phi. In the Poisson equation solver Phi = 1/k^2 S00 / N^3 / dx^2, where k = 2 pi n being n the integer vector identifying the Fourier mode.
The previous power spectrum corresponds to Pk(T00) and const k^4 Pk(Phi), and they don't match. The proportionality constant is correct, but their shape it just isn't the same. The issue is that Pk( f(k) A ) is different from f(k)^2 Pk(A).
In the PowerI4 documentation by Emiliano Sefussati, we find that the power spectrum of a field A is defined as P_k(A) = kf^3 / |\Omegak| \sum{q \in \Omega_k} |A_k|^2 where kf = 2 pi/L is the fundamental mode, \Omega_k is a set of modes such that q\in \Omega_k if and only if k-dk/2 <= |q| < k+dk/2, and dk is a bin size. With that in mind, given a mode q we can find which set \Omega_k it belongs: k/dk = \floor(1/2 + |q|/dk) where |q| =\sqrt(q_1^2 + q_2^2 + q_3^3)
Then consider f(|q|) a function of the modulo of the mode P_k(f A) = kf^3 / |\Omegak| \sum{q\in \Omega_k} |A_q|^2 |f(|q|)|^2 if |q| had the same value for every q \in \Omega_k, for instance when dk is very small, q\in\Omega_k implies |q|=k, then P_k(f A) = kf^3 / |\Omegak| |f(k)|^2 \sum{q\in\Omega_k} |A_q|^2.
However, if we have dk = kf, we are actually stating that the bins are shells thick enough as the distance between two consecutive nodes, and yet we will show that for the different values of q\in\Omega_k there will be non unique values of |q|.
Let us count for example some modes and their bins for dk = kf | q/kf | degeneracy | q^2/kf^2 | k/kf |
---|---|---|---|---|
(0,0,0) | 1 | 0 | 0 | |
(0,0,1) | 3 | 1 | 1 | |
(0,1,1) | 6 | 2 | 1 | |
(1,1,1) | 4 | 3 | 2 | |
(0,0,2) | 3 | 4 | 2 | |
(0,1,2) | 12 | 5 | 2 | |
(1,1,2) | 12 | 6 | 2 |
so several values of q^2 will correspond to the same value of k, and we cannot have P_k(f A) = f(k)^2 P_k(A)
One possible definition of the power spectrum that is bilinear could be the following, P_k(A,B) = 1/|Nk| \sum{q \in N_k} A*(q) B(q) where N_k = { q | q^2 = k^2 },
then P_k(A f, B g) = f*(k) g(k) P_k(A,B)
with this new definition, we obtain a more wigling power spectrum because there are less points on each bins, but P_k is bilinear and the power spectra of the source and the potential are related by a factor k^4
Smoothing the power spectra with bins of size 1 we obtain. This is a gevolution (1 step) simulation L=500 Mpc/h, N = 256^3 particles.
My guess is that with a bigger boxsize we will be able to start noticing deviations between T00 and k^2 Phi at large scales due to the Hubble terms in the Einstein equation for Phi. For example, for L=2000 Mpc/h, N=256^3 we find the following
pw-phi-plus.pdf