Scattered fixes

PerezHz · Dec 20, 2024 · 08c9712 · 08c9712
1 parent 4be8ec2
commit 08c9712
Show file tree

Hide file tree

Showing 4 changed files with 16 additions and 56 deletions.
diff --git a/src/integrator/cache.jl b/src/integrator/cache.jl
@@ -22,16 +22,16 @@ struct ScalarCache{TV, XV, PSOL} <: AbstractTaylorIntegrationCache
     psol::PSOL
 end
 
-function init_cache(::Type{ScalarCache}, dense::Val{D}, t0::T, x::Taylor1{U}, maxsteps::Int) where {D, U, T}
-    return ScalarCache(
+function init_cache(cachetype::Type{ScalarCache}, dense::Val{D}, t0::T, x::Taylor1{U}, maxsteps::Int) where {D, U, T}
+    return cachetype(
         Array{T}(undef, maxsteps + 1),
         Array{U}(undef, maxsteps + 1),
         init_psol(dense, maxsteps, 1, x))
 end
 
-function init_cache(::Type{ScalarCache}, ::Val{false}, trange::AbstractVector{T}, x::Taylor1{U}, maxsteps::Int) where {U, T}
+function init_cache(cachetype::Type{ScalarCache}, ::Val{false}, trange::AbstractVector{T}, x::Taylor1{U}, maxsteps::Int) where {U, T}
     nn = length(trange)
-    cache = ScalarCache(
+    cache = cachetype(
         trange,
         Array{U}(undef, nn),
         init_psol(Val(false), maxsteps, 1, x))
@@ -48,9 +48,9 @@ struct VectorCache{TV, XV, PSOL, XAUX} <: AbstractVectorCache
     xaux::XAUX
 end
 
-function init_cache(::Type{VectorCache}, dense::Val{D}, t0::T, x::Vector{Taylor1{U}}, maxsteps::Int) where {D, U, T}
+function init_cache(cachetype::Type{VectorCache}, dense::Val{D}, t0::T, x::Vector{Taylor1{U}}, maxsteps::Int) where {D, U, T}
     dof = length(x)
-    return VectorCache(
+    return cachetype(
         Array{T}(undef, maxsteps + 1),
         Array{U}(undef, dof, maxsteps + 1),
         init_psol(dense, maxsteps, dof, x),
@@ -68,10 +68,10 @@ struct VectorTRangeCache{TV, XV, PSOL, XAUX, X0, X1} <: AbstractVectorCache
     x1::X1
 end
 
-function init_cache(::Type{VectorTRangeCache}, ::Val{false}, trange::AbstractVector{T}, x::Vector{Taylor1{U}}, maxsteps::Int) where {U, T}
+function init_cache(cachetype::Type{VectorTRangeCache}, ::Val{false}, trange::AbstractVector{T}, x::Vector{Taylor1{U}}, maxsteps::Int) where {U, T}
     nn = length(trange)
     dof = length(x)
-    cache = VectorTRangeCache(
+    cache = cachetype(
         trange,
         Array{U}(undef, dof, nn),
         init_psol(Val(false), maxsteps, dof, x),
@@ -106,10 +106,10 @@ struct LyapunovSpectrumCache{TV, XV, PSOL, XAUX, X0, Λ, ΛTSUM, ΔX, DΔX, JAC,
     vⱼ::VJ
 end
 
-function init_cache(::Type{LyapunovSpectrumCache}, dense, t0::T, x::Vector{Taylor1{U}}, maxsteps::Int) where {U, T}
+function init_cache(cachetype::Type{LyapunovSpectrumCache}, dense, t0::T, x::Vector{Taylor1{U}}, maxsteps::Int) where {U, T}
     nx0 = length(x) # equals dof + dof^2
     dof = Int(sqrt(nx0 + 1/4) - 1/2)
-    cache = LyapunovSpectrumCache(
+    cache = cachetype(
         Array{T}(undef, maxsteps+1),
         Array{U}(undef, dof, maxsteps+1),
         nothing,
@@ -153,11 +153,11 @@ struct LyapunovSpectrumTRangeCache{TV, XV, PSOL, XAUX, X0, Q1, Λ, ΛTSUM, ΔX,
     vⱼ::VJ
 end
 
-function init_cache(::Type{LyapunovSpectrumTRangeCache}, dense, trange::AbstractVector{T}, x::Vector{Taylor1{U}}, maxsteps::Int) where {U, T}
+function init_cache(cachetype::Type{LyapunovSpectrumTRangeCache}, dense, trange::AbstractVector{T}, x::Vector{Taylor1{U}}, maxsteps::Int) where {U, T}
     nx0 = length(x) # equals dof + dof^2
     dof = Int(sqrt(nx0 + 1/4) - 1/2)
     nn = length(trange)
-    cache = LyapunovSpectrumTRangeCache(
+    cache = cachetype(
         trange,
         Array{U}(undef, dof, nn),
         nothing,

diff --git a/src/integrator/taylorinteg.jl b/src/integrator/taylorinteg.jl
@@ -62,8 +62,8 @@ function taylorinteg(f, x0::U, t0::T, tmax::T, order::Int, abstol::T, params = n
     parse_eqs, rv = _determine_parsing!(parse_eqs, f, t, x, params)
 
     # Re-initialize the Taylor1 expansions
-    t = t0 + Taylor1( T, order )
-    x = Taylor1( x0, order )
+    # t = t0 + Taylor1( T, order )
+    # x = Taylor1( x0, order )
     return _taylorinteg!(Val(dense), f, t, x, x0, t0, tmax, abstol, rv, cache, params; parse_eqs, maxsteps)
 end
 
@@ -124,10 +124,6 @@ function taylorinteg(f!, q0::Array{U,1}, t0::T, tmax::T, order::Int, abstol::T,
     # Determine if specialized jetcoeffs! method exists
     parse_eqs, rv = _determine_parsing!(parse_eqs, f!, t, x, dx, params)
 
-    # Re-initialize the Taylor1 expansions
-    t = t0 + Taylor1( T, order )
-    x .= Taylor1.( q0, order )
-    dx .= Taylor1.( zero.(q0), order)
     return _taylorinteg!(Val(dense), f!, t, x, dx, q0, t0, tmax, abstol, rv,
         cache, params; parse_eqs, maxsteps)
 end
@@ -284,9 +280,6 @@ function taylorinteg(f, x0::U, trange::AbstractVector{T},
     # Determine if specialized jetcoeffs! method exists
     parse_eqs, rv = _determine_parsing!(parse_eqs, f, t, x, params)
 
-    # Re-initialize the Taylor1 expansions
-    t = t0 + Taylor1( T, order )
-    x = Taylor1( x0, order )
     return _taylorinteg!(f, t, x, x0, trange, abstol, rv, cache, params; parse_eqs, maxsteps)
 end
 
@@ -358,10 +351,6 @@ function taylorinteg(f!, q0::Array{U,1}, trange::AbstractVector{T},
     # Determine if specialized jetcoeffs! method exists
     parse_eqs, rv = _determine_parsing!(parse_eqs, f!, t, x, dx, params)
 
-    # Re-initialize the Taylor1 expansions
-    t = t0 + Taylor1( T, order )
-    x .= Taylor1.( q0, order )
-    dx .= Taylor1.( zero.(q0), order )
     return _taylorinteg!(f!, t, x, dx, q0, trange, abstol, rv,
         cache, params; parse_eqs, maxsteps)
 end

diff --git a/src/lyapunovspectrum.jl b/src/lyapunovspectrum.jl
@@ -278,9 +278,6 @@ function lyap_taylorinteg(f!, q0::Array{U,1}, t0::T, tmax::T,
     parse_eqs, rv = _determine_parsing!(parse_eqs, f!, t,
                     view(x, 1:dof), view(dx, 1:dof), params)
 
-    # Re-initialize the Taylor1 expansions
-    t = t0 + Taylor1( T, order )
-    x .= Taylor1.( x0, order )
     return _lyap_taylorinteg!(f!, t, x, dx, q0, t0, tmax, abstol, jt, _dv, rv,
         cache, params, jacobianfunc!; maxsteps, parse_eqs)
 end
@@ -381,10 +378,6 @@ function lyap_taylorinteg(f!, q0::Array{U,1}, trange::AbstractVector{T},
     parse_eqs, rv = _determine_parsing!(parse_eqs, f!, t,
                     view(x, 1:dof), view(dx, 1:dof), params)
 
-    # Re-initialize the Taylor1 expansions
-    t = trange[1] + Taylor1( T, order )
-    x .= Taylor1.( x0, order )
-    tmpTaylor, arrTaylor = rv.v0, rv.v1
     return _lyap_taylorinteg!(f!, t, x, dx, q0, trange, abstol, jt, _dv, rv, cache,
         params, jacobianfunc!; parse_eqs, maxsteps)
 

diff --git a/src/rootfinding.jl b/src/rootfinding.jl
@@ -223,18 +223,11 @@ function taylorinteg(f!, g, q0::Array{U,1}, t0::T, tmax::T,
     end
 
     # Allocation
-    cache = VectorCache(
-        Array{T}(undef, maxsteps + 1),
-        Array{U}(undef, dof, maxsteps + 1),
-        init_psol(Val(dense), maxsteps, dof, x),
-        Array{Taylor1{U}}(undef, dof))
+    cache = init_cache(VectorCache, Val(dense), t0, x, maxsteps)
 
     # Determine if specialized jetcoeffs! method exists
     parse_eqs, rv = _determine_parsing!(parse_eqs, f!, t, x, dx, params)
 
-    # Re-initialize the Taylor1 expansions
-    t = t0 + Taylor1( T, order )
-    x .= Taylor1.( q0, order )
     return _taylorinteg!(Val(dense), f!, g, t, x, dx, q0, t0, tmax, abstol, rv, cache, params;
         parse_eqs, maxsteps, eventorder, newtoniter, nrabstol)
 end
@@ -328,26 +321,11 @@ function taylorinteg(f!, g, q0::Array{U,1}, trange::AbstractVector{T},
     end
 
     # Allocation
-    nn = length(trange)
-    dof = length(q0)
-    cache = VectorTRangeCache(
-        trange,
-        Array{U}(undef, dof, nn),
-        init_psol(Val(false), maxsteps, dof, x),
-        Array{Taylor1{U}}(undef, dof),
-        similar(q0),
-        similar(q0))
-    fill!(cache.x0, T(NaN))
-    for ind in 1:nn
-        @inbounds cache.xv[:,ind] .= cache.x0
-    end
+    cache = init_cache(VectorTRangeCache, Val(false), trange, x, maxsteps)
 
     # Determine if specialized jetcoeffs! method exists
     parse_eqs, rv = _determine_parsing!(parse_eqs, f!, t, x, dx, params)
 
-    # Re-initialize the Taylor1 expansions
-    t = t0 + Taylor1( T, order )
-    x .= Taylor1.(q0, order)
     return _taylorinteg!(f!, g, t, x, dx, q0, trange, abstol, rv, cache, params;
         parse_eqs, maxsteps, eventorder, newtoniter, nrabstol)
 end