CosmoBolognaLib
Free Software C++/Python libraries for cosmological calculations
RSD.cpp
Go to the documentation of this file.
1 /********************************************************************
2  * Copyright (C) 2010 by Federico Marulli *
3  * federico.marulli3@unibo.it *
4  * *
5  * This program is free software; you can redistribute it and/or *
6  * modify it under the terms of the GNU General Public License as *
7  * published by the Free Software Foundation; either version 2 of *
8  * the License, or (at your option) any later version. *
9  * *
10  * This program is distributed in the hope that it will be useful, *
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
13  * GNU General Public License for more details. *
14  * *
15  * You should have received a copy of the GNU General Public *
16  * License along with this program; if not, write to the Free *
17  * Software Foundation, Inc., *
18  * 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. *
19  ********************************************************************/
20 
36 #include "Cosmology.h"
37 
38 using namespace std;
39 
40 using namespace cbl;
41 
42 
43 // =====================================================================================
44 
45 
46 double cbl::cosmology::Cosmology::fsigma8 (const double redshift, const string method_Pk, const bool store_output, const string output_root, const bool NL, const double k_min, const double k_max, const double prec, const string file_par) const
47 {
48  return linear_growth_rate(redshift)*sigma8_Pk(method_Pk, redshift, store_output, output_root, NL, k_min, k_max, prec, file_par);
49 }
50 
51 
52 // =====================================================================================
53 
54 
55 double cbl::cosmology::Cosmology::beta (const double redshift, const double bias) const
56 {
57  return linear_growth_rate(redshift)/bias;
58 }
59 
60 
61 // =====================================================================================
62 
63 
64 double cbl::cosmology::Cosmology::error_beta (const double redshift, const double bias, const double err_bias) const
65 {
66  return linear_growth_rate(redshift)/pow(bias, 2)*err_bias;
67 }
68 
69 
70 // =====================================================================================
71 
72 
73 double cbl::cosmology::Cosmology::beta (const double Mass_min, const double Mass_max, const double redshift, const std::string model_bias, const std::string model_MF, const std::string method_SS, const bool store_output, const std::string output_root, const double Delta, const double kk, const std::string interpType, const int norm, const double k_min, const double k_max, const double prec, const std::string input_file, const bool is_parameter_file)
74 {
75  return linear_growth_rate(redshift)/bias_eff(Mass_min, Mass_max, redshift, model_bias, model_MF, method_SS, store_output, output_root, Delta, kk, interpType, norm, k_min, k_max, prec, input_file, is_parameter_file);
76 }
77 
78 
79 // =====================================================================================
80 
81 
82 double cbl::cosmology::Cosmology::error_beta (const double Mass_min, const double Mass_max, const double redshift, const std::string model_bias, const std::string model_MF, const std::string method_SS, const double err_bias, const bool store_output, const std::string output_root, const double Delta, const double kk, const std::string interpType, const int norm, const double k_min, const double k_max, const double prec, const std::string input_file, const bool is_parameter_file)
83 {
84  return linear_growth_rate(redshift)/pow(bias_eff(Mass_min, Mass_max, redshift, model_bias, model_MF, method_SS, store_output, output_root, Delta, kk, interpType, norm, k_min, k_max, prec, input_file, is_parameter_file),2)*err_bias;
85 }
86 
87 
88 // =====================================================================================
89 
90 
91 double cbl::cosmology::Cosmology::beta (const std::vector<double> MM, const std::vector<double> MF, const double redshift, const std::string model_bias, const std::string method_SS, const bool store_output, const std::string output_root, const double Delta, const double kk, const std::string interpType, const int norm, const double k_min, const double k_max, const double prec, const std::string input_file, const bool is_parameter_file)
92 {
93  return linear_growth_rate(redshift)/bias_eff(MM, MF, redshift, model_bias, method_SS, store_output, output_root, Delta, kk, interpType, norm, k_min, k_max, prec, input_file, is_parameter_file);
94 }
95 
96 
97 // =====================================================================================
98 
99 
100 double cbl::cosmology::Cosmology::error_beta (const std::vector<double> MM, const std::vector<double> MF, const double redshift, const std::string model_bias, const std::string method_SS, const double err_bias, const bool store_output, const std::string output_root, const double Delta, const double kk, const std::string interpType, const int norm, const double k_min, const double k_max, const double prec, const std::string input_file, const bool is_parameter_file)
101 {
102  return linear_growth_rate(redshift)/pow(bias_eff(MM, MF, redshift, model_bias, method_SS, store_output, output_root, Delta, kk, interpType, norm, k_min, k_max, prec, input_file, is_parameter_file),2)*err_bias;
103 }
104 
105 
106 // =====================================================================================
107 
108 
109 double cbl::cosmology::Cosmology::error_beta_measured (const double Volume, const double density, const double Mass_min, const double Mass_max, const double redshift, const std::string model_bias, const std::string model_MF, const std::string method_SS, const bool store_output, const std::string output_root, const double Delta, const double kk, const std::string interpType, const int norm, const double k_min, const double k_max, const double prec, const std::string input_file, const bool is_parameter_file)
110 { // from Eq. 20 of Bianchi et al. 2012
111 
112  double bias = bias_eff(Mass_min, Mass_max, redshift, model_bias, model_MF, method_SS, store_output, output_root, Delta, kk, interpType, norm, k_min, k_max, prec, input_file, is_parameter_file);
113 
114  return relative_error_beta (bias, Volume, density);
115 }
116 
117 
118 // =====================================================================================
119 
120 
121 double cbl::cosmology::Cosmology::quadrupole (const double Mass_min, const double Mass_max, const double redshift, const std::string model_bias, const std::string model_MF, const std::string method_SS, const bool store_output, const std::string output_root, const double Delta, const double kk, const std::string interpType, const int norm, const double k_min, const double k_max, const double prec, const std::string input_file, const bool is_parameter_file)
122 {
123  double Beta = beta(Mass_min, Mass_max, redshift, model_bias, model_MF, method_SS, store_output, output_root, Delta, kk, interpType, norm, k_min, k_max, prec, input_file, is_parameter_file);
124  return (4./3.*Beta+4./7.*Beta*Beta)/(1.+2./3.*Beta+1./5.*Beta*Beta);
125 }
126 
127 
128 // =====================================================================================
129 
130 
131 double cbl::cosmology::Cosmology::quadrupole (const std::vector<double> MM, const std::vector<double> MF, const double redshift, const std::string model_bias, const std::string method_SS, const bool store_output, const std::string output_root, const double Delta, const double kk, const std::string interpType, const int norm, const double k_min, const double k_max, const double prec, const std::string input_file, const bool is_parameter_file)
132 {
133  double Beta = beta(MM, MF, redshift, model_bias, method_SS, store_output, output_root, Delta, kk, interpType, norm, k_min, k_max, prec, input_file, is_parameter_file);
134  return (4./3.*Beta+4./7.*Beta*Beta)/(1.+2./3.*Beta+1./5.*Beta*Beta);
135 }
136 
137 
138 // =====================================================================================
139 
140 
141 double cbl::cosmology::Cosmology::Pk_DeltaDelta_fitting_function (const double kk, const std::string method_Pk, const double redshift, const std::string author, const bool store_output, const std::string output_root, const int norm, double k_min, double k_max, const double prec, const std::string file_par, const bool unit1)
142 {
143  double Pkdd = 0;
144  if(author == "Pezzotta" || author == "Bel")
145  Pkdd = Pk_matter({kk}, method_Pk, false, redshift, store_output, output_root, norm, k_min, k_max, prec, file_par, unit1)[0];
146  else WarningMsgCBL("the current implementation is not correct with author = " + author, "Pk_DeltaDelta_fitting_function", "RSD.cpp");
147  return Pkdd;
148 }
149 
150 
151 // =====================================================================================
152 
153 
154 double cbl::cosmology::Cosmology::Pk_DeltaTheta_fitting_function (const double kk, const std::string method_Pk, const double redshift, const std::string author, const bool store_output, const std::string output_root, const bool NL, const int norm, double k_min, double k_max, const double prec, const std::string file_par, const bool unit1)
155 {
156  double sigma8_z = sigma8_Pk(method_Pk, redshift, store_output, output_root, NL, k_min, k_max, prec, file_par);
157  double kd = 1./(-0.017 + 1.496*pow(sigma8_z, 2.));
158  double Pkdt = 0;
159  if (author == "Pezzotta"){
160  Pkdt = sqrt(Pk_matter({kk}, method_Pk, true, redshift, store_output, output_root, norm, k_min, k_max, prec, file_par, unit1)[0]*Pk_matter({kk}, method_Pk, false, redshift, store_output, output_root, norm, k_min, k_max, prec, file_par, unit1)[0])*exp(-kk/kd);
161  }
162  else if (author == "Bel"){
163  double b = 0.091 + 0.702*sigma8_z*sigma8_z;
164  Pkdt = sqrt(Pk_matter({kk}, method_Pk, true, redshift, store_output, output_root, norm, k_min, k_max, prec, file_par, unit1)[0]*Pk_matter({kk}, method_Pk, false, redshift, store_output, output_root, norm, k_min, k_max, prec, file_par, unit1)[0])*exp(-kk/kd-b*pow(kk,6.0));
165  }
166  else WarningMsgCBL("the current implementation is not correct with author = " + author, "Pk_DeltaTheta_fitting_function", "RSD.cpp");
167  return Pkdt;
168 }
169 
170 
171 // =====================================================================================
172 
173 
174 double cbl::cosmology::Cosmology::Pk_ThetaTheta_fitting_function (const double kk, const std::string method_Pk, const double redshift, const std::string author, const bool store_output, const std::string output_root, const bool NL, const int norm, double k_min, double k_max, const double prec, const std::string file_par, const bool unit1)
175 {
176  double sigma8_z = sigma8_Pk(method_Pk, redshift, store_output, output_root, NL, k_min, k_max, prec, file_par);
177  double Pktt = 0;
178  if (author == "Pezzotta"){
179  double kt = 1./(-0.048 + 1.917*sigma8_z*sigma8_z);
180  Pktt = Pk_matter({kk}, method_Pk, false, redshift, store_output, output_root, norm, k_min, k_max, prec, file_par, unit1)[0]*exp(-kk/kt);
181  }
182 
183  else if (author == "Bel"){
184  double a1 = -0.817 + 3.198*sigma8_z;
185  double a2 = 0.877 - 4.191*sigma8_z;
186  double a3 = -1.199 + 4.629*sigma8_z;
187  Pktt = Pk_matter({kk}, method_Pk, false, redshift, store_output, output_root, norm, k_min, k_max, prec, file_par, unit1)[0]*exp(-kk*(a1 + a2*kk + a3*kk*kk));
188  }
189 
190  else WarningMsgCBL("the current implementation is not correct with author = " + author, "Pk_ThetaTheta_fitting_function", "RSD.cpp");
191 
192  return Pktt;
193 }
194 
195 
196 // =====================================================================================
197 
198 
199 double cbl::cosmology::Cosmology::sigma_v (const double redshift, const std::string method_Pk, const bool store_output, const std::string output_root, const int norm, const double k_min, const double k_max, const int bin_k, const double prec, const std::string file_par, const bool unit1)
200 {
201  const vector<double> kk = logarithmic_bin_vector(bin_k, k_min, k_max);
202  const vector<double> Pk = Pk_matter(kk, method_Pk, false, redshift, store_output, output_root, norm, k_min, k_max, prec, file_par, unit1);
203 
204  auto interp_Pk = glob::FuncGrid(kk, Pk, "Spline");
205 
206  auto integrand = [&] (const double log_q)
207  {
208  const double qq = exp(log_q);
209  return qq*interp_Pk(qq);
210  };
211 
212  return sqrt(1./(6.*pow(par::pi, 2))*wrapper::gsl::GSL_integrate_cquad(integrand, log(k_min), log(k_max), 1.e-5, 0., 10000));
213 }
The class Cosmology.
double sigma_v(const double redshift=0., const std::string method_Pk="CAMB", const bool store_output=true, const std::string output_root="test", const int norm=-1, const double k_min=0.001, const double k_max=100., const int bin_k=512, const double prec=1.e-2, const std::string file_par=par::defaultString, const bool unit1=false)
the linear-order one-dimensional pairwise velocity dispersion,
Definition: RSD.cpp:199
double fsigma8(const double redshift, const std::string method_Pk, const bool store_output=true, const std::string output_root="test", const bool NL=false, const double k_min=0.001, const double k_max=100., const double prec=1.e-2, const std::string file_par=par::defaultString) const
f*σ8: the linear growth rate times the dark matter rms mass fluctuation within 8 Mpc/h
Definition: RSD.cpp:46
double quadrupole(const double Mass_min, const double Mass_max, const double redshift, const std::string model_bias, const std::string model_MF, const std::string method_SS, const bool store_output=true, const std::string output_root="test", const double Delta=200., const double kk=-1., const std::string interpType="Linear", const int norm=-1, const double k_min=0.001, const double k_max=100., const double prec=1.e-2, const std::string input_file=par::defaultString, const bool is_parameter_file=true)
the normalised quadrupole Q
Definition: RSD.cpp:121
double error_beta(const double redshift, const double bias, const double err_bias) const
the error on the specific growth rate β
Definition: RSD.cpp:64
double beta(const double redshift, const double bias) const
the specific growth rate β
Definition: RSD.cpp:55
double error_beta_measured(const double Volume, const double density, const double Mass_min, const double Mass_max, const double redshift, const std::string model_bias, const std::string model_MF, const std::string method_SS, const bool store_output=true, const std::string output_root="test", const double Delta=200., const double kk=-1., const std::string interpType="Linear", const int norm=-1, const double k_min=0.001, const double k_max=100., const double prec=1.e-2, const std::string input_file=par::defaultString, const bool is_parameter_file=true)
the error on the specific growth rate β from Bianchi et al. 2012
Definition: RSD.cpp:109
double Pk_ThetaTheta_fitting_function(const double kk, const std::string method_Pk, const double redshift, const std::string author, const bool store_output, const std::string output_root, const bool NL, const int norm, double k_min, double k_max, const double prec, const std::string file_par, const bool unit1)
the dark matter velocity divergence power spectrum
Definition: RSD.cpp:174
double Pk_DeltaTheta_fitting_function(const double kk, const std::string method_Pk, const double redshift, const std::string author, const bool store_output, const std::string output_root, const bool NL, const int norm, double k_min, double k_max, const double prec, const std::string file_par, const bool unit1)
the dark matter cross power spectrum
Definition: RSD.cpp:154
double Pk_DeltaDelta_fitting_function(const double kk, const std::string method_Pk, const double redshift, const std::string author, const bool store_output, const std::string output_root, const int norm, double k_min, double k_max, const double prec, const std::string file_par, const bool unit1)
the non-linear dark matter power spectrum using fitting functions given by Bel et....
Definition: RSD.cpp:141
The class FuncGrid.
Definition: FuncGrid.h:55
static const double pi
: the ratio of a circle's circumference to its diameter
Definition: Constants.h:199
double bias(const double Mmin, const double sigmalgM, const double M0, const double M1, const double alpha, const std::shared_ptr< void > inputs)
the mean galaxy bias
double GSL_integrate_cquad(gsl_function func, const double a, const double b, const double rel_err=1.e-3, const double abs_err=0, const int nevals=100)
integral, using the gsl cquad method
Definition: GSLwrapper.cpp:159
The global namespace of the CosmoBolognaLib
Definition: CAMB.h:38
std::vector< T > logarithmic_bin_vector(const size_t nn, const T min, const T max)
fill a std::vector with logarithmically spaced values
Definition: Kernel.h:1621
double relative_error_beta(const double bias, const double Volume, const double density)
estimated relative error on
Definition: Func.cpp:1072
void WarningMsgCBL(const std::string msg, const std::string functionCBL, const std::string fileCBL)
internal CBL warning message
Definition: Kernel.h:747