Package pycosmicstar :: Module structures
[hide private]
[frames] | no frames]

Source Code for Module pycosmicstar.structures

  1  #!/usr/bin/env python3 
  2  # *-* Coding: UTF-8 *-* 
  3  from __future__ import division, absolute_import 
  4   
  5  __author__ = "Eduardo dos Santos Pereira" 
  6  __email__ = "pereira.somoza@gmail.com" 
  7  __credits__ = ["Eduardo dos Santos Pereira"] 
  8  __license__ = "GPLV3" 
  9  __version__ = "1.0.1" 
 10  __maintainer__ = "Eduardo dos Santos Pereira" 
 11  __status__ = "Stable" 
 12   
 13   
 14  """Cosmological Dark Halos History 
 15  From the formalism of Reed et al (MNRAS, 346, 565-572, 2003) 
 16  it is calculated the mass fraction of dark matter halos. 
 17  The code obtain the mass density of dark halos and the fraction 
 18  of brions into structures as a function of the time. 
 19  Here is used the Transfer function from Efstathiou, Bond & White 
 20  -- (MNRAS, 258, 1P, 1992). 
 21  The current version it is assumed the normalization of WMAP (withou 
 22  gravitational waves) adapted from Eisenstein e Hu (ApJ 511, 5, 1999) that 
 23  in the way that return  sigma_8 = 0,84. 
 24  The fraction of mass of dark halos is obtained by the work of 
 25  Sheth e Tormen (MNRAS 308, 119, 1999). 
 26  All models consider Omega_Total = Omega_M + Omega_L = 1,0 
 27   
 28  This file is part of pystar. 
 29  copyright : Eduardo dos Santos Pereira 
 30   
 31  pystar is free software: you can redistribute it and/or modify 
 32  it under the terms of the GNU General Public License as published by 
 33  the Free Software Foundation, either version 3 of the License. 
 34  pystar is distributed in the hope that it will be useful, 
 35  but WITHOUT ANY WARRANTY; without even the implied warranty of 
 36  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the 
 37  GNU General Public License for more details. 
 38   
 39  You should have received a copy of the GNU General Public License 
 40  along with Foobar.  If not, see <http://www.gnu.org/licenses/>. 
 41   
 42  """ 
 43   
 44  from numpy import sqrt, pi, log, log10, exp, array, abs 
 45  import scipy.interpolate as spint 
 46  from scipy.interpolate import InterpolatedUnivariateSpline as spline 
 47  from scipy.special import gamma 
 48   
 49  from .structuresabstract import structuresabstract 
 50   
 51  import sys 
 52  pyversion = sys.version_info 
 53  if pyversion[0] >= 3: 
 54      from . import filedict 
 55  else: 
 56      print("Importing filedict for python2.7") 
 57      from . import filedict_old as filedict 
 58   
 59  import os 
 60  from .diferencial import dfridr, locate 
 61   
 62  from .paralleloverlist import parallel_list 
 63   
 64   
