1
2
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
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
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
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
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
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
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
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
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
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
303
304
305
306
307
308
309
310
311
312
313
315 return self._cosmology.rodm(z)[0]
316
319
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
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
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
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
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
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
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
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
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
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
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
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
535
536
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
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
601
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
612
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
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
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
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
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
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
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
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
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
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
776 if(self.__massFunctionType == "TK"):
777 self.__delta_halo = delta_halo
778 return True
779 else:
780 return False
781
783 if(self.__massFunctionType == "TK"):
784 return self.__delta_halo
785 else:
786 raise NameError("TinkerNotDefined")
787
789 """
790 Set the q value of dark haloes mass function derived from Burr
791 distribuction.
792 """
793 self.__qBurr = q
794
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
806 """
807 Add a new key and function in the dark haloes mass function dictionary
808 """
809
810 self.__massFunctionDict[key] = function
811