Improve specialization and other small fixes (#272)

* Make explicit NumberNotSeries and NumberNotSeriesN * Avoid Any as parametric types * Fix numbr2str for Complex{Interval} with tests * Tiny fix in evaluate * Improve specialization for Type, Function, and Vararg, following the manual
JuliaDiff · Apr 21, 2021 · dfb3f15 · dfb3f15
1 parent 176e16d
commit dfb3f15
Show file tree

Hide file tree

Showing 16 changed files with 191 additions and 195 deletions.
diff --git a/src/arithmetic.jl b/src/arithmetic.jl
@@ -12,7 +12,7 @@
 for T in (:Taylor1, :TaylorN)
 
     @eval begin
-        ==(a::$T{T}, b::$T{S}) where {T<:Number, S<:Number} = ==(promote(a,b)...)
+        ==(a::$T{T}, b::$T{S}) where {T<:Number,S<:Number} = ==(promote(a,b)...)
 
         function ==(a::$T{T}, b::$T{T}) where {T<:Number}
             if a.order != b.order
@@ -84,7 +84,7 @@ for (f, fc) in ((:+, :(add!)), (:-, :(subst!)))
 
     for T in (:Taylor1, :TaylorN)
         @eval begin
-            ($f)(a::$T{T}, b::$T{S}) where {T<:Number, S<:Number} =
+            ($f)(a::$T{T}, b::$T{S}) where {T<:Number,S<:Number} =
                 $f(promote(a,b)...)
 
             function $f(a::$T{T}, b::$T{T}) where {T<:Number}
@@ -102,7 +102,7 @@ for (f, fc) in ((:+, :(add!)), (:-, :(subst!)))
                 return $T(v, a.order)
             end
 
-            ($f)(a::$T{T}, b::S) where {T<:Number, S<:Number} =
+            ($f)(a::$T{T}, b::S) where {T<:Number,S<:Number} =
                 $f(promote(a,b)...)
 
             function $f(a::$T{T}, b::T) where {T<:Number}
@@ -111,7 +111,7 @@ for (f, fc) in ((:+, :(add!)), (:-, :(subst!)))
                 return $T(coeffs, a.order)
             end
 
-            ($f)(b::S, a::$T{T}) where {T<:Number, S<:Number} =
+            ($f)(b::S, a::$T{T}) where {T<:Number,S<:Number} =
                 $f(promote(b,a)...)
 
             function $f(b::T, a::$T{T}) where {T<:Number}
@@ -122,23 +122,23 @@ for (f, fc) in ((:+, :(add!)), (:-, :(subst!)))
             end
 
             ## add! and subst! ##
-            function ($fc)(v::$T{T}, a::$T{T}, k::Int) where {T}
+            function ($fc)(v::$T{T}, a::$T{T}, k::Int) where {T<:Number}
                 @inbounds v[k] = ($f)(a[k])
                 return nothing
             end
-            function ($fc)(v::$T, a::NumberNotSeries, k::Int)
-                @inbounds v[k] = k==0 ? ($f)(zero(v[0]),a) : zero(v[k])
+            function ($fc)(v::$T{T}, a::T, k::Int) where {T<:Number}
+                @inbounds v[k] = k==0 ? ($f)(a) : zero(a)
                 return nothing
             end
             function ($fc)(v::$T, a::$T, b::$T, k::Int)
                 @inbounds v[k] = ($f)(a[k], b[k])
                 return nothing
             end
-            function ($fc)(v::$T, a::$T, b::NumberNotSeries, k::Int)
+            function ($fc)(v::$T, a::$T, b::Number, k::Int)
                 @inbounds v[k] = k==0 ? ($f)(a[0], b) : ($f)(a[k], zero(b))
                 return nothing
             end
-            function ($fc)(v::$T, a::NumberNotSeries, b::$T, k::Int)
+            function ($fc)(v::$T, a::Number, b::$T, k::Int)
                 @inbounds v[k] = k==0 ? ($f)(a, b[0]) : ($f)(zero(a), b[k])
                 return nothing
             end
@@ -150,7 +150,7 @@ for (f, fc) in ((:+, :(add!)), (:-, :(subst!)))
             {T<:NumberNotSeriesN,S<:NumberNotSeriesN} = $f(promote(a,b)...)
 
         function $f(a::HomogeneousPolynomial{T}, b::HomogeneousPolynomial{T}) where
-                T<:NumberNotSeriesN
+                {T<:NumberNotSeriesN}
             @assert a.order == b.order
             v = similar(a.coeffs)
             @__dot__ v = $f(a.coeffs, b.coeffs)
@@ -164,7 +164,7 @@ for (f, fc) in ((:+, :(add!)), (:-, :(subst!)))
         end
 
         function ($f)(a::TaylorN{Taylor1{T}}, b::Taylor1{S}) where
-                {T<:NumberNotSeries, S<:NumberNotSeries}
+                {T<:NumberNotSeries,S<:NumberNotSeries}
             @inbounds aux = $f(a[0][1], b)
             R = eltype(aux)
             coeffs = Array{HomogeneousPolynomial{Taylor1{R}}}(undef, a.order+1)
@@ -215,26 +215,30 @@ end
 for T in (:Taylor1, :HomogeneousPolynomial, :TaylorN)
 
     @eval begin
-        *(a::Bool, b::$T) = *(convert(Int, a), b)
-
-        *(a::$T, b::Bool) = b * a
-
-        function *(a::T, b::$T) where {T<:NumberNotSeries}
+        function *(a::T, b::$T{S}) where {T<:NumberNotSeries,S<:NumberNotSeries}
             @inbounds aux = a * b.coeffs[1]
             v = Array{typeof(aux)}(undef, length(b.coeffs))
             @__dot__ v = a * b.coeffs
             return $T(v, b.order)
         end
 
-        *(b::$T, a::T) where {T<:NumberNotSeries} = a * b
+        *(b::$T{S}, a::T) where {T<:NumberNotSeries,S<:NumberNotSeries} = a * b
+
+        function *(a::T, b::$T{T}) where {T<:Number}
+            v = Array{T}(undef, length(b.coeffs))
+            @__dot__ v = a * b.coeffs
+            return $T(v, b.order)
+        end
+
+        *(b::$T{T}, a::T) where {T<:Number} = a * b
     end
 end
 
 for T in (:HomogeneousPolynomial, :TaylorN)
 
     @eval begin
         function *(a::Taylor1{T}, b::$T{Taylor1{S}}) where
-                {T<:NumberNotSeries, S<:NumberNotSeries}
+                {T<:NumberNotSeries,S<:NumberNotSeries}
 
             @inbounds aux = a * b.coeffs[1]
             R = typeof(aux)
@@ -246,22 +250,19 @@ for T in (:HomogeneousPolynomial, :TaylorN)
         *(b::$T{Taylor1{R}}, a::Taylor1{T}) where
             {T<:NumberNotSeries,R<:NumberNotSeries} = a * b
 
-        function *(a::$T{T}, b::Taylor1{$T{S}}) where {T<:NumberNotSeries, S<:NumberNotSeries}
+        function *(a::$T{T}, b::Taylor1{$T{S}}) where {T<:NumberNotSeries,S<:NumberNotSeries}
             @inbounds aux = a * b[0]
             R = typeof(aux)
             coeffs = Array{R}(undef, length(b.coeffs))
             @__dot__ coeffs = a * b.coeffs
             return Taylor1(coeffs, b.order)
         end
 
-        *(b::Taylor1{$T{S}}, a::$T{T}) where
-            {T<:NumberNotSeries,S<:NumberNotSeries} = a * b
+        *(b::Taylor1{$T{S}}, a::$T{T}) where {T<:NumberNotSeries,S<:NumberNotSeries} = a * b
     end
 end
 
 for (T, W) in ((:Taylor1, :Number), (:TaylorN, :NumberNotSeriesN))
-    @eval *(a::$T{T}, b::$T{S}) where {T<:$W, S<:$W} = *(promote(a,b)...)
-
     @eval function *(a::$T{T}, b::$T{T}) where {T<:$W}
         if a.order != b.order
             a, b = fixorder(a, b)
@@ -293,7 +294,7 @@ end
 
 # Internal multiplication functions
 for T in (:Taylor1, :TaylorN)
-    @eval @inline function mul!(c::$T{T}, a::$T{T}, b::$T{T}, k::Int) where {T}
+    @eval @inline function mul!(c::$T{T}, a::$T{T}, b::$T{T}, k::Int) where {T<:Number}
 
         if $T == Taylor1
             @inbounds c[k] = a[0] * b[k]
@@ -376,21 +377,27 @@ end
 
 
 ## Division ##
-function /(a::Taylor1{Rational{T}}, b::S) where {T<:Integer, S<:NumberNotSeries}
+function /(a::Taylor1{Rational{T}}, b::S) where {T<:Integer,S<:NumberNotSeries}
     R = typeof( a[0] // b)
     v = Array{R}(undef, a.order+1)
     @__dot__ v = a.coeffs // b
     return Taylor1(v, a.order)
 end
 
 for T in (:Taylor1, :HomogeneousPolynomial, :TaylorN)
-
-    @eval function /(a::$T{T}, b::S) where {T<:NumberNotSeries, S<:NumberNotSeries}
-        R = promote_type(T,S)
-        return convert($T{R}, a) * inv(convert(R, b))
+    @eval function /(a::$T{T}, b::S) where {T<:NumberNotSeries,S<:NumberNotSeries}
+        @inbounds aux = a.coeffs[1] / b
+        v = Array{typeof(aux)}(undef, length(a.coeffs))
+        @__dot__ v = a.coeffs / b
+        return $T(v, a.order)
     end
 
-    @eval /(a::$T, b::T) where {T<:NumberNotSeries} = a * inv(b)
+    @eval function /(a::$T{T}, b::T) where {T<:Number}
+        @inbounds aux = a.coeffs[1] / b
+        v = Array{typeof(aux)}(undef, length(a.coeffs))
+        @__dot__ v = a.coeffs / b
+        return $T(v, a.order)
+    end
 end
 
 for T in (:HomogeneousPolynomial, :TaylorN)
@@ -413,9 +420,7 @@ for T in (:HomogeneousPolynomial, :TaylorN)
     end
 end
 
-
-
-/(a::Taylor1{T}, b::Taylor1{S}) where {T<:Number, S<:Number} = /(promote(a,b)...)
+/(a::Taylor1{T}, b::Taylor1{S}) where {T<:Number,S<:Number} = /(promote(a,b)...)
 
 function /(a::Taylor1{T}, b::Taylor1{T}) where {T<:Number}
     if a.order != b.order
@@ -434,7 +439,7 @@ function /(a::Taylor1{T}, b::Taylor1{T}) where {T<:Number}
 end
 
 /(a::TaylorN{T}, b::TaylorN{S}) where
-    {T<:NumberNotSeriesN, S<:NumberNotSeriesN} = /(promote(a,b)...)
+    {T<:NumberNotSeriesN,S<:NumberNotSeriesN} = /(promote(a,b)...)
 
 function /(a::TaylorN{T}, b::TaylorN{T}) where {T<:NumberNotSeriesN}
     @assert !iszero(constant_term(b))
@@ -601,7 +606,7 @@ end
 # Specialize a method of `lu` for Matrix{Taylor1{T}}, which avoids pivoting,
 # since the polynomial field is not an ordered one.
 # We can't assume an ordered field so we first try without pivoting
-function lu(A::AbstractMatrix{Taylor1{T}}; check::Bool = true) where T
+function lu(A::AbstractMatrix{Taylor1{T}}; check::Bool = true) where {T<:Number}
     S = Taylor1{lutype(T)}
     F = lu!(copy_oftype(A, S), Val(false); check = false)
     if issuccess(F)

diff --git a/src/auxiliary.jl b/src/auxiliary.jl
@@ -161,8 +161,7 @@ getindex(a::TaylorN, c::Colon) = view(a.coeffs, c)
 getindex(a::TaylorN, u::StepRange{Int,Int}) where {T<:Number} =
     view(a.coeffs, u[:] .+ 1)
 
-function setindex!(a::TaylorN{T}, x::HomogeneousPolynomial{T}, n::Int) where
-        {T<:Number}
+function setindex!(a::TaylorN{T}, x::HomogeneousPolynomial{T}, n::Int) where {T<:Number}
     @assert x.order == n
     a.coeffs[n+1] = x
 end
@@ -275,7 +274,7 @@ end
 ## copyto! ##
 # Inspired from base/abstractarray.jl, line 665
 for T in (:Taylor1, :HomogeneousPolynomial, :TaylorN)
-    @eval function copyto!(dst::$T{T}, src::$T{T}) where {T}
+    @eval function copyto!(dst::$T{T}, src::$T{T}) where {T<:Number}
         length(dst) < length(src) && throw(ArgumentError(string("Destination has fewer elements than required; no copy performed")))
         destiter = eachindex(dst)
         y = iterate(destiter)

diff --git a/src/broadcasting.jl b/src/broadcasting.jl
@@ -12,13 +12,13 @@ import .Broadcast: BroadcastStyle, Broadcasted, broadcasted
 
 # BroadcastStyle definitions and basic precedence rules
 struct Taylor1Style{T} <: Base.Broadcast.AbstractArrayStyle{0} end
-Taylor1Style{T}(::Val{N}) where {T, N}= Base.Broadcast.DefaultArrayStyle{N}()
+Taylor1Style{T}(::Val{N}) where {T,N}= Base.Broadcast.DefaultArrayStyle{N}()
 BroadcastStyle(::Type{<:Taylor1{T}}) where {T} = Taylor1Style{T}()
 BroadcastStyle(::Taylor1Style{T}, ::Base.Broadcast.DefaultArrayStyle{0}) where {T} = Taylor1Style{T}()
 BroadcastStyle(::Taylor1Style{T}, ::Base.Broadcast.DefaultArrayStyle{1}) where {T} = Base.Broadcast.DefaultArrayStyle{1}()
 #
 struct HomogeneousPolynomialStyle{T} <: Base.Broadcast.AbstractArrayStyle{0} end
-HomogeneousPolynomialStyle{T}(::Val{N}) where {T, N}= Base.Broadcast.DefaultArrayStyle{N}()
+HomogeneousPolynomialStyle{T}(::Val{N}) where {T,N}= Base.Broadcast.DefaultArrayStyle{N}()
 BroadcastStyle(::Type{<:HomogeneousPolynomial{T}}) where {T} = HomogeneousPolynomialStyle{T}()
 BroadcastStyle(::HomogeneousPolynomialStyle{T}, ::Base.Broadcast.DefaultArrayStyle{0}) where {T} = HomogeneousPolynomialStyle{T}()
 BroadcastStyle(::HomogeneousPolynomialStyle{T}, ::Base.Broadcast.DefaultArrayStyle{1}) where {T} = Base.Broadcast.DefaultArrayStyle{1}()
@@ -91,7 +91,7 @@ Base.Broadcast.eltypes(t::Tuple{AbstractArray,TaylorN}) =
 # Adapted from Base.Broadcast.copyto!, base/broadcasting.jl, line 832
 for T in (:Taylor1, :HomogeneousPolynomial, :TaylorN)
     @eval begin
-        @inline function copyto!(dest::$T{T}, bc::Broadcasted) where {T}
+        @inline function copyto!(dest::$T{T}, bc::Broadcasted) where {T<:Number}
             axes(dest) == axes(bc) || Base.Broadcast.throwdm(axes(dest), axes(bc))
             # Performance optimization: broadcast!(identity, dest, A) is equivalent to copyto!(dest, A) if indices match
             if bc.f === identity && bc.args isa Tuple{$T{T}} # only a single input argument to broadcast!
@@ -109,9 +109,9 @@ end
 
 
 # Broadcasted extensions
-@inline broadcasted(::Taylor1Style{T}, ::Type{Float32}, a::Taylor1{T}) where {T} =
+@inline broadcasted(::Taylor1Style{T}, ::Type{Float32}, a::Taylor1{T}) where {T<:Number} =
     Taylor1(Float32.(a.coeffs), a.order)
-@inline broadcasted(::TaylorNStyle{T}, ::Type{Float32}, a::TaylorN{T}) where {T} =
+@inline broadcasted(::TaylorNStyle{T}, ::Type{Float32}, a::TaylorN{T}) where {T<:Number} =
     convert(TaylorN{Float32}, a)
 
 # # This prevents broadcasting being applied to the Taylor1/TaylorN params

diff --git a/src/calculus.jl b/src/calculus.jl
@@ -69,7 +69,7 @@ derivative!(p::Taylor1, a::Taylor1, k::Int) = k < a.order ? p[k] = (k+1)*a[k+1]
 Compute recursively the `Taylor1` polynomial of the n-th derivative of
 `a::Taylor1`.
 """
-function derivative(a::Taylor1{T}, n::Int) where {T <: Number}
+function derivative(a::Taylor1{T}, n::Int) where {T<:Number}
     @assert a.order ≥ n ≥ 0
     if n==0
         return a
@@ -99,7 +99,7 @@ end
 Return the integral of `a::Taylor1`. The constant of integration
 (0-th order coefficient) is set to `x`, which is zero if ommitted.
 """
-function integrate(a::Taylor1{T}, x::S) where {T<:Number, S<:Number}
+function integrate(a::Taylor1{T}, x::S) where {T<:Number,S<:Number}
     order = get_order(a)
     aa = a[0]/1 + zero(x)
     R = typeof(aa)
@@ -250,9 +250,7 @@ function jacobian(vf::Array{TaylorN{T},1}) where {T<:Number}
 
     return transpose(jac)
 end
-function jacobian(vf::Array{TaylorN{T},1}, vals::Array{S,1}) where
-        {T<:Number, S<:Number}
-
+function jacobian(vf::Array{TaylorN{T},1}, vals::Array{S,1}) where {T<:Number,S<:Number}
     R = promote_type(T,S)
     numVars = get_numvars()
     @assert length(vf) == numVars == length(vals)
@@ -319,7 +317,7 @@ Return the hessian matrix (jacobian of the gradient) of `f::TaylorN`,
 evaluated at the vector `vals`. If `vals` is ommited, it is evaluated at
 zero.
 """
-hessian(f::TaylorN{T}, vals::Array{S,1}) where {T<:Number, S<:Number} =
+hessian(f::TaylorN{T}, vals::Array{S,1}) where {T<:Number,S<:Number} =
     (R = promote_type(T,S); jacobian( gradient(f), vals::Array{R,1}) )
 
 hessian(f::TaylorN{T}) where {T<:Number} = hessian( f, zeros(T, get_numvars()) )
@@ -406,9 +404,9 @@ function integrate(a::TaylorN, r::Int, x0::TaylorN)
     res = integrate(a, r)
     return x0+res
 end
-integrate(a::TaylorN, r::Int, x0::NumberNotSeries) =
+integrate(a::TaylorN, r::Int, x0) =
     integrate(a,r,TaylorN(HomogeneousPolynomial([convert(eltype(a),x0)], 0)))
 
 integrate(a::TaylorN, s::Symbol) = integrate(a, lookupvar(s))
 integrate(a::TaylorN, s::Symbol, x0::TaylorN) = integrate(a, lookupvar(s), x0)
-integrate(a::TaylorN, s::Symbol, x0::NumberNotSeries) = integrate(a, lookupvar(s), x0)
+integrate(a::TaylorN, s::Symbol, x0) = integrate(a, lookupvar(s), x0)
diff --git a/src/constructors.jl b/src/constructors.jl
@@ -163,7 +163,7 @@ struct TaylorN{T<:Number} <: AbstractSeries{T}
     end
 end
 
-TaylorN(x::TaylorN{T}) where T<:Number = x
+TaylorN(x::TaylorN{T}) where {T<:Number} = x
 function TaylorN(x::Array{HomogeneousPolynomial{T},1}, order::Int) where {T<:Number}
     if order == 0
         order = maxorderH(x)
@@ -202,11 +202,10 @@ TaylorN(nv::Int; order::Int=get_order()) = TaylorN(Float64, nv, order=order)
 
 
 # A `Number` which is not an `AbstractSeries`
-const NumberNotSeries = Union{setdiff(subtypes(Number), [AbstractSeries])...}
+const NumberNotSeries = Union{Real,Complex}
 
 # A `Number` which is not `TaylorN` nor a `HomogeneousPolynomial`
-const NumberNotSeriesN =
-    Union{setdiff(subtypes(Number), [AbstractSeries])..., Taylor1}
+const NumberNotSeriesN = Union{Real,Complex,Taylor1}
 
 ## Additional Taylor1 outer constructor ##
-Taylor1{T}(x::S) where {T<:Number, S<:NumberNotSeries} = Taylor1([convert(T,x)], 0)
+Taylor1{T}(x::S) where {T<:Number,S<:NumberNotSeries} = Taylor1([convert(T,x)], 0)