Fix imports of special functions (#268)

* Fix special functions from Calculus

* Fix imports of special functions

* Add test

* Avoid overwriting gamma

* Relax signatures for gamma/loggamma

Author: jishnub

Project.toml

@@ -1,6 +1,6 @@
 name = "ApproxFunBase"
 uuid = "fbd15aa5-315a-5a7d-a8a4-24992e37be05"
-version = "0.7.34"
+version = "0.7.35"
 
 [deps]
 AbstractFFTs = "621f4979-c628-5d54-868e-fcf4e3e8185c"

src/ApproxFunBase.jl

@@ -1,42 +1,44 @@
 module ApproxFunBase
 using Base: AnyDict
-using Base, BlockArrays, BandedMatrices, BlockBandedMatrices, DomainSets, IntervalSets,
-            SpecialFunctions, AbstractFFTs, FFTW, SpecialFunctions, DSP, DualNumbers,
-            LinearAlgebra, SparseArrays, LowRankApprox, FillArrays, InfiniteArrays, InfiniteLinearAlgebra #, Arpack
-import StaticArrays, Calculus
-
-import DomainSets: Domain, indomain, UnionDomain, ProductDomain, FullSpace, Point, elements, DifferenceDomain,
-            Interval, ChebyshevInterval, boundary, ∂, rightendpoint, leftendpoint,
-            dimension, WrappedDomain, VcatDomain, component, components, ncomponents
-
-
+using Base, BlockArrays, BandedMatrices, BlockBandedMatrices, DomainSets,
+              IntervalSets, SpecialFunctions, AbstractFFTs, FFTW,
+              SpecialFunctions, DSP, DualNumbers, LinearAlgebra, SparseArrays,
+              LowRankApprox, FillArrays, InfiniteArrays, InfiniteLinearAlgebra
+              # Arpack
 
-import AbstractFFTs: Plan, fft, ifft
-import FFTW: plan_r2r!, fftwNumber, REDFT10, REDFT01, REDFT00, RODFT00, R2HC, HC2R,
-                r2r!, r2r,  plan_fft, plan_ifft, plan_ifft!, plan_fft!
-
-
-import Base: values, convert, getindex, setindex!, *, +, -, ==, <, <=, >, |, !, !=, eltype, iterate,
-                >=, /, ^, \, ∪, transpose, size, tail, broadcast, broadcast!, copyto!, copy, to_index, (:),
-                similar, map, vcat, hcat, hvcat, show, summary, stride, sum, cumsum, sign, imag, conj, inv,
-                complex, reverse, exp, sqrt, abs, abs2, sign, issubset, values, in, first, last, rand, intersect, setdiff,
-                isless, union, angle, join, isnan, isapprox, isempty, sort, merge, promote_rule,
-                minimum, maximum, extrema, argmax, argmin, findmax, findmin, isfinite,
-                zeros, zero, one, promote_rule, repeat, length, resize!, isinf,
-                getproperty, findfirst, unsafe_getindex, fld, cld, div,
-                @_inline_meta, eachindex, firstindex, lastindex, keys, isreal, OneTo,
-                Array, Vector, Matrix, view, ones, @propagate_inbounds, print_array,
-                split, iszero, permutedims, rad2deg, deg2rad
+import StaticArrays, Calculus
 
-import Base.Broadcast: BroadcastStyle, Broadcasted, AbstractArrayStyle, broadcastable,
-                        DefaultArrayStyle, broadcasted
+import DomainSets: Domain, indomain, UnionDomain, ProductDomain, Point, ∂,
+              elements, DifferenceDomain, Interval, ChebyshevInterval, boundary,
+              rightendpoint, leftendpoint, dimension, WrappedDomain, VcatDomain,
+              component, components, ncomponents
+
+using AbstractFFTs: Plan
+
+import Base: values, convert, getindex, setindex!, *, +, -, ==, <, <=, >, |, !,
+              !=, eltype, iterate, /, ^, \,
+              transpose, size, tail, broadcast, broadcast!, copyto!, copy,
+              to_index, (:), similar, map, vcat, hcat, hvcat, show, summary,
+              stride, sum, cumsum, imag, conj, inv, complex, reverse, exp,
+              sqrt, abs, abs2, sign, issubset, in, first, last, rand, intersect,
+              setdiff, isless, union, angle, join, isnan, isapprox, isempty,
+              sort, merge, minimum, maximum, extrema, argmax,
+              argmin, findmax, findmin, isfinite,
+              zeros, zero, one, promote_rule, repeat, length, resize!, isinf,
+              getproperty, findfirst, unsafe_getindex, fld, div,
+              eachindex, firstindex, lastindex, isreal,
+              OneTo, Array, Vector, Matrix, view, ones, @propagate_inbounds,
+              print_array, split, iszero, permutedims, rad2deg, deg2rad
+
+import Base.Broadcast: BroadcastStyle, Broadcasted, AbstractArrayStyle,
+              broadcastable, DefaultArrayStyle, broadcasted
 
 import Statistics: mean
 
-import LinearAlgebra: BlasInt, BlasFloat, norm, ldiv!, mul!, det, eigvals, cross,
-                        qr, qr!, rank, isdiag, istril, istriu, issymmetric, ishermitian,
-                        Tridiagonal, diagm, diagm_container, factorize, nullspace,
-                        Hermitian, Symmetric, adjoint, transpose, char_uplo
+import LinearAlgebra: BlasInt, BlasFloat, norm, ldiv!, mul!, det, cross,
+              qr, qr!, rank, isdiag, istril, istriu, issymmetric,
+              Tridiagonal, diagm, diagm_container, factorize,
+              nullspace, Hermitian, Symmetric, adjoint, transpose, char_uplo
 
 import SparseArrays: blockdiag
 
@@ -45,35 +47,38 @@ import SparseArrays: blockdiag
 # we need to import all special functions to use Calculus.symbolic_derivatives_1arg
 # we can't do importall Base as we replace some Base definitions
 import SpecialFunctions: sinpi, cospi, airy, besselh,
-                    asinh, acosh,atanh, erfcx, dawson, erf, erfi,
-                    sin, cos, sinh, cosh, airyai, airybi, airyaiprime, airybiprime,
-                    hankelh1, hankelh2, besselj, besselj0, bessely, besseli, besselk,
-                    besselkx, hankelh1x, hankelh2x, exp2, exp10, log2, log10,
-                    tan, tanh, csc, asin, acsc, sec, acos, asec,
-                    cot, atan, acot, sinh, csch, asinh, acsch,
-                    sech, acosh, asech, tanh, coth, atanh, acoth,
-                    expm1, log1p, lfact, sinc, cosc, erfinv, erfcinv, beta, lbeta,
-                    eta, zeta, gamma,  lgamma, polygamma, invdigamma, digamma, trigamma,
-                    abs, sign, log, expm1, tan, abs2, sqrt, angle, max, min, cbrt, log,
-                    atan, acos, asin, erfc, inv
-
-import StaticArrays: SVector
-
-import BandedMatrices: bandrange, bandshift,
-                        inbands_getindex, inbands_setindex!, bandwidth, AbstractBandedMatrix,
-                        colstart, colstop, colrange, rowstart, rowstop, rowrange,
-                        bandwidths, _BandedMatrix, BandedMatrix
+              sin, cos, cosh, exp2, exp10, log2, log10, csc, acsc, sec,
+              asec, cot, acot, sinh, csch, asinh, acsch,
+              sech, acosh, asech, tanh, coth, atanh, acoth,
+              log1p, lfact, sinc, cosc, beta, lbeta,
+              eta, zeta, polygamma, logabsgamma, loggamma,
+              abs, sign, log, expm1, tan, abs2, sqrt, angle,
+              max, min, cbrt, atan, acos, asin, inv,
+              besselj, bessely, besseli, besselk, besselkx,
+              hankelh1, hankelh2, hankelh1x, hankelh2x,
+              # functions from Calculus.symbolic_derivatives_1arg
+              erf, erfinv, erfc, erfcinv, erfi, gamma, lgamma,
+              digamma, invdigamma, trigamma,
+              airyai, airybi, airyaiprime, airybiprime,
+              besselj0, besselj1, bessely0, bessely1,
+              erfcx, dawson
+
+import BandedMatrices: bandrange, inbands_setindex!, bandwidth,
+              colstart, colstop, colrange, rowstart, rowstop, rowrange,
+              bandwidths, _BandedMatrix, BandedMatrix
 
 import BlockArrays: blocksize, block, blockaxes, blockindex
-import BlockBandedMatrices: blockbandwidth, blockbandwidths, blockcolstop, blockcolrange,
-                            blockcolstart, blockrowstop, blockrowstart, blockrowrange,
-                            subblockbandwidth, subblockbandwidths, _BlockBandedMatrix,
-                            _BandedBlockBandedMatrix, BandedBlockBandedMatrix, BlockBandedMatrix,
-                            isblockbanded, isbandedblockbanded, bb_numentries, BlockBandedSizes
+import BlockBandedMatrices: blockbandwidth, blockbandwidths, blockcolstop,
+              blockcolrange, blockcolstart, blockrowstop, blockrowstart,
+              subblockbandwidth, subblockbandwidths, _BlockBandedMatrix,
+              _BandedBlockBandedMatrix, BandedBlockBandedMatrix,
+              BlockBandedMatrix, isblockbanded, isbandedblockbanded,
+              bb_numentries, BlockBandedSizes
 
 import FillArrays: AbstractFill, getindex_value
 import LazyArrays: cache, CachedVector, cacheddata
-import InfiniteArrays: PosInfinity, InfRanges, AbstractInfUnitRange, OneToInf, InfiniteCardinal
+import InfiniteArrays: PosInfinity, InfRanges, AbstractInfUnitRange,
+              OneToInf, InfiniteCardinal
 
 
 # convenience for 1-d block ranges
@@ -93,7 +98,8 @@ export pad!, pad, chop!, sample,
        complexroots, roots,
        reverseorientation
 
-export .., Interval, ChebyshevInterval, leftendpoint, rightendpoint, endpoints, cache
+export .., Interval, ChebyshevInterval, leftendpoint, rightendpoint
+export endpoints, cache
 
 import Base: oneto
 

src/Fun.jl

@@ -245,7 +245,9 @@ iterate(x::ScalarFun) = (x, nothing)
 iterate(x::ScalarFun, ::Any) = nothing
 isempty(x::ScalarFun) = false
 
-iterate(A::ArrayFun, i=1) = (@_inline_meta; (i % UInt) - 1 < length(A) ? (@inbounds A[i], i + 1) : nothing)
+@inline function iterate(A::ArrayFun, i=1)
+    (i % UInt) - 1 < length(A) ? (@inbounds A[i], i + 1) : nothing
+end
 
 in(x::ScalarFun, y::ScalarFun) = x == y
 

src/specialfunctions.jl

@@ -467,6 +467,31 @@ for (funsym, exp) in Calculus.symbolic_derivatives_1arg()
     end
 end
 
+# Other special functions
+for f in [:logabsgamma]
+    @eval function $f(z::Fun{<:ConstantSpace, <:Real})
+        t = $f(Number(z))
+        Fun(t[1], space(z)), t[2]
+    end
+end
+function loggamma(z::Fun{<:ConstantSpace})
+    t = loggamma(Number(z))
+    Fun(t, space(z))
+end
+for f in [:gamma, :loggamma]
+    @eval begin
+        function $f(a, z::Fun{<:ConstantSpace})
+            t = $f(a, Number(z))
+            Fun(t, space(z))
+        end
+    end
+end
+
+for f in [:besselj, :besselk, :besselkx, :bessely, :besseli,
+            :hankelh1x, :hankelh2x, :hankelh1, :hankelh2]
+    @eval $f(nu, x::Fun{<:ConstantSpace}) = Fun($f(nu, Number(x)), space(x))
+end
+
 # Roots
 
 for op in (:(argmax),:(argmin))

test/runtests.jl

@@ -1,6 +1,7 @@
 using ApproxFunBase, LinearAlgebra, Random, Test
 using ApproxFunBase: ∞
 using Aqua
+using SpecialFunctions
 
 @testset "Project quality" begin
     Aqua.test_all(ApproxFunBase, ambiguities=false)
@@ -389,5 +390,29 @@ end
     end
 end
 
+@testset "Special functions" begin
+    pt = 0.5
+    x = Fun(pt, ConstantSpace(1..2))
+    for f in [erf, erfinv, erfc, erfcinv, erfi, gamma,
+                digamma, invdigamma, trigamma, loggamma,
+                airyai, airybi, airyaiprime, airybiprime,
+                besselj0, besselj1, bessely0, bessely1,
+                erfcx, dawson]
+
+        @test Number(f(x)) ≈ f(pt)
+        if f ∉ [erfinv, erfcinv, invdigamma]
+            @test Number(f(x*im)) ≈ f(pt*im)
+        end
+    end
+    for f in [besselj, bessely, besseli, besselk, besselkx,
+              hankelh1, hankelh2, hankelh1x, hankelh2x]
+        @test Number(f(1, x)) ≈ f(1, pt)
+    end
+    @test Number(logabsgamma(x)[1]) ≈ logabsgamma(pt)[1]
+    @test logabsgamma(x)[2] ≈ logabsgamma(pt)[2]
+    @test Number(gamma(2, x)) ≈ gamma(2, pt)
+    @test Number(gamma(2, x*im)) ≈ gamma(2, pt*im)
+end
+
 @time include("ETDRK4Test.jl")
 include("show.jl")