Merge remote-tracking branch 'origin/master' into jp/evaluate

Author: PerezHz

Project.toml

@@ -1,6 +1,6 @@
 name = "TaylorSeries"
 uuid = "6aa5eb33-94cf-58f4-a9d0-e4b2c4fc25ea"
-version = "0.22.1"
+version = "0.22.4"
 repo = "https://github.com/JuliaDiff/TaylorSeries.jl.git"
 
 [deps]

ext/TaylorSeriesIAExt.jl

@@ -318,15 +318,16 @@ function TS.sqrt!(c::Taylor1{Interval{T}}, a::Taylor1{Interval{T}},
     @inbounds for i = imin+1:imax
         c[k] += c[i] * c[k+k0-i]
     end
+    intvl2 = interval(T(2))
     if k+k0 ≤ a_order
-        @inbounds aux = a[k+k0] - interval(T(2)) * c[k]
+        @inbounds aux = a[k+k0] - intvl2 * c[k]
     else
-        @inbounds aux = - interval(T(2)) * c[k]
+        @inbounds aux = - intvl2 * c[k]
     end
     if kodd == 0
         @inbounds aux = aux - c[kend+k0+1]^2
     end
-    @inbounds c[k] = aux / (interval(T(2)) * c[k0])
+    @inbounds c[k] = aux / (intvl2 * c[k0])
     return nothing
 end
 
@@ -347,12 +348,13 @@ function TS.sqrt!(c::TaylorN{Interval{T}}, a::TaylorN{Interval{T}},
         # @inbounds c[k] <- c[k] - (c[kend+1])^2
         @inbounds TS.mul_scalar!(c[k], -interval(T(1)), c[kend+1], c[kend+1])
     end
+    intvl2 = interval(T(2))
     @inbounds for i = 1:kend
         # c[k] <- c[k] - 2*c[i]*c[k-i]
-        TS.mul_scalar!(c[k], -interval(T(2)), c[i], c[k-i])
+        TS.mul_scalar!(c[k], -intvl2, c[i], c[k-i])
     end
     # @inbounds c[k] <- c[k] / (2*c[0])
-    TS.div!(c[k], c[k], interval(T(2))*constant_term(c))
+    TS.div!(c[k], c[k], intvl2*constant_term(c))
 
     return nothing
 end
@@ -469,11 +471,12 @@ function TS.exp!(c::Taylor1{Interval{T}}, a::Taylor1{Interval{T}}, k::Int) where
         @inbounds c[0] = exp(constant_term(a))
         return nothing
     end
-    @inbounds c[k] = interval(T(k)) * a[k] * c[0]
+    intvlk = interval(T(k))
+    @inbounds c[k] = intvlk * a[k] * c[0]
     @inbounds for i = 1:k-1
         c[k] += interval(T(k-i)) * a[k-i] * c[i]
     end
-    @inbounds c[k] = c[k] / interval(T(k))
+    @inbounds c[k] = c[k] / intvlk
     return nothing
 end
 
@@ -482,11 +485,12 @@ function TS.exp!(c::TaylorN{Interval{T}}, a::TaylorN{Interval{T}}, k::Int) where
         @inbounds c[0] = exp(constant_term(a))
         return nothing
     end
-    @inbounds TS.mul!(c[k], interval(T(k)) * a[k], c[0])
+    intvlk = interval(T(k))
+    @inbounds TS.mul!(c[k], intvlk * a[k], c[0])
     @inbounds for i = 1:k-1
         TS.mul!(c[k], interval(T(k-i)) * a[k-i], c[i])
     end
-    @inbounds c[k] = c[k] / interval(T(k))
+    @inbounds c[k] = c[k] / intvlk
     return nothing
 end
 
@@ -496,11 +500,12 @@ function TS.expm1!(c::Taylor1{Interval{T}}, a::Taylor1{Interval{T}}, k::Int) whe
         return nothing
     end
     c0 = c[0]+one(c[0])
-    @inbounds c[k] = interval(T(k)) * a[k] * c0
+    intvlk = interval(T(k))
+    @inbounds c[k] = intvlk * a[k] * c0
     @inbounds for i = 1:k-1
         c[k] += interval(T(k-i)) * a[k-i] * c[i]
     end
-    @inbounds c[k] = c[k] / interval(T(k))
+    @inbounds c[k] = c[k] / intvlk
     return nothing
 end
 
@@ -510,11 +515,12 @@ function TS.expm1!(c::TaylorN{Interval{T}}, a::TaylorN{Interval{T}}, k::Int) whe
         return nothing
     end
     c0 = c[0]+one(c[0])
-    @inbounds TS.mul!(c[k], interval(T(k)) * a[k], c0)
+    intvlk = interval(T(k))
+    @inbounds TS.mul!(c[k], intvlk * a[k], c0)
     @inbounds for i = 1:k-1
         TS.mul!(c[k], interval(T(k-i)) * a[k-i], c[i])
     end
-    @inbounds c[k] = c[k] / interval(T(k))
+    @inbounds c[k] = c[k] / intvlk
     return nothing
 end
 
@@ -601,8 +607,9 @@ function TS.sincos!(s::Taylor1{Interval{T}}, c::Taylor1{Interval{T}},
         s[k] += x * c[k-i]
         c[k] -= x * s[k-i]
     end
-    @inbounds s[k] = s[k] / interval(T(k))
-    @inbounds c[k] = c[k] / interval(T(k))
+    intvlk = interval(T(k))
+    @inbounds s[k] = s[k] / intvlk
+    @inbounds c[k] = c[k] / intvlk
     return nothing
 end
 
@@ -621,8 +628,9 @@ function TS.sincos!(s::TaylorN{Interval{T}}, c::TaylorN{Interval{T}},
         TS.mul!(s[k], x, c[k-i])
         TS.mul!(c[k], -x, s[k-i])
     end
-    @inbounds s[k] = s[k] / interval(T(k))
-    @inbounds c[k] = c[k] / interval(T(k))
+    intvlk = interval(T(k))
+    @inbounds s[k] = s[k] / intvlk
+    @inbounds c[k] = c[k] / intvlk
     return nothing
 end
 
@@ -634,11 +642,12 @@ function TS.tan!(c::Taylor1{Interval{T}}, a::Taylor1{Interval{T}},
         @inbounds c2[0] = aux^2
         return nothing
     end
-    @inbounds c[k] = interval(T(k)) * a[k] * c2[0]
+    intvlk = interval(T(k))
+    @inbounds c[k] = intvlk * a[k] * c2[0]
     @inbounds for i = 1:k-1
         c[k] += interval(T(k-i)) * a[k-i] * c2[i]
     end
-    @inbounds c[k] = a[k] + c[k]/interval(T(k))
+    @inbounds c[k] = a[k] + c[k]/intvlk
     TS.sqr!(c2, c, zero(c[0]), k)
     return nothing
 end
@@ -651,11 +660,12 @@ function TS.tan!(c::TaylorN{Interval{T}}, a::TaylorN{Interval{T}},
         @inbounds c2[0] = aux^2
         return nothing
     end
-    @inbounds TS.mul!(c[k], interval(T(k)) * a[k], c2[0])
+    intvlk = interval(T(k))
+    @inbounds TS.mul!(c[k], intvlk * a[k], c2[0])
     @inbounds for i = 1:k-1
         TS.mul!(c[k], interval(T(k-i)) * a[k-i], c2[i])
     end
-    @inbounds c[k] = a[k] + c[k]/interval(T(k))
+    @inbounds c[k] = a[k] + c[k]/intvlk
     TS.sqr!(c2, c, zero(c[0][1]), k)
     return nothing
 end
@@ -777,8 +787,9 @@ function TS.sinhcosh!(s::Taylor1{Interval{T}}, c::Taylor1{Interval{T}},
         s[k] += x * c[k-i]
         c[k] += x * s[k-i]
     end
-    s[k] = s[k] / interval(T(k))
-    c[k] = c[k] / interval(T(k))
+    intvlk = interval(T(k))
+    s[k] = s[k] / intvlk
+    c[k] = c[k] / intvlk
     return nothing
 end
 
@@ -797,8 +808,9 @@ function TS.sinhcosh!(s::TaylorN{Interval{T}}, c::TaylorN{Interval{T}},
         TS.mul!(s[k], x, c[k-i])
         TS.mul!(c[k], x, s[k-i])
     end
-    s[k] = s[k] / interval(T(k))
-    c[k] = c[k] / interval(T(k))
+    intvlk = interval(T(k))
+    s[k] = s[k] / intvlk
+    c[k] = c[k] / intvlk
     return nothing
 end
 
@@ -1045,7 +1057,7 @@ function TS._evaluate(a::HomogeneousPolynomial{T},
     @inbounds suma = a[1]*interval(zero(T))
     for (i, a_coeff) in enumerate(a.coeffs)
         TS._isthinzero(a_coeff) && continue
-        @inbounds tmp = prod(dx .^ ct[i])
+        @inbounds tmp = prod(Base.literal_pow.(^, dx, Val.(ct[i])))
         suma += a_coeff * tmp
     end
     return suma
@@ -1112,7 +1124,7 @@ function TS._evaluate(a::HomogeneousPolynomial{T},
     suma = zero(a[1])*vals[1]
     for (i, a_coeff) in enumerate(a.coeffs)
         TS._isthinzero(a_coeff) && continue
-        @inbounds tmp = prod( vals .^ ct[i] )
+        @inbounds tmp = prod( Base.literal_pow.(^, vals, Val.(ct[i])) )
         suma += a_coeff * tmp
     end
     return suma

src/TaylorSeries.jl

@@ -36,7 +36,7 @@ import Base: iterate, size, eachindex, firstindex, lastindex,
 import Base: zero, one, zeros, ones, isinf, isnan, iszero, isless,
     convert, promote_rule, promote, show, sum!,
     real, imag, conj, adjoint,
-    rem, mod, mod2pi, abs, abs2,
+    round, rem, mod, mod2pi, abs, abs2,
     sqrt, exp, expm1, log, log1p,
     sin, cos, sincos, sinpi, cospi, sincospi, tan,
     asin, acos, atan, sinh, cosh, tanh, atanh, asinh, acosh,

src/evaluate.jl

@@ -224,10 +224,10 @@ function _evaluate(a::HomogeneousPolynomial{T}, vals::NTuple) where {T}
     order(a) == 0 && return a[1]*one(vals[1])
     ct = a.space.coeff_table[order(a)+1]
     suma = zero(a[1])*vals[1]
-    vv = vals .^ ct[1]
+    vv = Base.literal_pow.(^, vals, Val.(ct[1]))
     for (i, a_coeff) in enumerate(a.coeffs)
         TS._isthinzero(a_coeff) && continue
-        @inbounds vv .= vals .^ ct[i]
+        @inbounds vv .= Base.literal_pow.(^, vals, Val.(ct[i]))
         tmp = prod( vv )
         suma += a_coeff * tmp
     end
@@ -441,7 +441,7 @@ function _evaluate!(suma::TaylorN{T}, a::HomogeneousPolynomial{T}, ind::Int,
         suma[0][1] = a[1]*one(val)
         return nothing
     end
-    vv = val .^ (0:order)
+    vv = Base.literal_pow.(^, val, Val.(0:order))
     vct = zero(a.space.coeff_table[order+1][1])
     zct = zero(a.space.coeff_table[order+1][1])
     for (i, a_coeff) in enumerate(a.coeffs)

src/other_functions.jl

@@ -23,6 +23,11 @@ for T in (:Taylor1, :HomogeneousPolynomial, :TaylorN)
     @eval isnan(a::$T) = any(isnan, a.coeffs)
 end
 
+# Rounding
+round(::Type{T}, x::AbstractSeries; kwargs...) where {T <: NumberNotSeries} =
+    round(T, x, RoundNearest; kwargs...)
+round(::Type{T}, x::AbstractSeries, r::RoundingMode; kwargs...) where {T <: NumberNotSeries} =
+    round(T, constant_term(x), r; kwargs...)
 
 ## Division functions: rem and mod ##
 for op in (:mod, :rem)
@@ -264,7 +269,7 @@ end
 
 #update! function for TaylorN
 function update!(a::TaylorN{T}, vals::Vector{T}) where {T<:Number}
-    a.coeffs .= evaluate(a, variables(order(a)) .+ vals).coeffs
+    a.coeffs .= evaluate(a, variables(a.space; order=order(a)) .+ vals).coeffs
     return nothing
 end
 function update!(a::TaylorN{T}, vals::Vector{S}) where {T<:Number, S<:Number}

test/manyvariables.jl

@@ -172,7 +172,7 @@ end
     @test get_variable_symbols() == [:x, :y]
     @test get_variable_symbols() == [:x, :y]
     @test TS.lookupvar(:x) == 1
-    @test TS.lookupvar(:α) == 0
+    @test TS.lookupvar(TS.default_space[], :α) == 0
     @test TS.get_variable_names() == ["x", "y"]
     @test x == HomogeneousPolynomial(Float64, 1)
     @test x == HomogeneousPolynomial(1)
@@ -364,6 +364,7 @@ end
     @test (1/(1-xT))[3] == HomogeneousPolynomial([1.0],3)
     @test xH^20 == HomogeneousPolynomial([0], order())
     @test (yT/(1-xT))[4] == xH^3 * yH
+    @test round(Int, 1.0 + xT) === round(Int, 1.0 + yT) === round(Int, 1.0)
     @test mod(1+xT,1) == +xT
     @test (rem(1+xT,1))[0] == 0
     @test mod(1+xT,1.0) == +xT

test/mixtures.jl

@@ -126,6 +126,7 @@ using Test
         @test -2*tN1 < -tN1^2 ≤ 0
     end
 
+    @test round(Int, tN1 + 1.0) === round(Int, t1N + 1.0) === round(Int, 1.0)
     @test mod(tN1+1,1.0) == 0+tN1
     @test mod(tN1-1.125,2) == 0.875+tN1
     @test (rem(tN1+1.125,1.0))[0][1] == 0.125 + t

test/onevariable.jl

@@ -237,6 +237,7 @@ Base.iszero(::SymbNumber) = false
     @test trational^3/complex(7,1) == Taylor1([0,0,0,complex(7//50,-1//50)],15)
     @test sqrt(zero(t)) == zero(t)
 
+    @test round(Int, 4.1 + t) === round(Int, 4.1)
     @test isapprox( rem(4.1 + t,4)[0], 0.1 )
     @test isapprox( mod(4.1 + t,4)[0], 0.1 )
     @test isapprox( rem(1+Taylor1(Int,4),4.0)[0], 1.0 )