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

Added possibility to do partial normalization

Browse Source
This commit is contained in:
Guilhem Lavaux 2015-02-24 14:48:08 +01:00
parent 535f25d710
commit b7392204ad
3 changed files with 10 additions and 6 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

6
python/_cosmo_power.pyx
View File

@ -44,7 +44,7 @@ cdef extern from "cosmopower.hpp" namespace "CosmoTool":
void setFunction(CosmoFunction) void setFunction(CosmoFunction)
void updateCosmology() void updateCosmology()
void updatePhysicalCosmology() void updatePhysicalCosmology()
void normalize() void normalize(double)
void setNormalization(double) void setNormalization(double)
double power(double) double power(double)
@ -75,7 +75,7 @@ cdef class CosmologyPower:
self.power.updateCosmology() self.power.updateCosmology()
def normalize(self,s8): def normalize(self,s8,k_max=-1):
"""normalize(self, sigma8) """normalize(self, sigma8)
Compute the normalization of the power spectrum using sigma8. Compute the normalization of the power spectrum using sigma8.
@ -84,7 +84,7 @@ cdef class CosmologyPower:
sigma8 (float): standard deviation of density field smoothed at 8 Mpc/h sigma8 (float): standard deviation of density field smoothed at 8 Mpc/h
""" """
self.power.SIGMA8 = s8 self.power.SIGMA8 = s8
self.power.normalize() self.power.normalize(k_max)
def setFunction(self,funcname): def setFunction(self,funcname):

8
src/cosmopower.cpp
View File

@ -262,12 +262,16 @@ double CosmoPower::integrandNormalize(double x)
return power(k)*k*k*f*f/(x*x); return power(k)*k*k*f*f/(x*x);
} }
void CosmoPower::normalize() void CosmoPower::normalize(double k_max)
{ {
double normVal = 0; double normVal = 0;
double abserr; double abserr;
gsl_integration_workspace *w = gsl_integration_workspace_alloc(NUM_ITERATION); gsl_integration_workspace *w = gsl_integration_workspace_alloc(NUM_ITERATION);
gsl_function f; gsl_function f;
double x_min = 0;
if (k_max > 0)
x_min = 1/(1+k_max);
f.function = gslPowSpecNorm; f.function = gslPowSpecNorm;
f.params = this; f.params = this;
@ -282,7 +286,7 @@ void CosmoPower::normalize()
} }
// gsl_integration_qagiu(&f, 0, 0, TOLERANCE, NUM_ITERATION, w, &normVal, &abserr); // gsl_integration_qagiu(&f, 0, 0, TOLERANCE, NUM_ITERATION, w, &normVal, &abserr);
gsl_integration_qag(&f, 0, 1, 0, TOLERANCE, NUM_ITERATION, GSL_INTEG_GAUSS61, w, &normVal, &abserr); gsl_integration_qag(&f, x_min, 1, 0, TOLERANCE, NUM_ITERATION, GSL_INTEG_GAUSS61, w, &normVal, &abserr);
gsl_integration_workspace_free(w); gsl_integration_workspace_free(w);
normVal /= (2*M_PI*M_PI); normVal /= (2*M_PI*M_PI);

2
src/cosmopower.hpp
View File

@ -92,7 +92,7 @@ namespace CosmoTool {
void updateCosmology(); void updateCosmology();
void updatePhysicalCosmology(); void updatePhysicalCosmology();
void normalize(); void normalize(double k_max = -1);
void setNormalization(double A_K); void setNormalization(double A_K);
void updateHuWigglesConsts(); void updateHuWigglesConsts();