From 9b84b16c3f17ef31ec74e219c371ed6b4513ec03 Mon Sep 17 00:00:00 2001 From: Alex Arslan Date: Thu, 14 Dec 2017 17:45:19 -0800 Subject: [PATCH] Fix some 0.7 deprecations --- src/bessel.jl | 56 +++++++++++++++++++++++------------------------ src/deprecated.jl | 8 +++---- src/erf.jl | 8 +++---- src/gamma.jl | 6 ++--- test/runtests.jl | 28 ++++++++++++------------ 5 files changed, 53 insertions(+), 53 deletions(-) diff --git a/src/bessel.jl b/src/bessel.jl index 343d6f31..15261cea 100644 --- a/src/bessel.jl +++ b/src/bessel.jl @@ -26,7 +26,7 @@ function Base.showerror(io::IO, ex::AmosException) end ## Airy functions -function _airy(z::Complex128, id::Int32, kode::Int32) +function _airy(z::Complex{Float64}, id::Int32, kode::Int32) ai1, ai2 = Ref{Float64}(), Ref{Float64}() ae1, ae2 = Ref{Int32}(), Ref{Int32}() @@ -43,7 +43,7 @@ function _airy(z::Complex128, id::Int32, kode::Int32) end end -function _biry(z::Complex128, id::Int32, kode::Int32) +function _biry(z::Complex{Float64}, id::Int32, kode::Int32) ai1, ai2 = Ref{Float64}(), Ref{Float64}() ae1 = Ref{Int32}() @@ -67,7 +67,7 @@ end Airy function of the first kind ``\\operatorname{Ai}(x)``. """ function airyai end -airyai(z::Complex128) = _airy(z, Int32(0), Int32(1)) +airyai(z::Complex{Float64}) = _airy(z, Int32(0), Int32(1)) """ airyaiprime(x) @@ -75,7 +75,7 @@ airyai(z::Complex128) = _airy(z, Int32(0), Int32(1)) Derivative of the Airy function of the first kind ``\\operatorname{Ai}'(x)``. """ function airyaiprime end -airyaiprime(z::Complex128) = _airy(z, Int32(1), Int32(1)) +airyaiprime(z::Complex{Float64}) = _airy(z, Int32(1), Int32(1)) """ airybi(x) @@ -83,7 +83,7 @@ airyaiprime(z::Complex128) = _airy(z, Int32(1), Int32(1)) Airy function of the second kind ``\\operatorname{Bi}(x)``. """ function airybi end -airybi(z::Complex128) = _biry(z, Int32(0), Int32(1)) +airybi(z::Complex{Float64}) = _biry(z, Int32(0), Int32(1)) """ airybiprime(x) @@ -91,7 +91,7 @@ airybi(z::Complex128) = _biry(z, Int32(0), Int32(1)) Derivative of the Airy function of the second kind ``\\operatorname{Bi}'(x)``. """ function airybiprime end -airybiprime(z::Complex128) = _biry(z, Int32(1), Int32(1)) +airybiprime(z::Complex{Float64}) = _biry(z, Int32(1), Int32(1)) """ airyaix(x) @@ -100,7 +100,7 @@ Scaled Airy function of the first kind ``\\operatorname{Ai}(x) e^{\\frac{2}{3} x \\sqrt{x}}``. Throws `DomainError` for negative `Real` arguments. """ function airyaix end -airyaix(z::Complex128) = _airy(z, Int32(0), Int32(2)) +airyaix(z::Complex{Float64}) = _airy(z, Int32(0), Int32(2)) """ airyaiprimex(x) @@ -109,7 +109,7 @@ Scaled derivative of the Airy function of the first kind ``\\operatorname{Ai}'(x e^{\\frac{2}{3} x \\sqrt{x}}``. Throws `DomainError` for negative `Real` arguments. """ function airyaiprimex end -airyaiprimex(z::Complex128) = _airy(z, Int32(1), Int32(2)) +airyaiprimex(z::Complex{Float64}) = _airy(z, Int32(1), Int32(2)) """ airybix(x) @@ -117,7 +117,7 @@ airyaiprimex(z::Complex128) = _airy(z, Int32(1), Int32(2)) Scaled Airy function of the second kind ``\\operatorname{Bi}(x) e^{- \\left| \\operatorname{Re} \\left( \\frac{2}{3} x \\sqrt{x} \\right) \\right|}``. """ function airybix end -airybix(z::Complex128) = _biry(z, Int32(0), Int32(2)) +airybix(z::Complex{Float64}) = _biry(z, Int32(0), Int32(2)) """ airybiprimex(x) @@ -125,14 +125,14 @@ airybix(z::Complex128) = _biry(z, Int32(0), Int32(2)) Scaled derivative of the Airy function of the second kind ``\\operatorname{Bi}'(x) e^{- \\left| \\operatorname{Re} \\left( \\frac{2}{3} x \\sqrt{x} \\right) \\right|}``. """ function airybiprimex end -airybiprimex(z::Complex128) = _biry(z, Int32(1), Int32(2)) +airybiprimex(z::Complex{Float64}) = _biry(z, Int32(1), Int32(2)) for afn in (:airyai, :airyaiprime, :airybi, :airybiprime, :airyaix, :airyaiprimex, :airybix, :airybiprimex) @eval begin $afn(z::Complex) = $afn(float(z)) $afn(z::Complex{<:AbstractFloat}) = throw(MethodError($afn,(z,))) - $afn(z::Complex64) = Complex64($afn(Complex128(z))) + $afn(z::Complex{Float32}) = Complex{Float32}($afn(Complex{Float64}(z))) end if afn in (:airyaix, :airyaiprimex) @eval $afn(x::Real) = x < 0 ? throw(DomainError(x, "`x` must be nonnegative.")) : real($afn(complex(float(x)))) @@ -172,7 +172,7 @@ for jy in ("j","y"), nu in (0,1) end -function _besselh(nu::Float64, k::Int32, z::Complex128, kode::Int32) +function _besselh(nu::Float64, k::Int32, z::Complex{Float64}, kode::Int32) ai1, ai2 = Ref{Float64}(), Ref{Float64}() ae1, ae2 = Ref{Int32}(), Ref{Int32}() @@ -189,7 +189,7 @@ function _besselh(nu::Float64, k::Int32, z::Complex128, kode::Int32) end end -function _besseli(nu::Float64, z::Complex128, kode::Int32) +function _besseli(nu::Float64, z::Complex{Float64}, kode::Int32) ai1, ai2 = Ref{Float64}(), Ref{Float64}() ae1, ae2 = Ref{Int32}(), Ref{Int32}() @@ -206,7 +206,7 @@ function _besseli(nu::Float64, z::Complex128, kode::Int32) end end -function _besselj(nu::Float64, z::Complex128, kode::Int32) +function _besselj(nu::Float64, z::Complex{Float64}, kode::Int32) ai1, ai2 = Ref{Float64}(), Ref{Float64}() ae1, ae2 = Ref{Int32}(), Ref{Int32}() @@ -223,7 +223,7 @@ function _besselj(nu::Float64, z::Complex128, kode::Int32) end end -function _besselk(nu::Float64, z::Complex128, kode::Int32) +function _besselk(nu::Float64, z::Complex{Float64}, kode::Int32) ai1, ai2 = Ref{Float64}(), Ref{Float64}() ae1, ae2 = Ref{Int32}(), Ref{Int32}() @@ -240,7 +240,7 @@ function _besselk(nu::Float64, z::Complex128, kode::Int32) end end -function _bessely(nu::Float64, z::Complex128, kode::Int32) +function _bessely(nu::Float64, z::Complex{Float64}, kode::Int32) ai1, ai2 = Ref{Float64}(), Ref{Float64}() ae1, ae2 = Ref{Int32}(), Ref{Int32}() wrk1, wrk2 = Ref{Float64}(), Ref{Float64}() @@ -268,7 +268,7 @@ selecting [`hankelh1`](@ref) or [`hankelh2`](@ref), respectively. """ function besselh end -function besselh(nu::Float64, k::Integer, z::Complex128) +function besselh(nu::Float64, k::Integer, z::Complex{Float64}) if nu < 0 s = (k == 1) ? 1 : -1 return _besselh(-nu,Int32(k),z,Int32(1)) * complex(cospi(nu),-s*sinpi(nu)) @@ -291,7 +291,7 @@ exponential factor (analytically), so it avoids these problems. """ function besselhx end -function besselhx(nu::Float64, k::Integer, z::Complex128) +function besselhx(nu::Float64, k::Integer, z::Complex{Float64}) if nu < 0 s = (k == 1) ? 1 : -1 return _besselh(-nu,Int32(k),z,Int32(2)) * complex(cospi(nu),-s*sinpi(nu)) @@ -299,7 +299,7 @@ function besselhx(nu::Float64, k::Integer, z::Complex128) return _besselh(nu,Int32(k),z,Int32(2)) end -function besseli(nu::Float64, z::Complex128) +function besseli(nu::Float64, z::Complex{Float64}) if nu < 0 if isinteger(nu) return _besseli(-nu,z,Int32(1)) @@ -311,7 +311,7 @@ function besseli(nu::Float64, z::Complex128) end end -function besselix(nu::Float64, z::Complex128) +function besselix(nu::Float64, z::Complex{Float64}) if nu < 0 if isinteger(nu) return _besseli(-nu,z,Int32(2)) @@ -323,7 +323,7 @@ function besselix(nu::Float64, z::Complex128) end end -function besselj(nu::Float64, z::Complex128) +function besselj(nu::Float64, z::Complex{Float64}) if nu < 0 if isinteger(nu) return _besselj(-nu,z,Int32(1))*cospi(nu) @@ -339,7 +339,7 @@ besselj(nu::Cint, x::Float64) = ccall((:jn, libm), Float64, (Cint, Float64), nu, besselj(nu::Cint, x::Float32) = ccall((:jnf, libm), Float32, (Cint, Float32), nu, x) -function besseljx(nu::Float64, z::Complex128) +function besseljx(nu::Float64, z::Complex{Float64}) if nu < 0 if isinteger(nu) return _besselj(-nu,z,Int32(2))*cospi(nu) @@ -351,9 +351,9 @@ function besseljx(nu::Float64, z::Complex128) end end -besselk(nu::Float64, z::Complex128) = _besselk(abs(nu), z, Int32(1)) +besselk(nu::Float64, z::Complex{Float64}) = _besselk(abs(nu), z, Int32(1)) -besselkx(nu::Float64, z::Complex128) = _besselk(abs(nu), z, Int32(2)) +besselkx(nu::Float64, z::Complex{Float64}) = _besselk(abs(nu), z, Int32(2)) function bessely(nu::Cint, x::Float64) if x < 0 @@ -368,7 +368,7 @@ function bessely(nu::Cint, x::Float32) ccall((:ynf, libm), Float32, (Cint, Float32), nu, x) end -function bessely(nu::Float64, z::Complex128) +function bessely(nu::Float64, z::Complex{Float64}) if nu < 0 return _bessely(-nu,z,Int32(1))*cospi(nu) - _besselj(-nu,z,Int32(1))*sinpi(nu) else @@ -376,7 +376,7 @@ function bessely(nu::Float64, z::Complex128) end end -function besselyx(nu::Float64, z::Complex128) +function besselyx(nu::Float64, z::Complex{Float64}) if nu < 0 return _bessely(-nu,z,Int32(2))*cospi(nu) - _besselj(-nu,z,Int32(2))*sinpi(nu) else @@ -500,7 +500,7 @@ for f in ("i", "ix", "j", "jx", "k", "kx", "y", "yx") $bfn(Tf(nu), Complex{Tf}(z)) end $bfn(k::T, z::Complex{T}) where {T<:AbstractFloat} = throw(MethodError($bfn,(k,z))) - $bfn(nu::Float32, x::Complex64) = Complex64($bfn(Float64(nu), Complex128(x))) + $bfn(nu::Float32, x::Complex{Float32}) = Complex{Float32}($bfn(Float64(nu), Complex{Float64}(x))) end end @@ -517,7 +517,7 @@ for bfn in (:besselh, :besselhx) end $bfn(nu::T, k::Integer, z::Complex{T}) where {T<:AbstractFloat} = throw(MethodError($bfn,(nu,k,z))) - $bfn(nu::Float32, k::Integer, x::Complex64) = Complex64($bfn(Float64(nu), k, Complex128(x))) + $bfn(nu::Float32, k::Integer, x::Complex{Float32}) = Complex{Float32}($bfn(Float64(nu), k, Complex{Float64}(x))) end end diff --git a/src/deprecated.jl b/src/deprecated.jl index d0e0d368..c8853cad 100644 --- a/src/deprecated.jl +++ b/src/deprecated.jl @@ -4,7 +4,7 @@ @deprecate airyx(z::Number) airyaix(z) @deprecate airyprime(z::Number) airyaiprime(z) -function _airy(k::Integer, z::Complex128) +function _airy(k::Integer, z::Complex{Float64}) Base.depwarn("`airy(k,x)` is deprecated, use `airyai(x)`, `airyaiprime(x)`, `airybi(x)` or `airybiprime(x)` instead.",:airy) id = Int32(k==1 || k==3) if k == 0 || k == 1 @@ -15,7 +15,7 @@ function _airy(k::Integer, z::Complex128) throw(DomainError(k, "`k` must be between 0 and 3.")) end end -function _airyx(k::Integer, z::Complex128) +function _airyx(k::Integer, z::Complex{Float64}) Base.depwarn("`airyx(k,x)` is deprecated, use `airyaix(x)`, `airyaiprimex(x)`, `airybix(x)` or `airybiprimex(x)` instead.",:airyx) id = Int32(k==1 || k==3) if k == 0 || k == 1 @@ -31,7 +31,7 @@ for afn in (:airy,:airyx) _afn = Symbol("_"*string(afn)) suf = string(afn)[5:end] @eval begin - function $afn(k::Integer, z::Complex128) + function $afn(k::Integer, z::Complex{Float64}) afn = $(QuoteNode(afn)) suf = $(QuoteNode(suf)) Base.depwarn("`$afn(k,x)` is deprecated, use `airyai$suf(x)`, `airyaiprime$suf(x)`, `airybi$suf(x)` or `airybiprime$suf(x)` instead.",$(QuoteNode(afn))) @@ -40,7 +40,7 @@ for afn in (:airy,:airyx) $afn(k::Integer, z::Complex) = $afn(k, float(z)) $afn(k::Integer, z::Complex{<:AbstractFloat}) = throw(MethodError($afn,(k,z))) - $afn(k::Integer, z::Complex64) = Complex64($afn(k, Complex128(z))) + $afn(k::Integer, z::Complex{Float32}) = Complex{Float32}($afn(k, Complex{Float64}(z))) $afn(k::Integer, x::Real) = $afn(k, float(x)) $afn(k::Integer, x::AbstractFloat) = real($afn(k, complex(x))) end diff --git a/src/erf.jl b/src/erf.jl index fc9e17c4..22ba3c57 100644 --- a/src/erf.jl +++ b/src/erf.jl @@ -9,7 +9,7 @@ for f in (:erf, :erfc) ($f)(x::Float32) = ccall(($(string(f,"f")),libm), Float32, (Float32,), x) ($f)(x::Real) = ($f)(float(x)) ($f)(a::Float16) = Float16($f(Float32(a))) - ($f)(a::Complex32) = Complex32($f(Complex64(a))) + ($f)(a::Complex{Float16}) = Complex{Float16}($f(Complex{Float32}(a))) function ($f)(x::BigFloat) z = BigFloat() ccall(($(string(:mpfr_,f)), :libmpfr), Int32, (Ref{BigFloat}, Ref{BigFloat}, Int32), z, x, ROUNDING_MODE[]) @@ -22,9 +22,9 @@ end for f in (:erf, :erfc, :erfcx, :erfi, :Dawson) fname = (f === :Dawson) ? :dawson : f @eval begin - ($fname)(z::Complex128) = Complex128(ccall(($(string("Faddeeva_",f)),openspecfun), Complex{Float64}, (Complex{Float64}, Float64), z, zero(Float64))) - ($fname)(z::Complex64) = Complex64(ccall(($(string("Faddeeva_",f)),openspecfun), Complex{Float64}, (Complex{Float64}, Float64), Complex128(z), Float64(eps(Float32)))) - ($fname)(z::Complex) = ($fname)(Complex128(z)) + ($fname)(z::Complex{Float64}) = Complex{Float64}(ccall(($(string("Faddeeva_",f)),openspecfun), Complex{Float64}, (Complex{Float64}, Float64), z, zero(Float64))) + ($fname)(z::Complex{Float32}) = Complex{Float32}(ccall(($(string("Faddeeva_",f)),openspecfun), Complex{Float64}, (Complex{Float64}, Float64), Complex{Float64}(z), Float64(eps(Float32)))) + ($fname)(z::Complex) = ($fname)(Complex{Float64}(z)) end end diff --git a/src/gamma.jl b/src/gamma.jl index 1301bdbb..932eafda 100644 --- a/src/gamma.jl +++ b/src/gamma.jl @@ -100,13 +100,13 @@ function cotderiv_q(m::Int) q₋ = cotderiv_q(m-1) d = length(q₋) - 1 # degree of q₋ if isodd(m-1) - q = Array{Float64}(length(q₋)) + q = Vector{Float64}(uninitialized, length(q₋)) q[end] = d * q₋[end] * 2/m for i = 1:length(q)-1 q[i] = ((i-1)*q₋[i] + i*q₋[i+1]) * 2/m end else # iseven(m-1) - q = Array{Float64}(length(q₋) + 1) + q = Vector{Float64}(uninitialized, length(q₋) + 1) q[1] = q₋[1] / m q[end] = (1 + 2d) * q₋[end] / m for i = 2:length(q)-1 @@ -486,7 +486,7 @@ end for T in (Float16, Float32, Float64) - @eval f64(x::Complex{$T}) = Complex128(x) + @eval f64(x::Complex{$T}) = Complex{Float64}(x) @eval f64(x::$T) = Float64(x) end diff --git a/test/runtests.jl b/test/runtests.jl index 276f09c0..65932f1d 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -27,7 +27,7 @@ relerrc(z, x) = max(relerr(real(z),real(x)), relerr(imag(z),imag(x))) @test SF.erfc(Float16(1)) ≈ 0.15729920705028513066 @test SF.erfcx(1) ≈ 0.42758357615580700442 @test SF.erfcx(Float32(1)) ≈ 0.42758357615580700442 - @test SF.erfcx(Complex64(1)) ≈ 0.42758357615580700442 + @test SF.erfcx(Complex{Float32}(1)) ≈ 0.42758357615580700442 @test SF.erfi(1) ≈ 1.6504257587975428760 @test SF.erfinv(0.84270079294971486934) ≈ 1 @test SF.erfcinv(0.15729920705028513066) ≈ 1 @@ -104,13 +104,13 @@ end @test_throws SF.AmosException SF.airyai(200im) @test_throws SF.AmosException SF.airybi(200) - for T in [Float32, Float64, Complex64,Complex128] + for T in [Float32, Float64, Complex{Float32},Complex{Float64}] @test SF.airyai(T(1.8)) ≈ 0.0470362168668458052247 @test SF.airyaiprime(T(1.8)) ≈ -0.0685247801186109345638 @test SF.airybi(T(1.8)) ≈ 2.595869356743906290060 @test SF.airybiprime(T(1.8)) ≈ 2.98554005084659907283 end - for T in [Complex64, Complex128] + for T in [Complex{Float32}, Complex{Float64}] z = convert(T,1.8 + 1.0im) @test SF.airyaix(z) ≈ SF.airyai(z) * exp(2/3 * z * sqrt(z)) @test SF.airyaiprimex(z) ≈ SF.airyaiprime(z) * exp(2/3 * z * sqrt(z)) @@ -172,13 +172,13 @@ end @test SF.besseli(-3,3) ≈ true_i33 @test SF.besseli(3,-3) ≈ -true_i33 @test SF.besseli(-3,-3) ≈ -true_i33 - @test SF.besseli(Float32(-3),Complex64(-3,0)) ≈ -true_i33 + @test SF.besseli(Float32(-3),Complex{Float32}(-3,0)) ≈ -true_i33 true_im3p1_3 = 0.84371226532586351965 @test SF.besseli(-3.1,3) ≈ true_im3p1_3 for i in [-5 -3 -1 1 3 5] @test SF.besseli(i,0) == 0.0 @test SF.besseli(i,Float32(0)) == 0 - @test SF.besseli(i,Complex64(0)) == 0 + @test SF.besseli(i,Complex{Float32}(0)) == 0 end @testset "Error throwing" begin @test_throws SF.AmosException SF.besseli(1,1000) @@ -194,7 +194,7 @@ end for i in [-5 -3 -1 1 3 5] @test SF.besselj(i,0) == 0 @test SF.besselj(i,Float32(0)) == 0 - @test SF.besselj(i,Complex64(0)) == 0 + @test SF.besselj(i,Complex{Float32}(0)) == 0 end j33 = SF.besselj(3,3.) @@ -280,7 +280,7 @@ end end @testset "besselhx" begin - for elty in [Complex64,Complex128] + for elty in [Complex{Float32},Complex{Float64}] z = convert(elty, 1.0 + 1.9im) @test SF.besselhx(1.0, 1, z) ≈ convert(elty,-0.5949634147786144 - 0.18451272807835967im) @test SF.besselhx(Float32(1.0), 1, z) ≈ convert(elty,-0.5949634147786144 - 0.18451272807835967im) @@ -294,7 +294,7 @@ end end @testset "scaled bessel[ijky] and hankelh[12]" begin for x in (1.0, 0.0, -1.0), y in (1.0, 0.0, -1.0), nu in (1.0, 0.0, -1.0) - z = Complex128(x + y * im) + z = Complex{Float64}(x + y * im) z == zero(z) || @test SF.hankelh1x(nu, z) ≈ SF.hankelh1(nu, z) * exp(-z * im) z == zero(z) || @test SF.hankelh2x(nu, z) ≈ SF.hankelh2(nu, z) * exp(z * im) (nu < 0 && z == zero(z)) || @test SF.besselix(nu, z) ≈ SF.besseli(nu, z) * exp(-abs(real(z))) @@ -308,8 +308,8 @@ end @test SF.besselix(i,0) == 0 @test SF.besseljx(i,Float32(0)) == 0 @test SF.besselix(i,Float32(0)) == 0 - @test SF.besseljx(i,Complex64(0)) == 0 - @test SF.besselix(i,Complex64(0)) == 0 + @test SF.besseljx(i,Complex{Float32}(0)) == 0 + @test SF.besselix(i,Complex{Float32}(0)) == 0 end @testset "Error throwing" begin @test_throws SF.AmosException SF.hankelh1x(1, 0) @@ -325,8 +325,8 @@ end end @testset "issue #6653" begin @testset "$f" for f in (SF.besselj,SF.bessely,SF.besseli,SF.besselk,SF.hankelh1,SF.hankelh2) - @test f(0,1) ≈ f(0,Complex128(1)) - @test f(0,1) ≈ f(0,Complex64(1)) + @test f(0,1) ≈ f(0,Complex{Float64}(1)) + @test f(0,1) ≈ f(0,Complex{Float32}(1)) end end end @@ -389,14 +389,14 @@ end @test SF.eta(1) ≈ log(2) @test SF.eta(2) ≈ pi^2/12 @test SF.eta(Float32(2)) ≈ SF.eta(2) - @test SF.eta(Complex64(2)) ≈ SF.eta(2) + @test SF.eta(Complex{Float32}(2)) ≈ SF.eta(2) end end @testset "zeta" begin @test SF.zeta(0) ≈ -0.5 @test SF.zeta(2) ≈ pi^2/6 - @test SF.zeta(Complex64(2)) ≈ SF.zeta(2) + @test SF.zeta(Complex{Float32}(2)) ≈ SF.zeta(2) @test SF.zeta(4) ≈ pi^4/90 @test SF.zeta(1,Float16(2.)) ≈ SF.zeta(1,2.) @test SF.zeta(1.,Float16(2.)) ≈ SF.zeta(1,2.)