Fixes and additions in evaluate, including some in-line methods (#357)

* Minor fixes in sqr_orderzero! and power_by_squaring

* Fix methods of evaluate for Taylor1{TaylorN{T}}
with tests. Some methods have been commented for now

* Run Aqua tests only once at the end

* Improvements when dealing with intervals

* Add evaluate(::Taylor1{T}, ::TaylorN{T}), and other fixes

* Add evaluate(::Taylor1{T}, ::Taylor1{TaylorN{T})

* Add _horner! as in-line evaluate method for certain mixtures

* Add `fixorder` when needed for evaluate

* Fix typo

* Another fix (in fixorder)

* Add methods for evaluate of a TaylorN by a value or TaylorN, including in-line methods

* Fixes and tests

* More tests a minor improvement

* Avoid using `@isonethread` in evaluate, differentiate and integrate

* Bump patch version

Author: lbenet

Project.toml

@@ -1,7 +1,7 @@
 name = "TaylorSeries"
 uuid = "6aa5eb33-94cf-58f4-a9d0-e4b2c4fc25ea"
 repo = "https://github.com/JuliaDiff/TaylorSeries.jl.git"
-version = "0.17.6"
+version = "0.17.7"
 
 [deps]
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"

ext/TaylorSeriesIAExt.jl

@@ -159,7 +159,7 @@ for T in (:Taylor1, :TaylorN)
 end
 
 
-function evaluate(a::Taylor1{T}, dx::Interval{S}) where {T<:Real, S<:Real}
+function evaluate(a::Taylor1, dx::Interval{S}) where {S<:Real}
     order = a.order
     uno = one(dx)
     dx2 = dx^2
@@ -183,34 +183,13 @@ function evaluate(a::TaylorN, dx::IntervalBox{N,T}) where {T<:Real,N}
     @assert N == get_numvars()
     a_length = length(a)
     suma = zero(constant_term(a)) + Interval{T}(0, 0)
-    @inbounds for homPol in length(a):-1:1
-        suma += evaluate(a.coeffs[homPol], dx)
+    @inbounds for homPol in reverse(eachindex(a))
+        suma += evaluate(a[homPol], dx)
     end
 
     return suma
 end
 
-function evaluate(a::Taylor1{TaylorN{T}}, dx::Interval{S}) where {T<:Real, S<:Real}
-    order = a.order
-    uno = one(dx)
-    dx2 = dx^2
-    if iseven(order)
-        kend = order-2
-        @inbounds sum_even = a[end]*uno
-        @inbounds sum_odd = a[end-1]*zero(dx)
-    else
-        kend = order-3
-        @inbounds sum_odd = a[end]*uno
-        @inbounds sum_even = a[end-1]*uno
-    end
-    @inbounds for k in kend:-2:0
-        sum_odd = sum_odd*dx2 + a[k+1]
-        sum_even = sum_even*dx2 + a[k]
-    end
-    return sum_even + sum_odd*dx
-end
-
-
 function evaluate(a::HomogeneousPolynomial, dx::IntervalBox{N,T}) where {T<:Real,N}
     @assert N == get_numvars()
     dx == IntervalBox(-1..1, Val(N)) && return _evaluate(a, dx, Val(true))
@@ -226,16 +205,16 @@ function _evaluate(a::HomogeneousPolynomial, dx::IntervalBox{N,T}, ::Val{true} )
     @inbounds suma = a[1]*Interval{T}(0,0)
 
     Ieven = Interval{T}(0,1)
-    for (i,a_coeff) in enumerate(a.coeffs)
+    for (i, a_coeff) in enumerate(a.coeffs)
         iszero(a_coeff) && continue
         if isodd(sum(ct[i]))
-            suma += sum(a_coeff) * dx[1]
-            continue
-        end
-        @inbounds tmp = iseven(ct[i][1]) ? Ieven : dx[1]
-        for n in 2:N
-            @inbounds vv = iseven(ct[i][n]) ? Ieven : dx[1]
-            tmp *= vv
+            tmp = dx[1]
+        else
+            tmp = Ieven
+            for n in eachindex(ct[i])
+                iseven(ct[i][n]) && continue
+                tmp *= dx[1]
+            end
         end
         suma += a_coeff * tmp
     end

src/auxiliary.jl

@@ -259,7 +259,7 @@ end
 ## fixorder ##
 for T in (:Taylor1, :TaylorN)
     @eval begin
-        @inline function fixorder(a::$T{T}, b::$T{T}) where {T<:Number}
+        @inline function fixorder(a::$T, b::$T)
             a.order == b.order && return a, b
             minorder = _minorder(a, b)
             return $T(copy(a.coeffs), minorder), $T(copy(b.coeffs), minorder)

src/calculus.jl

@@ -141,9 +141,10 @@ function differentiate(a::HomogeneousPolynomial, r::Int)
     coeffs = zeros(T, num_coeffs)
     @inbounds posTb = pos_table[a.order]
     @inbounds num_coeffs = size_table[a.order+1]
-
+    ct = deepcopy(coeff_table[a.order+1])
     @inbounds for i = 1:num_coeffs
-        iind = @isonethread coeff_table[a.order+1][i]
+        # iind = @isonethread coeff_table[a.order+1][i]
+        iind = ct[i]
         n = iind[r]
         n == 0 && continue
         iind[r] -= 1
@@ -370,9 +371,10 @@ function integrate(a::HomogeneousPolynomial, r::Int)
 
     T = promote_type(TS.numtype(a), TS.numtype(a[1]/1))
     coeffs = zeros(T, size_table[a.order+2])
-
+    ct = deepcopy(coeff_table[a.order+1])
     @inbounds for i = 1:num_coeffs
-        iind = @isonethread coeff_table[a.order+1][i]
+        # iind = @isonethread coeff_table[a.order+1][i]
+        iind = ct[i]
         n = iind[r]
         n == order_max && continue
         iind[r] += 1