guilhem_lavaux/cosmotool
0
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Merged in AdamAndrewsSU/cosmotool/AA_Primordial (pull request #1)

Browse Source 
additions from PNG-project
This commit is contained in:
Adam Andrews 2019-07-11 19:16:33 +00:00 committed by Guilhem Lavaux
commit df29233852
2 changed files with 84 additions and 38 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

44
src/cosmopower.cpp
View File

@ -212,6 +212,43 @@ double CosmoPower::powerHuWiggles(double k)
return normPower * pow(k,n) * T_k * T_k;
}
double CosmoPower::primordialPowerSpectrum(double k)
{
//Primordial power spectrum, needed for PNG
return normPower * pow(k,n);
}
double CosmoPower::matterTransferFunctionHu(double k)
{
// EISENSTEIN ET HU (1998)
// FULL POWER SPECTRUM WITH BARYONS AND WIGGLES
double k_silk = ehu.k_silk;
double s = ehu.s;
double k_eq = ehu.k_eq;
double alpha_c = ehu.alpha_c;
double beta_c = ehu.beta_c;
double alpha_b = ehu.alpha_b;
double beta_b = ehu.beta_b;
double xx = k * s;
double s_tilde = s * pow(1 + pow(ehu.beta_node / (xx), 3), -1./3);
double f = 1 / (1 + pow(xx / 5.4, 4));
double q = k / (13.41 * k_eq);
double T_c = f * T_tilde_0(q, 1, beta_c) + (1 - f) * T_tilde_0(q, alpha_c, beta_c);
double T_b = (
T_tilde_0(q, 1, 1) / (1 + pow(xx / 5.2, 2)) +
alpha_b / (1 + pow(beta_b / xx, 3)) * exp(-pow(k/k_silk, 1.4))
) * j_0(k * s_tilde);
double T_k = OMEGA_B/OMEGA_0 * T_b + OMEGA_C/OMEGA_0 * T_c;
return T_k;
}
double CosmoPower::powerHuBaryons(double k)
{
double s = 44.5 * log(9.83 / OmegaEff) / (sqrt(1 + 10 * pow(OMEGA_B * h * h, 0.75)));
@ -383,6 +420,12 @@ void CosmoPower::setFunction(CosmoFunction f)
updateHuWigglesConsts();
eval = &CosmoPower::powerHuWiggles;
break;
case PRIMORDIAL_PS:
eval = &CosmoPower::primordialPowerSpectrum;
break;
case MATTER_TK:
eval = &CosmoPower::matterTransferFunctionHu;
break;
case HU_WIGGLES_ORIGINAL:
updateHuWigglesOriginal();
eval = &CosmoPower::powerHuWigglesOriginal;
@ -414,4 +457,3 @@ void CosmoPower::setNormalization(double A_K)
{
normPower = A_K;///power(0.002);
}

4
src/cosmopower.hpp
View File

@ -81,6 +81,8 @@ namespace CosmoTool {
enum CosmoFunction {
POWER_EFSTATHIOU,
HU_WIGGLES,
PRIMORDIAL_PS,
MATTER_TK,
HU_BARYON,
OLD_POWERSPECTRUM,
POWER_BARDEEN,
@ -111,6 +113,8 @@ namespace CosmoTool {
double powerEfstathiou(double k);
double powerHuWiggles(double k);
double primordialPowerSpectrum(double k);
double matterTransferFunctionHu(double k);
double powerHuBaryons(double k);
double powerOld(double k);
double powerBardeen(double k);