math-functions 0.2.0.0 → 0.2.0.1
raw patch · 4 files changed
+16/−13 lines, 4 filesPVP ok
version bump matches the API change (PVP)
API changes (from Hackage documentation)
Files
- Numeric/SpecFunctions/Internal.hs +7/−7
- changelog.md +5/−0
- math-functions.cabal +1/−1
- tests/Tests/SpecFunctions.hs +3/−5
Numeric/SpecFunctions/Internal.hs view
@@ -352,7 +352,7 @@ -- | Compute the natural logarithm of the beta function. -- -- \[--- B(a,b) = \int_0^1 t^{a-1}(1-t)^{1-b}\,dt = \frac{\Gamma{a}\Gamma{b}}{\Gamma{a+b}}+-- B(a,b) = \int_0^1 t^{a-1}(1-t)^{1-b}\,dt = \frac{\Gamma(a)\Gamma(b)}{\Gamma(a+b)} -- \] logBeta :: Double -- ^ /a/ > 0@@ -490,12 +490,10 @@ | a < 0 || a > 1 = modErr $ printf "invIncompleteBeta x must be in [0,1]. p=%g q=%g a=%g" p q a | a == 0 || a == 1 = a- | a > 0.5 = 1 - invIncompleteBetaWorker (logBeta p q) q p (1 - a)- | otherwise = invIncompleteBetaWorker (logBeta p q) p q a+ | otherwise = invIncompleteBetaWorker (logBeta p q) p q a invIncompleteBetaWorker :: Double -> Double -> Double -> Double -> Double--- NOTE: p <= 0.5. invIncompleteBetaWorker beta a b p = loop (0::Int) (invIncBetaGuess beta a b p) where a1 = a - 1@@ -562,14 +560,14 @@ in x -- If both a and b larger or equal that 1 but not too big we use -- same approximation as above but calculate it a bit differently- | a+b <= 6 && a>=1 && b>=1 =+ | a+b <= 6 && a>1 && b>1 = let x_infl = (a - 1) / (a + b - 2) p_infl = incompleteBeta a b x_infl x | p < p_infl = exp ((log(p * a) + beta) / a) | otherwise = 1 - exp((log((1-p) * b) + beta) / b) in x -- For small a and not too big b we use approximation from boost.- | b < 5 && a < 1 =+ | b < 5 && a <= 1 = let x | p**(1/a) < 0.5 = (p * a * exp beta) ** (1/a) | otherwise = 1 - (1 - p ** (b * exp beta))**(1/b) in x@@ -829,7 +827,9 @@ -- | Calculate the error term of the Stirling approximation. This is -- only defined for non-negative values. ----- > stirlingError @n@ = @log(n!) - log(sqrt(2*pi*n)*(n/e)^n)+-- \[+-- \operatorname{stirlingError}(n) = \log(n!) - \log(\sqrt{2\pi n}\frac{n}{e}^n)+-- \] stirlingError :: Double -> Double stirlingError n | n <= 15.0 = case properFraction (n+n) of
changelog.md view
@@ -1,5 +1,10 @@ Changes in 0.2.0.0 + * Bug fixes and documentation tweaks+++Changes in 0.2.0.0+ * `logGamma` now uses Lancsoz approximation and same as `logGammaL`. Old implementation of `logGamma` moved to `Numeric.SpecFunctions.Extra.logGammaAS245`.
math-functions.cabal view
@@ -1,5 +1,5 @@ name: math-functions-version: 0.2.0.0+version: 0.2.0.1 cabal-version: >= 1.10 license: BSD3 license-file: LICENSE
tests/Tests/SpecFunctions.hs view
@@ -50,9 +50,8 @@ $ and [ eq 1e-15 (logGammaL (fromIntegral n)) (logFactorial (n-1)) | n <- [3..10000::Int]]- -- FIXME: Too low! , testAssertion "logGammaL is expected to be precise at 1e-10 level [fractional points]"- $ and [ eq 1e-10 (logGammaL x) lg | (x,lg) <- tableLogGamma ]+ $ and [ eq (64*m_epsilon) (logGammaL x) lg | (x,lg) <- tableLogGamma ] -- FIXME: loss of precision when logBeta p q ≈ 0. -- Relative error doesn't work properly in this case. , testAssertion "logBeta is expected to be precise at 1e-6 level"@@ -66,9 +65,8 @@ -- Relative precision is lost when digamma(x) ≈ 0 , testAssertion "digamma is expected to be precise at 1e-12" $ and [ eq 1e-12 r (digamma x) | (x,r) <- tableDigamma ]- -- FIXME: Why 1e-8? Is it due to poor precision of logBeta?- , testAssertion "incompleteBeta is expected to be precise at 1e-8 level"- $ and [ eq 1e-8 (incompleteBeta p q x) ib | (p,q,x,ib) <- tableIncompleteBeta ]+ , testAssertion "incompleteBeta is expected to be precise at 32*m_epsilon level"+ $ and [ eq (32 * m_epsilon) (incompleteBeta p q x) ib | (p,q,x,ib) <- tableIncompleteBeta ] , testAssertion "incompleteBeta with p > 3000 and q > 3000" $ and [ eq 1e-11 (incompleteBeta p q x) ib | (x,p,q,ib) <- [ (0.495, 3001, 3001, 0.2192546757957825068677527085659175689142653854877723)