diff --git a/Project.toml b/Project.toml index 979b3743..3b451311 100644 --- a/Project.toml +++ b/Project.toml @@ -5,12 +5,14 @@ version = "1.6.1" [deps] ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" LogExpFunctions = "2ab3a3ac-af41-5b50-aa03-7779005ae688" +OpenLibm_jll = "05823500-19ac-5b8b-9628-191a04bc5112" OpenSpecFun_jll = "efe28fd5-8261-553b-a9e1-b2916fc3738e" [compat] ChainRulesCore = "0.9.44, 0.10, 1" ChainRulesTestUtils = "0.6.8, 0.7, 1" LogExpFunctions = "0.2, 0.3" +OpenLibm_jll = "0.7, 0.8" OpenSpecFun_jll = "0.5" julia = "1.3" diff --git a/src/SpecialFunctions.jl b/src/SpecialFunctions.jl index 3de7ec81..c199ab84 100644 --- a/src/SpecialFunctions.jl +++ b/src/SpecialFunctions.jl @@ -3,6 +3,7 @@ module SpecialFunctions import ChainRulesCore import LogExpFunctions +using OpenLibm_jll using OpenSpecFun_jll export diff --git a/src/bessel.jl b/src/bessel.jl index ba0afd2a..e398b9da 100644 --- a/src/bessel.jl +++ b/src/bessel.jl @@ -189,14 +189,14 @@ for jy in ("j","y"), nu in (0,1) bjynu = Symbol("bessel",jy,nu) if jy == "y" @eval begin - $bjynu(x::Float64) = nan_dom_err(ccall(($jynu,libm), Float64, (Float64,), x), x) - $bjynu(x::Float32) = nan_dom_err(ccall(($jynuf,libm), Float32, (Float32,), x), x) + $bjynu(x::Float64) = nan_dom_err(ccall(($jynu,libopenlibm), Float64, (Float64,), x), x) + $bjynu(x::Float32) = nan_dom_err(ccall(($jynuf,libopenlibm), Float32, (Float32,), x), x) $bjynu(x::Float16) = Float16($bjynu(Float32(x))) end else @eval begin - $bjynu(x::Float64) = ccall(($jynu,libm), Float64, (Float64,), x) - $bjynu(x::Float32) = ccall(($jynuf,libm), Float32, (Float32,), x) + $bjynu(x::Float64) = ccall(($jynu,libopenlibm), Float64, (Float64,), x) + $bjynu(x::Float32) = ccall(($jynuf,libopenlibm), Float32, (Float32,), x) $bjynu(x::Float16) = Float16($bjynu(Float32(x))) end end @@ -389,8 +389,8 @@ function besselj(nu::Float64, z::Complex{Float64}) end end -besselj(nu::Cint, x::Float64) = ccall((:jn, libm), Float64, (Cint, Float64), nu, x) -besselj(nu::Cint, x::Float32) = ccall((:jnf, libm), Float32, (Cint, Float32), nu, x) +besselj(nu::Cint, x::Float64) = ccall((:jn, libopenlibm), Float64, (Cint, Float64), nu, x) +besselj(nu::Cint, x::Float32) = ccall((:jnf, libopenlibm), Float32, (Cint, Float32), nu, x) function besseljx(nu::Float64, z::Complex{Float64}) @@ -413,13 +413,13 @@ function bessely(nu::Cint, x::Float64) if x < 0 throw(DomainError(x, "`x` must be nonnegative.")) end - ccall((:yn, libm), Float64, (Cint, Float64), nu, x) + ccall((:yn, libopenlibm), Float64, (Cint, Float64), nu, x) end function bessely(nu::Cint, x::Float32) if x < 0 throw(DomainError(x, "`x` must be nonnegative.")) end - ccall((:ynf, libm), Float32, (Cint, Float32), nu, x) + ccall((:ynf, libopenlibm), Float32, (Cint, Float32), nu, x) end function bessely(nu::Float64, z::Complex{Float64}) diff --git a/src/erf.jl b/src/erf.jl index f5e0f155..1b44c0e3 100644 --- a/src/erf.jl +++ b/src/erf.jl @@ -1,12 +1,12 @@ # This file contains code that was formerly a part of Julia. License is MIT: http://julialang.org/license -using Base.Math: @horner, libm +using Base.Math: @horner using Base.MPFR: ROUNDING_MODE for f in (:erf, :erfc) @eval begin - ($f)(x::Float64) = ccall(($(string(f)),libm), Float64, (Float64,), x) - ($f)(x::Float32) = ccall(($(string(f,"f")),libm), Float32, (Float32,), x) + ($f)(x::Float64) = ccall(($(string(f)),libopenlibm), Float64, (Float64,), x) + ($f)(x::Float32) = ccall(($(string(f,"f")),libopenlibm), Float32, (Float32,), x) ($f)(x::Real) = ($f)(float(x)) ($f)(a::Float16) = Float16($f(Float32(a))) ($f)(a::Complex{Float16}) = Complex{Float16}($f(Complex{Float32}(a))) diff --git a/src/gamma.jl b/src/gamma.jl index f75a972a..04507fa6 100644 --- a/src/gamma.jl +++ b/src/gamma.jl @@ -570,8 +570,8 @@ export gamma, loggamma, logabsgamma, beta, logbeta, logabsbeta, logfactorial, lo ## from base/special/gamma.jl -gamma(x::Float64) = nan_dom_err(ccall((:tgamma, libm), Float64, (Float64,), x), x) -gamma(x::Float32) = nan_dom_err(ccall((:tgammaf, libm), Float32, (Float32,), x), x) +gamma(x::Float64) = nan_dom_err(ccall((:tgamma, libopenlibm), Float64, (Float64,), x), x) +gamma(x::Float32) = nan_dom_err(ccall((:tgammaf, libopenlibm), Float32, (Float32,), x), x) gamma(x::Float16) = Float16(gamma(Float32(x))) gamma(x::AbstractFloat) = throw(MethodError(gamma, x)) @@ -619,12 +619,12 @@ gamma(x::Number) = gamma(float(x)) function logabsgamma(x::Float64) signp = Ref{Int32}() - y = ccall((:lgamma_r,libm), Float64, (Float64, Ptr{Int32}), x, signp) + y = ccall((:lgamma_r,libopenlibm), Float64, (Float64, Ptr{Int32}), x, signp) return y, Int(signp[]) end function logabsgamma(x::Float32) signp = Ref{Int32}() - y = ccall((:lgammaf_r,libm), Float32, (Float32, Ptr{Int32}), x, signp) + y = ccall((:lgammaf_r,libopenlibm), Float32, (Float32, Ptr{Int32}), x, signp) return y, Int(signp[]) end logabsgamma(x::Real) = logabsgamma(float(x)) diff --git a/test/gamma.jl b/test/gamma.jl index 26f287f6..6ff346d1 100644 --- a/test/gamma.jl +++ b/test/gamma.jl @@ -59,11 +59,7 @@ @test eta(Complex{Float32}(2)) ≈ eta(2) end @testset "gamma, loggamma, logabsgamma (complex argument)" begin - if Base.Math.libm == "libopenlibm" - @test gamma.(Float64[1:25;]) == gamma.(1:25) - else - @test gamma.(Float64[1:25;]) ≈ gamma.(1:25) - end + @test gamma.(Float64[1:25;]) ≈ gamma.(1:25) for elty in (Float32, Float64) @test gamma(convert(elty,1/2)) ≈ convert(elty,sqrt(π)) @test gamma(convert(elty,-1/2)) ≈ convert(elty,-2sqrt(π))