mirror/monofonIC
1
0
Fork 0
mirror of https://github.com/cosmo-sims/monofonIC.git synced 2024-09-19 17:03:45 +02:00
Code Releases Activity

working commit, decaying mode integration for single fluid

Browse source
This commit is contained in:
Oliver Hahn 2019-11-27 22:19:39 +01:00
parent a5253bcace
commit 093363791e
2 changed files with 67 additions and 10 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

7
example.conf
View file

@ -4,7 +4,7 @@ GridRes = 128
# length of the box in Mpc/h
BoxLength = 250
# starting redshift
zstart = 49.0
zstart = 100.0
# order of the LPT to be used (1,2 or 3)
LPTorder = 3
# also do baryon ICs?
@ -15,8 +15,9 @@ DoFixing = no
ParticleLoad = sc
[cosmology]
#transfer = CLASS
transfer = eisenstein
transfer = CLASS
ztarget = 100.0
#transfer = eisenstein
Omega_m = 0.302
Omega_b = 0.045
Omega_L = 0.698

70
src/plugins/transfer_CLASS.cc
View file

@ -25,7 +25,13 @@ private:
std::vector<double> tab_lnk_, tab_dtot_, tab_dc_, tab_db_, tab_ttot_, tab_tc_, tab_tb_;
gsl_interp_accel *gsl_ia_dtot_, *gsl_ia_dc_, *gsl_ia_db_, *gsl_ia_ttot_, *gsl_ia_tc_, *gsl_ia_tb_;
gsl_spline *gsl_sp_dtot_, *gsl_sp_dc_, *gsl_sp_db_, *gsl_sp_ttot_, *gsl_sp_tc_, *gsl_sp_tb_;
double Omega_m_, Omega_b_, N_ur_, zstart_, ztarget_, kmax_, kmin_, h_;
// single fluid growing/decaying mode decomposition
gsl_interp_accel *gsl_ia_Cplus_, *gsl_ia_Cminus_;
gsl_spline *gsl_sp_Cplus_, *gsl_sp_Cminus_;
std::vector<double> tab_Cplus_, tab_Cminus_;
double Omega_m_, Omega_b_, N_ur_, zstart_, ztarget_, kmax_, kmin_, h_, astart_, atarget_;
void ClassEngine_get_data( void ){
std::vector<double> d_ncdm, t_ncdm, phi, psi;
@ -37,25 +43,46 @@ private:
pars.add("extra metric transfer functions", "yes");
pars.add("z_pk",ztarget_);
pars.add("P_k_max_h/Mpc", kmax_);
pars.add("h",h_);
pars.add("Omega_b",Omega_b_);
// pars.add("Omega_k",0.0);
// pars.add("Omega_ur",0.0);
pars.add("N_ur",N_ur_);
pars.add("Omega_cdm",Omega_m_-Omega_b_);
pars.add("Omega_Lambda",1.0-Omega_m_);
// pars.add("Omega_fld",0.0);
// pars.add("Omega_scf",0.0);
pars.add("Omega_k",0.0);
// pars.add("Omega_Lambda",1.0-Omega_m_);
pars.add("Omega_fld",0.0);
pars.add("Omega_scf",0.0);
pars.add("A_s",2.42e-9);
pars.add("n_s",.96); // tnis doesn't matter for TF
pars.add("n_s",.961); // this doesn't matter for TF
pars.add("output","dTk,vTk");
pars.add("YHe",0.248);
pars.add("lensing","no");
pars.add("alpha_s",0.0);
pars.add("P_k_ini type","analytic_Pk");
pars.add("gauge","synchronous");
pars.add("k_per_decade_for_pk",100);
pars.add("k_per_decade_for_bao",100);
pars.add("k_per_decade_for_pk",50);
pars.add("k_per_decade_for_bao",50);
pars.add("compute damping scale","yes");
pars.add("z_reio",-1.0); // make sure reionisation is not included
pars.add("tol_perturb_integration",1.e-8);
pars.add("tol_background_integration",1e-9);
// high precision options from cl_permille.pre:
// precision file to be passed as input in order to achieve at least percent precision on scalar Cls
pars.add("hyper_flat_approximation_nu", 7000. );
pars.add("transfer_neglect_delta_k_S_t0", 0.17 );
pars.add("transfer_neglect_delta_k_S_t1", 0.05 );
pars.add("transfer_neglect_delta_k_S_t2", 0.17 );
pars.add("transfer_neglect_delta_k_S_e", 0.13 );
pars.add("delta_l_max", 1000 );
std::unique_ptr<ClassEngine> CE = std::make_unique<ClassEngine>(pars, false);
CE->getTk(ztarget_, tab_lnk_, tab_dc_, tab_db_, d_ncdm, tab_dtot_,
@ -74,7 +101,9 @@ public:
Omega_b_ = pcf_->GetValue<double>("cosmology","Omega_b");
N_ur_ = pcf_->GetValueSafe<double>("cosmology","N_ur", 3.046);
ztarget_ = pcf_->GetValueSafe<double>("cosmology","ztarget",0.0);
atarget_ = 1.0/(1.0+ztarget_);
zstart_ = pcf_->GetValue<double>("setup","zstart");
astart_ = 1.0/(1.0+zstart_);
double lbox = pcf_->GetValue<double>("setup","BoxLength");
int nres = pcf_->GetValue<double>("setup","GridRes");
kmax_ = 2.0*M_PI/lbox * nres/2 * sqrt(3) * 2.0; // 120% of spatial diagonal
@ -102,6 +131,33 @@ public:
gsl_spline_init(gsl_sp_tc_, &tab_lnk_[0], &tab_tc_[0], tab_lnk_.size());
gsl_spline_init(gsl_sp_tb_, &tab_lnk_[0], &tab_tb_[0], tab_lnk_.size());
//--------------------------------------------------------------------------
// single fluid growing/decaying mode decomposition
//--------------------------------------------------------------------------
gsl_ia_Cplus_ = gsl_interp_accel_alloc();
gsl_ia_Cminus_ = gsl_interp_accel_alloc();
gsl_sp_Cplus_ = gsl_spline_alloc(gsl_interp_cspline, tab_lnk_.size());
gsl_sp_Cminus_ = gsl_spline_alloc(gsl_interp_cspline, tab_lnk_.size());
tab_Cplus_.assign(tab_lnk_.size(),0);
tab_Cminus_.assign(tab_lnk_.size(),0);
std::ofstream ofs("grow_decay.txt");
for( size_t i=0; i<tab_lnk_.size(); ++i ){
tab_Cplus_[i] = (3.0/5.0 * tab_dtot_[i]/atarget_ - 2.0/5.0*tab_ttot_[i]/atarget_);
tab_Cminus_[i] = (2.0/5.0 * std::pow(atarget_, 1.5) * ( tab_dtot_[i] + tab_ttot_[i] ));
ofs << std::exp(tab_lnk_[i]) << " " << tab_Cplus_[i] << " " << tab_Cminus_[i] << " " << tab_dtot_[i] << " " << tab_ttot_[i] << std::endl;
}
gsl_spline_init(gsl_sp_Cplus_, &tab_lnk_[0], &tab_Cplus_[0], tab_lnk_.size());
gsl_spline_init(gsl_sp_Cminus_, &tab_lnk_[0], &tab_Cminus_[0], tab_lnk_.size());
//--------------------------------------------------------------------------
kmin_ = std::exp(tab_lnk_[0]);
tf_distinct_ = true;