@@ -13,21 +13,6 @@ export dnPl
13
13
export collectdnPl
14
14
export collectdnPl!
15
15
16
- # cache values of logfactorial to avoid expensive calls
17
- const lfdictInt = Vector {Dict{Int, Float64}} ()
18
- const lfdictBigInt = Vector {Dict{BigInt, BigFloat}} ()
19
-
20
- function __init__ ()
21
- empty! (lfdictInt)
22
- empty! (lfdictBigInt)
23
- resize! (lfdictInt, Threads. nthreads ())
24
- resize! (lfdictBigInt, Threads. nthreads ())
25
- for i in 1 : Threads. nthreads ()
26
- lfdictInt[i] = Dict {Int, Float64} ()
27
- lfdictBigInt[i] = Dict {BigInt, BigFloat} ()
28
- end
29
- end
30
-
31
16
checkdomain (x) = true
32
17
function checkdomain (x:: Number )
33
18
f = abs (x) > 1
84
69
convert (T, P_n_l)
85
70
end
86
71
87
- # Allow memoization
88
- _logfactorial (n) = logfactorial (n)
89
- _logfactorial (n:: Int ) = get! (lfdictInt[Threads. threadid ()], n) do
90
- logfactorial (n)
91
- end
92
- _logfactorial (n:: BigInt ) = get! (lfdictBigInt[Threads. threadid ()], n) do
93
- logfactorial (n)
94
- end
95
-
96
72
function logplm_norm (l, m)
97
73
T = promote_type (typeof (l), typeof (m))
98
- (log (T (2 )) - log (2 T (l)+ 1 ) + _logfactorial (l + m) - _logfactorial (l - m))/ 2
74
+ (log (T (2 )) - log (2 T (l)+ 1 ) + logfactorial (l + m) - logfactorial (l - m))/ 2
99
75
end
100
76
function _maybebigexp (t)
101
77
if t < log (prevfloat (typemax (t)))
111
87
function logabspll_prefactor (l, T = typeof (l))
112
88
lT = T (l)
113
89
a = (logfactorial (2 lT- 1 ) + log (2 lT+ 1 ) - log (lT))/ 2
114
- b = l* log (T (2 )) + _logfactorial (lT- 1 )
90
+ b = l* log (T (2 )) + logfactorial (lT- 1 )
115
91
a - b
116
92
end
117
93
neg1pow (l) = iseven (l) ? 1 : - 1
118
94
function pll_prefactor (l, T = typeof (l); csphase:: Bool = true )
119
- l == 0 && return sqrt (oftype (_logfactorial (T (l)), 1 / 2 ))
95
+ l == 0 && return sqrt (oftype (logfactorial (T (l)), 1 / 2 ))
120
96
t = logabspll_prefactor (l, T)
121
97
(csphase ? neg1pow (l) : 1 ) * exp (t)
122
98
end
0 commit comments