wiggle parameterization

This commit is contained in:
Eleni 2022-12-01 18:52:18 +01:00
parent cd926a0834
commit 17698dab5b
3 changed files with 45 additions and 1 deletions

View File

@ -269,6 +269,19 @@ double CosmoPower::matterTransferFunctionHu(double k)
return T_k;
}
double CosmoPower::noWiggleTk(double k)
{
double s = 44.5 * log(9.83 / OmegaEff) / (sqrt(1 + 10 * pow(OMEGA_B * h * h, 0.75)));
double alpha_Gamma = 1 - 0.328 * log(431 * OmegaEff) * OMEGA_B / OMEGA_0 + 0.38 * log(22.3 * OmegaEff) * pow(OMEGA_B / OMEGA_0, 2);
double GammaEff = OMEGA_0 * h * (alpha_Gamma + (1 - alpha_Gamma)/(1 + pow(0.43 * k * s, 4)));
double q = k/(h*GammaEff) * pow(Theta_27, 2);
double L_0 = log(2 * M_E + 1.8 * q);
double C_0 = 14.2 + 731 / (1 + 62.5 * q);
double T0 = L_0 / (L_0 + C_0 * q * q);
return T0;
}
double CosmoPower::powerHuBaryons(double k)
{
double s = 44.5 * log(9.83 / OmegaEff) / (sqrt(1 + 10 * pow(OMEGA_B * h * h, 0.75)));
@ -473,6 +486,9 @@ void CosmoPower::setFunction(CosmoFunction f)
case POWER_TEST:
eval = &CosmoPower::powerTest;
break;
case NOWIGGLE_TK:
eval = &CosmoPower::noWiggleTk;
break;
default:
abort();
}

View File

@ -90,7 +90,8 @@ namespace CosmoTool {
POWER_BDM,
POWER_TEST,
HU_WIGGLES_ORIGINAL,
SAMPLE_WIGGLES
SAMPLE_WIGGLES,
NOWIGGLE_TK
};
CosmoPower();
@ -124,6 +125,7 @@ namespace CosmoTool {
double powerTest(double k);
double powerHuWigglesOriginal(double k);
double sample_BAO(double k);
double noWiggleTk(double k);
};
};

View File

@ -58,6 +58,7 @@ using namespace std;
#define POWER_TEST 8
#define POWER_SPECTRUM POWER_EFSTATHIOU
#define SAMPLE_WIGGLES 9
namespace Cosmology {
@ -206,6 +207,31 @@ double powG(double y)
*/
double powerSpectrum(double k, double normPower)
{
#if POWER_SPECTRUM == SAMPLE_WIGGLES
// BAO wiggle parameterization for reconstruction
// Babic et al. 2022, https://arxiv.org/abs/2203.06177
// No-wiggle transfer function
double s = 44.5 * log(9.83 / OmegaEff) / (sqrt(1 + 10 * pow(OMEGA_B * h * h, 0.75)));
double alpha_Gamma = 1 - 0.328 * log(431 * OmegaEff) * OMEGA_B / OMEGA_0 + 0.38 * log(22.3 * OmegaEff) * pow(OMEGA_B / OMEGA_0, 2);
double GammaEff = OMEGA_0 * h * (alpha_Gamma + (1 - alpha_Gamma)/(1 + pow(0.43 * k * s, 4)));
double q = k/(h*GammaEff) * pow(Theta_27, 2);
double L_0 = log(2 * M_E + 1.8 * q);
double C_0 = 14.2 + 731 / (1 + 62.5 * q);
double T0 = L_0 / (L_0 + C_0 * q * q);
// Wiggle parameterization
double A = 0;
double r_s = 10;
double k_D = 2 * M_PI / 100;
double param = 1 + A * sin(k * r_s) * exp(- k / k_D);
return normPower * pow(k, n) * T0 * T0 * param;
#endif
#if POWER_SPECTRUM == POWER_EFSTATHIOU
double a = 6.4/Gamma0;
double b = 3/Gamma0;