Some temporary updates related to JuliaDiff/TaylorSeries.jl#361

Author: PerezHz

src/parse_eqs.jl

@@ -35,11 +35,12 @@ mutable struct BookKeeping
     v_array4 :: Vector{Symbol}
     v_preamb :: Vector{Union{Symbol,Expr}}
     retvar   :: Symbol
+    numaux   :: Int
 
     function BookKeeping()
         return new(Dict{Symbol, Expr}(), Dict{Union{Symbol,Expr}, Number}(),
             Dict{Symbol, Expr}(), Symbol[], Symbol[], Symbol[], Symbol[], Symbol[], Symbol[],
-            Union{Symbol,Expr}[], :nothing)
+            Union{Symbol,Expr}[], :nothing, 0)
     end
 end
 
@@ -343,7 +344,7 @@ function _preamble_body(fnbody, fnargs)
     bkkeep = BookKeeping()
 
     # Rename vars to have the body in non-indexed form; bkkeep has different entries
-    # for bookkeeping variables/symbolds, including indexed ones
+    # for bookkeeping variables/symbols, including indexed ones
     fnbody, bkkeep.d_indx = _rename_indexedvars(fnbody)
 
     # Create `newfnbody` which corresponds to `fnbody`, cleaned (without irrelevant comments)
@@ -372,7 +373,7 @@ function _preamble_body(fnbody, fnargs)
     newfnbody = subs(newfnbody, bkkeep.d_indx)
 
     # Define retvar; for scalar eqs is the last entry included in v_newvars
-    bkkeep.retvar = length(fnargs) == 3 ? subs(bkkeep.v_newvars[end], bkkeep.d_indx) : fnargs[1]
+    bkkeep.retvar = length(fnargs) == 3 ? subs(bkkeep.v_newvars[end-bkkeep.numaux], bkkeep.d_indx) : fnargs[1]
 
     return defspreamble, defsprealloc, newfnbody, bkkeep
 end
@@ -836,10 +837,11 @@ function _replacecalls!(bkkeep::BookKeeping, fnold::Expr, newvar::Symbol)
         end
 
     elseif ll == 3
-        # Binary call; no auxiliary expressions needed
+        # Binary call
         newarg2 = fnold.args[3]
 
         # Replacements
+        # @show TaylorSeries._dict_binary_calls[dcall]
         fnexpr, def_fnexpr, aux_fnexpr = TaylorSeries._dict_binary_calls[dcall]
         fnexpr = subs(fnexpr,
             Dict(:_res => newvar, :_arg1 => newarg1, :_arg2 => newarg2, :_k => :ord))
@@ -861,6 +863,28 @@ function _replacecalls!(bkkeep::BookKeeping, fnold::Expr, newvar::Symbol)
         #     Dict(:_res => newvar, :_arg1 => :(constant_term($(newarg1))),
         #         :_arg2 => :(constant_term($(newarg2))), :_k => :ord))
 
+        # Auxiliary expression
+        if aux_fnexpr.head != :nothing
+            newaux = genname()
+
+            aux_alloc = :( _res = Taylor1($(aux_fnexpr.args[2]), order) )
+            aux_alloc = subs(aux_alloc,
+                Dict(:_res => newaux, :_arg1 => :(constant_term($(newarg1))), :_aux => newaux))
+
+            aux_fnexpr = Expr(:block,
+                :(TaylorSeries.zero!(_res)),
+                :(_res.coeffs[1] = $(aux_fnexpr.args[2])) )
+            aux_fnexpr = subs(aux_fnexpr,
+                Dict(:_res => newaux, :_arg1 => :(constant_term($(newarg1))), :_aux => newaux))
+
+            fnexpr = subs(fnexpr, Dict(:_aux => newaux))
+            if newvar ∈ bkkeep.v_arraydecl
+                push!(bkkeep.v_arraydecl, newaux)
+            else
+                bkkeep.numaux += 1
+                push!(bkkeep.v_newvars, newaux)
+            end
+        end
     else
         # Recognized call, but not a unary or binary call; copy expression
         fnexpr = :($newvar = $fnold)

test/runtests.jl

@@ -1,6 +1,7 @@
 # This file is part of the TaylorIntegration.jl package; MIT licensed
 
 testfiles = (
+    "solution.jl",
     "one_ode.jl",
     "many_ode.jl",
     "complex.jl",

test/taylorize.jl

@@ -246,7 +246,7 @@ import Logging: Warn
     end
 
 
-    # Pendulum integrtf = 100.0
+    # # Pendulum integration
     @testset "Integration of the pendulum and DiffEqs interface" begin
         @taylorize function pendulum!(dx::Array{T,1}, x::Array{T,1}, p, t) where {T}
             dx[1] = x[2]
@@ -1303,7 +1303,7 @@ import Logging: Warn
             _tmax = tf
             _q0 = q0TN
             _params = nothing
-            _maxsteps = 3000
+            _maxsteps = 100
             _T = eltype(_t0)
             _U = eltype(_q0)
             # Initialize the vector of Taylor1 expansions
@@ -1317,12 +1317,13 @@ import Logging: Warn
             end
             # Determine if specialized jetcoeffs! method exists
             __parse_eqs, __rv = TaylorIntegration._determine_parsing!(true, kepler1!, _t, _x, _dx, _params);
+            @test __parse_eqs
             # Re-initialize the Taylor1 expansions
             _t = _t0 + Taylor1( _T, _order )
             _x .= Taylor1.( _q0, _order )
             _dx .= Taylor1.( zero.(_q0), _order)
-            solTN = @inferred TaylorIntegration._taylorinteg!(Val(true), kepler1!, _t, _x, _dx, _q0, _t0, _tmax, _abstol, __rv, _params; parse_eqs=__parse_eqs, maxsteps=_maxsteps)
-            solTN2 = @inferred TaylorIntegration._taylorinteg!(Val(false), kepler1!, _t, _x, _dx, _q0, _t0, _tmax, _abstol, __rv, _params; parse_eqs=__parse_eqs, maxsteps=_maxsteps)
+            solTN = @test_logs((Warn, max_iters_reached()), @inferred TaylorIntegration._taylorinteg!(Val(true), kepler1!, _t, _x, _dx, _q0, _t0, _tmax, _abstol, __rv, _params; parse_eqs=__parse_eqs, maxsteps=_maxsteps) )
+            solTN2 = @test_logs((Warn, max_iters_reached()), @inferred TaylorIntegration._taylorinteg!(Val(false), kepler1!, _t, _x, _dx, _q0, _t0, _tmax, _abstol, __rv, _params; parse_eqs=__parse_eqs, maxsteps=_maxsteps) )
             @test solTN isa TaylorSolution{typeof(_t0), eltype(_q0), ndims(solTN.x), typeof(solTN.t), typeof(solTN.x), typeof(solTN.p), Nothing, Nothing, Nothing}
             @test solTN2 isa TaylorSolution{typeof(_t0), eltype(_q0), ndims(solTN.x), typeof(solTN.t), typeof(solTN.x), Nothing, Nothing, Nothing, Nothing}
         end