Added peak-background split bias parameters

Author: komatsu5147

Project.toml

@@ -1,12 +1,14 @@
 name = "HaloMF"
 uuid = "0e139a23-195b-4aec-b8fc-e9aa4361bc8d"
-authors = ["Eiichiro Komatsu <eiichirokomatsu@me.com>"]
 version = "0.5.0"
+authors = ["Eiichiro Komatsu <eiichirokomatsu@me.com>"]
 
 [deps]
 Dierckx = "39dd38d3-220a-591b-8e3c-4c3a8c710a94"
+ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210"
 
 [compat]
+ForwardDiff = "1.3.2"
 julia = "1"
 
 [extras]

README.md

@@ -36,13 +36,23 @@ Note:
 
 ## Classic multiplicity functions
 
-The package als contains some of the classic halo multiplicity functions
+The package also contains some of the classic halo multiplicity functions
 - `psMF(lnν)`: [Press & Schechter, 187, 425 (1974)](http://articles.adsabs.harvard.edu/pdf/1974ApJ...187..425P)
 - `stMF(lnν)`: Equation (2) of [Sheth & Tormen, MNRAS, 329, 61 (2002)](https://academic.oup.com/mnras/article/329/1/61/1112679)
 - `jenkinsMF(lnν)`: Equation (B3) of [Jenkins et al., MNRAS, 321, 372 (2001)](https://academic.oup.com/mnras/article/321/2/372/980658)
 
 These multiplicity functions are assumed to be "universal", in the sense that they depend only on ν and do not depend explicitly on `z`. This assumption was challenged by Tinker et al. (2008), hence the explicit dependence on `z` (see `tinker08MF` and `tinker10MF`).
 
+## Halo bias parameters from the peak-background split approximation
+
+The halo bias parameters can be calculated from a halo multiplicity function using the so-called peak-background split (PBS) approximation. See [Desjacques, Jeong & Schmdit, Phys. Rept., 733, 1 (2018)](https://www.sciencedirect.com/science/article/pii/S0370157317304192).
+
+The package contains 
+- `pbsBias(lnν, MF)` and  `pbsBias1(lnν, MF)`: The linear bias parameter from the PBS.
+- `pbsBias2(lnν, MF)`: The second-order bias parameter from the PBS.
+
+Here, `MF(lnν)` is any of the halo multiplicity functions from the above list. For example, ``pbsBias(lnν, stMF)`` and ``pbsBias(lnν, x -> tinker10MF(x, z, Δm))``.
+
 ## On normalization
 
 Some of the multiplicity functions (`tinker10MF` for `z=0`, `psMF`, `stMF`) are normalized such that
@@ -55,7 +65,7 @@ This normalization is equivalent to saying that all the mass in the Universe is
 
 where ρm is the mean mass density of the (present-day) Universe. This normalization is convenient mathematically but is not necessarily physical; thus, you do not have to pay too much attention to this. It is certainly useful for checking the code when the MF is normalized.
 
-The linear bias parameters (`tinker10Bias` and `stBias`) are also normalized such that
+The linear bias parameters (`tinker10Bias`, `stBias`) are also normalized such that
 
 - ``∫_-∞^∞ dlnν tinker10Bias(lnν, Δm) * tinker10MF(lnν, 0, Δm) = 1``
 - ``∫_-∞^∞ dlnν stBias(lnν) * stMF(lnν) = 1``
@@ -66,6 +76,8 @@ The nonlinear bias parameters (`stBias2` and `stBias3`) satisfy
 - ``∫_-∞^∞ dlnν stBias2(lnν) * stMF(lnν) = 0``
 - ``∫_-∞^∞ dlnν stBias3(lnν) * stMF(lnν) = 0``
 
+The halo bias parameters from the peak-background split approximation also satisfy the above relationships by construction.
+
 ## Relation to the halo mass function
 
 The halo mass function, dn/dM, can be computed from the halo multiplicity function, `MF`, in the following way.

src/HaloMF.jl

@@ -1,13 +1,15 @@
 module HaloMF
-using Dierckx
+using Dierckx, ForwardDiff
 export tinker08MF, tinker10MF
 export bocquetMFhy, bocquetMFdm
 export psMF, stMF, jenkinsMF
 export stBias, stBias1, stBias2, stBias3
 export tinker10Bias
+export pbsBias, pbsBias1, pbsBias2
 include("tinkerMF.jl")
 include("bocquetMF.jl")
 include("classicMF.jl")
 include("stBias.jl")
 include("tinkerBias.jl")
+include("pbsBias.jl")
 end # module

src/pbsBias.jl

@@ -0,0 +1,50 @@
+"""
+    pbsBias(lnν, MF)
+    pbsBias1(lnν, MF)
+
+Linear halo bias, b1(lnν), from the peak-background split (PBS) approximation.
+
+*Reference*: Section 3.3 of Desjacques, Jeong & Schmidt, Phys. Rept., 733, 1 (2018)
+
+# Arguments
+- `lnν::Real`: natural logarithm of a threshold, ν, i.e., `lnν` = log(ν), defined by ν ≡ [δc/σ(R,z)]^2. Here, δc = 1.6865 and σ(R,z) is the r.m.s. mass fluctuation within
+    a top-hat smoothing of scale R at a redshift `z`.
+- `MF`(lnν): a function which returns a halo multiplicity function with the argument lnν.
+
+The linear PBS bias is normalized as
+
+``∫_-∞^∞ dlnν MF(lnν) pbsBias(MF) = 1``
+"""
+function pbsBias(lnν::Real, MF)
+    δc = 1.6865
+    b1 = 1 - ForwardDiff.derivative(MF, lnν) * 2 / δc / MF(lnν)
+    return b1
+end
+const pbsBias1 = pbsBias
+
+"""
+    pbsBias2(lnν, MF)
+
+Second-order PBS halo bias, b2(lnν), from the peak-background split (PBS) approximation.
+
+*Reference*: Section 3.3 of Desjacques, Jeong & Schmidt, Phys. Rept., 733, 1 (2018)
+
+# Arguments
+- `lnν::Real`: natural logarithm of a threshold, ν (see `pbsBias`).
+- `MF`(lnν): a function which returns a halo multiplicity function with the argument lnν.
+
+The second-order PBS bias satisfies
+
+``∫_-∞^∞ dlnν MF(lnν) pbsBias2(MF) = 0``
+"""
+function pbsBias2(lnν::Real, MF)
+    δc = 1.6865
+    a2 = -17 / 21
+    dMFdlnν(x) = ForwardDiff.derivative(MF, x)
+    b2 = (
+        -2 * (1 + a2) * dMFdlnν(lnν) * 2 / δc / MF(lnν) +
+        ForwardDiff.derivative(dMFdlnν, lnν) * 4 / δc^2 / MF(lnν) -
+        dMFdlnν(lnν) * 2 / δc^2 / MF(lnν)
+    )
+    return b2
+end