65 -class structures(structuresabstract):
66 """This class was contructed based in the like Press-Schechter formalism 67 that provides characteristis like numerical density of dark matter halos 68 into the range m_h, m_h + dm_h, the fraction of barionic matter, 69 and, the accretion rate of barions into structures and the total number 70 of dark halos. 71 72 The models used to develop this class was presented for the first time 73 in the article of Pereira and Miranda (2010) - (MNRAS, 401, 1924, 2010). 74 75 The cosmologic background model is passed as a instance parameter: 76 cosmology 77 78 Keyword arguments: 79 lmin -- (default 6.0) log10 of the minal mass of the dark halo 80 where it is possible to have star formation. 81 82 zmax -- (defaul 20.0) - the maximum redshift to be considered 83 84 omegam -- (default 0.24) - The dark matter parameter 85 86 omegab -- (default 0.04) - The barionic parameter 87 88 omegal -- (default 0.73) - The dark energy parameter 89 90 h -- (default 0.7) - The h of the Hubble constant (H = h * 100) 91 92 massFunctionType: 93 (Dark Haloes Mass Function) 94 default 'ST' - Sheth et al. (2001) - z=[0,2] 95 'TK' - Tinker et al. (2008) - z=[0,2.5] 96 'PS' - Press and Schechter (1974) - z=- 97 'JK' - Jenkins et al. (2001) z=[0,5] 98 'W' - Warren et al. (2006) z=0 99 'WT1' - Watson et al. (2013) - Tinker Modified - z=[0,30] 100 'WT2' - Watson et al. (2013) - Gamma times times Tinker Modified 101 z=[0,30] 102 'B' - Burr Distribuction. Marassi and Lima (2006) - Press Schechter 103 modified. 104 qBurr: 105 (default 1) - The q value of Burr Distribuction. 106 107 """ 108
109 - def __init__(self, cosmology, **kwargs):
110 111 listParameters = ["lmin", "zmax", 112 "omegam", "omegab", "omegal", "h", 113 "cacheDir", "cacheFile", "massFunctionType", 114 "delta_halo", "qBurr"] 115 116 testeKeysArgs = [Ki for Ki in list(kwargs.keys()) 117 if Ki not in listParameters] 118 119 if(len(testeKeysArgs) >= 1): 120 nameError = "The key args are not defined:" 121 for i in range(len(testeKeysArgs)): 122 nameError += testeKeysArgs[i] 123 raise NameError(nameError) 124 125 #lmin=6.0, zmax=20.0, 126 127 if 'lmin' in list(kwargs.keys()): 128 lmin = kwargs['lmin'] 129 else: 130 lmin = 6.0 131 132 if 'zmax' in list(kwargs.keys()): 133 zmax = kwargs['zmax'] 134 else: 135 zmax = 20.0 136 137 #omegam=0.24, omegab=0.04, omegal=0.73, h=0.7, 138 139 if 'omegam' in list(kwargs.keys()): 140 omegam = kwargs['omegam'] 141 else: 142 omegam = 0.24 143 144 if 'omegab' in list(kwargs.keys()): 145 omegab = kwargs['omegab'] 146 else: 147 omegab = 0.04 148 149 if 'omegal' in list(kwargs.keys()): 150 omegal = kwargs['omegal'] 151 else: 152 omegal = 0.73 153 154 if 'h' in list(kwargs.keys()): 155 h = kwargs['h'] 156 else: 157 h = 0.7 158 159 #cacheDir=None, cacheFile=None, 160 161 if 'cacheDir' in list(kwargs.keys()): 162 cacheDir = kwargs['cacheDir'] 163 else: 164 cacheDir = None 165 166 if 'cacheFile' in list(kwargs.keys()): 167 cacheFile = kwargs['cacheFile'] 168 else: 169 cacheFile = None 170 171 #massFunctionType="ST", delta_halo=200, qBurr=1 172 173 if 'massFunctionType' in list(kwargs.keys()): 174 massFunctionType = kwargs['massFunctionType'] 175 else: 176 massFunctionType = "ST" 177 178 if 'delta_halo' in list(kwargs.keys()): 179 delta_halo = kwargs['delta_halo'] 180 else: 181 delta_halo = 200 182 183 if 'qBurr' in list(kwargs.keys()): 184 qBurr = kwargs['qBurr'] 185 else: 186 qBurr = 1 187 188 self._cosmology = cosmology(omegam, omegab, omegal, h) 189 190 if(cacheDir is None): 191 self._cacheDir = self. __creatCachDiretory()[0] 192 else: 193 self._cacheDir = cacheDir 194 195 if(cacheFile is None): 196 if(massFunctionType == "TK"): 197 198 cacheFile = str(self._cacheDir) + "/structures_cache_"\ 199 + massFunctionType + str(delta_halo) + "_" + "_" +\ 200 str(omegab) + "_" \ 201 + str(omegam) + "_" +\ 202 str(omegal) + "_ " \ 203 + str(h) + "_" + str(lmin) + "_" + str(zmax) 204 205 else: 206 cacheFile = str(self._cacheDir) + "/structures_cache_"\ 207 + massFunctionType + "_" + str(omegab) + "_" \ 208 + str(omegam) + "_" + str(omegal) + "_ " \ 209 + str(h) + "_" + str(lmin) + "_" + str(zmax) 210 else: 211 cacheFile = str(self._cacheDir) + cacheFile 212 213 self._cacheFIle = cacheFile 214 215 self._cache_dict = filedict.FileDict(filename=cacheFile + ".cache") 216 217 self.__mmin = 1.0e+4 218 self.__mmax = 1.0e+18 219 self.__lmax = log10(self.__mmax / 10.0) 220 self._zmax = zmax 221 self.__lmin = lmin 222 self.__deltac = self._cosmology.getDeltaC() 223 self.__pst = 0.3 224 225 h2 = h * h 226 h2om = h2 * omegam 227 #h2br = h2 * omegab 228 self.__ut = 1.0 / 3.0 229 self.__nr = 14000 230 self.__ct0 = 4.0 * pi 231 self.__ct1 = self.__ct0 * 2.76e+11 / 3.0 232 self.__ct2 = self.__ct1 * h2om 233 self.__ast1 = 0.322 234 self.__ast2 = 0.707 235 self.__pst = 0.3 236 self.__tilt2 = self._cosmology.getTilt() / log(10.0) 237 self.__ctst = self.__ast1 * sqrt(2.0 * self.__ast2 / pi) 238 239 #Burr q coeficiente 240 self.__qBurr = qBurr 241 242 self.__massFunctionType = massFunctionType 243 self.__delta_halo = delta_halo 244 245 self.__massFunctionDict = {"ST": self.__massFunctionST, 246 "TK": self.__massFunctionTinker, 247 "PS": self.__massFunctionPressSchechter, 248 "JK": self.__massFunctionJenkins, 249 "W": self.__massFunctionW, 250 "WT1": self.__massFunctionWT1, 251 "WT2": self.__massFunctionWT2, 252 "B": self.__massFunctionBurr 253 } 254 255 self.__startingSigmaAccretion()
256
257 - def __creatCachDiretory(self):
258 HOME = os.path.expanduser('~') 259 if not os.path.exists(HOME + '/.cosmicstarformation'): 260 print(('Creating .cosmicstarformation cache diretory in %s' 261 % HOME)) 262 os.makedirs(HOME + '/.cosmicstarformation') 263 return os.path.expanduser('~') + '/.cosmicstarformation', True
264
265 - def __ifSigmaNotInCache(self):
266 """Calculate the values necessaries to initialize the 267 numerical function of sigma 268 """ 269 numk = 1000.0 270 kscale = self.__mmax / self.__mmin 271 kls = log10(kscale) 272 numk = numk * kls 273 kls1 = kls / numk 274 deltaz = self._zmax / (numk) 275 276 def CalculaKm(i): 277 kmass = (10.0 ** ((i + 1) * kls1)) * self.__mmin 278 return kmass
279 280 def CalculaScale(i): 281 scale = (CalculaKm(i) / self.__ct2) ** self.__ut 282 return scale
283 284 self.__kmass = array([CalculaKm(i) for i in range(int(numk))]) 285 self.__scale = array([CalculaScale(i) for i in range(int(numk))]) 286 self.__zred = array([self._zmax - i * deltaz 287 for i in range(int(numk))]) 288 289 e, f = self._cosmology.sigma(self.__kmass) 290 291 self.__km = array([ei for ei in e]) 292 self.__sg = array([FI for FI in f]) 293 294 self.__t_z = parallel_list(self._cosmology.age, self.__zred) 295 296 self.__d_c2 = parallel_list(self.__deltaCz, self.__zred) 297 298 self.__rdm2 = parallel_list(self.__rodmz, self.__zred) 299 300 self.__rbr2 = parallel_list(self._cosmology.robr, self.__zred) 301 302 #self.__t_z = array([self._cosmology.age(zi) for zi in self.__zred]) 303 304 #self.__d_c2 = array([ 305 #self.__deltac / self._cosmology.growthFunction(zi) 306 #for zi in self.__zred 307 #]) 308 #self.__rdm2 = array([self._cosmology.rodm(zi)[0] 309 #for zi in self.__zred 310 #]) 311 #self.__rbr2 = array([self._cosmology.robr(zi) 312 #for zi in self.__zred]) 313
314 - def __rodmz(self, z):
315 return self._cosmology.rodm(z)[0]
316
317 - def __deltaCz(self, z):
318 return self.__deltac / self._cosmology.growthFunction(z)
319
320 - def __startingSigmaAccretion(self):
321 """ 322 Verify if the values are in cache 323 """ 324 325 try: 326 327 self.__km = self._cache_dict['km'] 328 self.__scale = self._cache_dict['scale'] 329 self.__zred = self._cache_dict['zred'] 330 self.__sg = self._cache_dict['sg'] 331 self.__t_z = self._cache_dict['t_z'] 332 self.__d_c2 = self._cache_dict['d_c2'] 333 self.__rdm2 = self._cache_dict['rdm2'] 334 self.__rbr2 = self._cache_dict['rbr2'] 335 self._abt2 = self._cache_dict['abt2'] 336 self._ascale = self._cache_dict['ascale'] 337 self._tck_ab = self._cache_dict['tck_ab'] 338 339 print("\nStructures Data in Cache\n") 340 341 except: 342 self.__ifSigmaNotInCache() 343 self.__startBarionicAccretionRate() 344 self.__cachingAtribut()
345
346 - def __cachingAtribut(self):
347 """Caching the values 348 """ 349 self._cache_dict['km'] = self.__km 350 self._cache_dict['scale'] = self.__scale 351 self._cache_dict['zred'] = self.__zred 352 self._cache_dict['sg'] = self.__sg 353 self._cache_dict['t_z'] = self.__t_z 354 self._cache_dict['d_c2'] = self.__d_c2 355 self._cache_dict['rdm2'] = self.__rdm2 356 self._cache_dict['rbr2'] = self.__rbr2 357 self._cache_dict['abt2'] = self._abt2 358 self._cache_dict['ascale'] = self._ascale 359 self._cache_dict['tck_ab'] = self._tck_ab
360
361 - def massFunction(self, lm, z):
362 """Return the mass function of dark halos. 363 364 Keyword arguments: 365 lm -- log10 of the mass of the dark halo 366 z -- redshift 367 """ 368 try: 369 return self.__massFunctionDict[self.__massFunctionType](lm, z) 370 except: 371 raise NameError("No Defined Mass Function")
372
373 - def __massFunctionJenkins(self, lm, z):
374 """Return the mass function of Jenkins et al. (2003). 375 Keyword arguments: 376 lm -- log10 of the mass of the dark halo 377 z -- redshift 378 """ 379 380 gte = self._cosmology.growthFunction(z) 381 rdmt, drdmt = self._cosmology.rodm(z) 382 step = lm / 2.0e+1 383 kmsgm = lm 384 kmass = 10.0 ** (kmsgm) 385 sgm = self.fstm(lm) 386 dsgm_dlgm = dfridr(self.fstm, lm, step, err=0.0) 387 sigma1 = self.__deltac / (sgm * gte) 388 #sigma2 = sigma1 ** 2.0 389 390 dsgm_dlgm = dfridr(self.fstm, lm, step, err=0.0) 391 392 fst = 0.315 * exp(- abs(log(sigma1) + 0.61) ** 3.8) 393 394 frst = (rdmt / kmass ** 2.0) * fst * abs(dsgm_dlgm) / sgm 395 dn_dm = frst 396 return dn_dm
397
398 - def __massFunctionPressSchechter(self, lm, z):
399 """Return the value of Press-Schechter (1974) mass function. 400 Keyword arguments: 401 lm -- log10 of the mass of the dark halo 402 z -- redshift 403 """ 404 405 gte = self._cosmology.growthFunction(z) 406 rdmt, drdmt = self._cosmology.rodm(z) 407 step = lm / 2.0e+1 408 kmsgm = lm 409 kmass = 10.0 ** (kmsgm) 410 sgm = self.fstm(lm) 411 dsgm_dlgm = dfridr(self.fstm, lm, step, err=0.0) 412 sigma1 = self.__deltac / (sgm * gte) 413 sigma2 = sigma1 ** 2.0 414 fst = sqrt(2.0 / pi) * (sigma1) * exp(-0.5 * sigma2) 415 frst = (rdmt / kmass ** 2.0) * fst * abs(dsgm_dlgm) / sgm 416 dn_dm = frst 417 return dn_dm
418
419 - def _masFunctionWT0(self, lm, z, A, a, b, c):
420 rdmt, drdmt = self._cosmology.rodm(z) 421 step = lm / 2.0e+1 422 kmsgm = lm 423 kmass = 10.0 ** (kmsgm) 424 sgm = self.fstm(lm) 425 dsgm_dlgm = dfridr(self.fstm, lm, step, err=0.0) 426 427 dsgm_dlgm = dfridr(self.fstm, lm, step, err=0.0) 428 429 fst = A * ((sgm / b) ** (-a) + 1) * exp(-c / sgm ** 2.0) 430 431 frst = (rdmt / kmass ** 2.0) * fst * abs(dsgm_dlgm) / sgm 432 dn_dm = frst 433 return dn_dm
434
435 - def __massFunctionWT1(self, lm, z):
436 """ 437 Return the value of Watson et al. (2013) (- Tinker Modified - z=[0,30]) 438 mass function 439 Keyword arguments: 440 lm -- log10 of the mass of the dark halo 441 z -- redshift 442 """ 443 A = 0.282 444 a = 2.163 445 b = 1.406 446 c = 1.21 447 return self._masFunctionWT0(lm, z, A, a, b, c)
448
449 - def __massFunctionWT2(self, lm, z):
450 """ 451 Return the value of Watson et al. (2013) (Gamma times 452 Tinker Modified z=[0,30]) mass function 453 Keyword arguments: 454 lm -- log10 of the mass of the dark halo 455 z -- redshift 456 """ 457 458 raise NameError("Not implemented yet") 459 460 if(z < 0): 461 raise NameError("z lower than zero.") 462 463 if(z == 0): 464 A = 0.194 465 a = 2.267 466 b = 1.805 467 c = 1.287 468 elif(z > 6): 469 A = 0.563 470 a = 0.874 471 b = 3.874 472 c = 1.453 473 474 else: 475 pass 476 477 return self._masFunctionWT0(lm, z, A, a, b, c)
478
479 - def _burrBq(self):
480 if(self.__qBurr > 0.0 and self.__qBurr < 1.0): 481 Bq = ((1.0 - self.__qBurr) ** 0.5) * ((3.0 - self.__qBurr) / 2.0)\ 482 * gamma(0.5 + 1.0 / (1.0 - self.__qBurr))\ 483 / gamma(1.0 / (1.0 - self.__qBurr)) 484 elif(self.__qBurr >= 1.0 and self.__qBurr < 2.0): 485 Bq = ((self.__qBurr - 1.0) ** 0.5)\ 486 * gamma(1.0 / (self.__qBurr - 1.0))\ 487 / gamma(1.0 / (self.__qBurr - 1.0) - 0.5) 488 else: 489 raise NameError('q of Burr function out of valide range') 490 491 return Bq
492
493 - def __massFunctionBurr(self, lm, z):
494 """ 495 Return the value of the Burr Distribution (Marassi and Lima (2006)) 496 - Press Schechter modified, mass function. 497 Keyword arguments: 498 lm -- log10 of the mass of the dark halo 499 z -- redshift 500 """ 501 if(self.__qBurr is None): 502 raise NameError('The Burr coeficient is None.') 503 504 gte = self._cosmology.growthFunction(z) 505 rdmt, drdmt = self._cosmology.rodm(z) 506 step = lm / 2.0e+1 507 kmsgm = lm 508 kmass = 10.0 ** (kmsgm) 509 sgm = self.fstm(lm) 510 dsgm_dlgm = dfridr(self.fstm, lm, step, err=0.0) 511 sigma1 = self.__deltac / (sgm * gte) 512 sigma2 = sigma1 ** 2.0 513 fst = self._burrBq() * sqrt(2.0 / pi) * (sigma1) * ( 514 1.0 - (1.0 - self.__qBurr) * 0.5 * sigma2 515 ) ** (1.0 / (1.0 - self.__qBurr)) 516 517 frst = (rdmt / kmass ** 2.0) * fst * abs(dsgm_dlgm) / sgm 518 dn_dm = frst 519 return dn_dm
520
521 - def __massFunctionTinker(self, lm, z):
522 """Return the mass function of dark halos of 523 Tinker mass function (Tinker et al. 2008) 524 525 This function was adapted from the work of: 526 S.G. Murray et al. 2013. Astronomy and Computing. 3-4. 23-34. 527 source of the original (https://github.com/steven-murray/hmf) 528 529 Keyword arguments: 530 lm -- log10 of the mass of the dark halo 531 z -- redshift 532 """ 533 534 #The Tinker function is a bit tricky - we use the code from 535 #http://cosmo.nyu.edu/~tinker/massfunction/MF_code.tar 536 #to aid us. 537 delta_virs = array([200, 300, 400, 600, 800, 1200, 1600, 2400, 3200]) 538 A_array = array([1.858659e-01, 539 1.995973e-01, 540 2.115659e-01, 541 2.184113e-01, 542 2.480968e-01, 543 2.546053e-01, 544 2.600000e-01, 545 2.600000e-01, 546 2.600000e-01]) 547 548 a_array = array([1.466904e+00, 549 1.521782e+00, 550 1.559186e+00, 551 1.614585e+00, 552 1.869936e+00, 553 2.128056e+00, 554 2.301275e+00, 555 2.529241e+00, 556 2.661983e+00]) 557 558 b_array = array([2.571104e+00, 559 2.254217e+00, 560 2.048674e+00, 561 1.869559e+00, 562 1.588649e+00, 563 1.507134e+00, 564 1.464374e+00, 565 1.436827e+00, 566 1.405210e+00]) 567 568 c_array = array([1.193958e+00, 569 1.270316e+00, 570 1.335191e+00, 571 1.446266e+00, 572 1.581345e+00, 573 1.795050e+00, 574 1.965613e+00, 575 2.237466e+00, 576 2.439729e+00]) 577 A_func = spline(delta_virs, A_array) 578 a_func = spline(delta_virs, a_array) 579 b_func = spline(delta_virs, b_array) 580 c_func = spline(delta_virs, c_array) 581 582 A_0 = A_func(self.__delta_halo) 583 a_0 = a_func(self.__delta_halo) 584 b_0 = b_func(self.__delta_halo) 585 c_0 = c_func(self.__delta_halo) 586 587 A = A_0 * (1 + z) ** (-0.14) 588 a = a_0 * (1 + z) ** (-0.06) 589 alpha = exp(-(0.75 / log(self.__delta_halo / 75)) ** 1.2) 590 b = b_0 * (1 + z) ** (-alpha) 591 c = c_0 592 593 #gte = self._cosmology.growthFunction(z) 594 rdmt, drdmt = self._cosmology.rodm(z) 595 step = lm / 2.0e+1 596 kmsgm = lm 597 kmass = 10.0 ** (kmsgm) 598 sgm = self.fstm(lm) 599 dsgm_dlgm = dfridr(self.fstm, lm, step, err=0.0) 600 #sigma1 = self.__deltac / (sgm * gte) 601 #sigma2 = sigma1 ** 2.0 602 603 dsgm_dlgm = dfridr(self.fstm, lm, step, err=0.0) 604 605 fst = A * ((sgm / b) ** (-a) + 1) * exp(-c / sgm ** 2.0) 606 607 frst = (rdmt / kmass ** 2.0) * fst * abs(dsgm_dlgm) / sgm 608 dn_dm = frst 609 return dn_dm
610
611 - def __massFunctionW(self, lm, z):
612 # LANL fitting function - Warren et al. 2005, astro-ph/0506395, eqtn. 5 613 A = 0.7234 614 a = 1.625 615 b = 0.2538 616 c = 1.1982 617 618 gte = self._cosmology.growthFunction(z) 619 rdmt, drdmt = self._cosmology.rodm(z) 620 step = lm / 2.0e+1 621 kmsgm = lm 622 kmass = 10.0 ** (kmsgm) 623 sgm = self.fstm(lm) 624 dsgm_dlgm = dfridr(self.fstm, lm, step, err=0.0) 625 sigma1 = self.__deltac / (sgm * gte) 626 sigma2 = sigma1 ** 2.0 627 628 dsgm_dlgm = dfridr(self.fstm, lm, step, err=0.0) 629 630 fst = A * ((sigma1 ** (-a)) + b) * exp(-c / sigma2) 631 632 frst = (rdmt / kmass ** 2.0) * fst * abs(dsgm_dlgm) / sgm 633 dn_dm = frst 634 return dn_dm
635
636 - def __massFunctionST(self, lm, z):
637 """Return the mass function of dark halos of 638 Sheth e Tormen (MNRAS 308, 119, 1999). 639 640 Keyword arguments: 641 lm -- log10 of the mass of the dark halo 642 z -- redshift 643 """ 644 gte = self._cosmology.growthFunction(z) 645 #gte2 = self._cosmology.dgrowth_dt(z) 646 rdmt, drdmt = self._cosmology.rodm(z) 647 step = lm / 2.0e+1 648 kmsgm = lm 649 kmass = 10.0 ** (kmsgm) 650 sgm = self.fstm(lm) 651 dsgm_dlgm = dfridr(self.fstm, lm, step, err=0.0) 652 sigma1 = self.__deltac / (sgm * gte) 653 sigma2 = sigma1 ** 2.0 654 expn = exp(-self.__ast2 * sigma2 / 2.0) 655 fst = self.__ctst * sigma1 * \ 656 (1.0 + (1.0 / (sigma2 * self.__ast2)) ** self.__pst) * expn 657 frst = (rdmt / kmass ** 2.0) * fst * abs(dsgm_dlgm) / sgm 658 dn_dm = frst 659 return dn_dm
660
661 - def fstm(self, lm):
662 '''Numerical function that return the value of sigm that 663 will be used by dfridr to calculate d_sigma_dlog10(m). 664 665 Keyword arguments: 666 lm -- log10 of the mass of dark halo 667 ''' 668 j = locate(self.__km, len(self.__km) - 1, lm) 669 return self.__sg[j]
670
671 - def __fmassM(self, lm, z):
672 """Return the mass function of dark halos multiplied by Mass - 673 Sheth e Tormen (MNRAS 308, 119, 1999). 674 """ 675 kmsgm = lm 676 kmass = 10.0 ** (kmsgm) 677 frst = self.massFunction(lm, z) * kmass 678 kmassa2 = self.__tilt2 * kmass 679 mdn_dm = kmassa2 * frst 680 return mdn_dm
681
682 - def halos_n(self, z):
683 """Return the integral of the mass function of dark halos multiplied 684 by mass in the range of log(M_min) a log(M_max) 685 686 Keyword arguments: 687 z -- redshift 688 """ 689 690 fmassM = lambda lm: self.__fmassM(lm, z) 691 692 deltal = (self.__lmax - self.__lmin) / 50.0 693 694 Lm = [self.__lmin + i * deltal for i in range(50)] 695 696 Fm = [fmassM(lm) for lm in Lm] 697 698 tck = spint.splrep(Lm, Fm) 699 Inte = spint.splint(self.__lmin, self.__lmax, tck) 700 return Inte
701
702 - def fbstruc(self, z):
703 """Return the faction of barions into structures 704 705 Keyword arguments: 706 z -- redshift 707 """ 708 rdm, drdm_dt = self._cosmology.rodm(z) 709 fb = self.halos_n(z) / rdm 710 return fb
711
712 - def numerical_density_halos(self, z):
713 """Return the numerial density of dark halos 714 within the comove volume 715 716 Keyword arguments: 717 z- redshift 718 """ 719 720 deltal = (self.__lmax - self.__lmin) / 50.0 721 Lm = [self.__lmin + i * deltal for i in range(50)] 722 723 Fm = [self.massFunction(lm, z) for lm in Lm] 724 725 tck = spint.splrep(Lm, Fm) 726 Inte = spint.splint(self.__lmin, self.__lmax, tck) 727 return Inte
728
729 - def abt(self, a):
730 """Return the accretion rate of barionic matter, as 731 a function of scala factor, into strutures. 732 733 Keyword arguments: 734 a -- scala factor (1.0 / (1.0 + z)) 735 """ 736 i = locate(self._ascale, len(self._ascale) - 1, a) 737 return self._abt2[i]
738
739 - def __startBarionicAccretionRate(self):
740 741 np = 1000 742 deltaz = self._zmax / float(np) 743 744 z = [self._zmax - i * deltaz for i in range(np)] 745 z.append(0) 746 z = array(z) 747 fbt2 = array([self.fbstruc(zi) for zi in z]) 748 ascale = array([1.0 / (1.0 + zi) for zi in z]) 749 self._ascale = ascale 750 751 tck = spint.splrep(ascale, fbt2) 752 ab3 = spint.splev(ascale, tck, der=1) 753 754 def a5(z, i): 755 a = 1.0 / (1.0 + z) 756 a2 = a * a 757 a3 = -1.0 * ab3[i] * a2 758 a4 = a3 759 a5 = self._cosmology.getRobr0() * abs(a4) \ 760 / self._cosmology.dt_dz(z) 761 return a5
762 763 self._abt2 = array([a5(z[i], i) for i in range(z.size)]) 764 self._tck_ab = spint.splrep(self._ascale, self._abt2) 765
766 - def getCacheDir(self):
767 """Return True and cache name if the cache directory existe 768 and false else. 769 """ 770 if(self._cacheDir is not None): 771 return True, self._cacheDir 772 else: 773 return False
774
775 - def setDeltaHTinker(self, delta_halo):
776 if(self.__massFunctionType == "TK"): 777 self.__delta_halo = delta_halo 778 return True 779 else: 780 return False
781
782 - def getDeltaHTinker(self):
783 if(self.__massFunctionType == "TK"): 784 return self.__delta_halo 785 else: 786 raise NameError("TinkerNotDefined")
787
788 - def setQBurrFunction(self, q):
789 """ 790 Set the q value of dark haloes mass function derived from Burr 791 distribuction. 792 """ 793 self.__qBurr = q
794
795 - def getmassFunctionDict(self):
796 """ 797 Return a list with key and function of implemented dark haloes 798 mass function 799 """ 800 mydict = [] 801 for key, value in list(self.__massFunctionDict.items()): 802 mydict.append([key, value]) 803 return mydict
804
805 - def setMassFunctionDict(self, key, function):
806 """ 807 Add a new key and function in the dark haloes mass function dictionary 808 """ 809 810 self.__massFunctionDict[key] = function
811