Skip to content
This repository has been archived by the owner on Jul 19, 2023. It is now read-only.

(WIP) Safer handling of DerivativeOperator.coefficients #244

Open
wants to merge 12 commits into
base: master
Choose a base branch
from

Conversation

ghost
Copy link

@ghost ghost commented May 6, 2020

This fixes bugs related to referencing undefined memory in the coefficients field of DerivativeOperator (see #240).

I found that the code / docs / tests had inconsistent notions on the type of coefficients. There were also differences between UpwindDifference and CenteredDifference that didn't have obvious motivation, and weren't mentioned in the docs. I hope we can agree on something consistent and include it in this PR.

Currently, I'm hitting a breaking issue dispatching on convolve_BC_left! (breaks the tests). For a MWE, change examples/heat_equation.jl to use UpwindDifference instead of CenteredDifference.

Based on the stacktrace, I can't figure out why it won't choose the method at convolutions.jl:106: function convolve_BC_left!(x_temp::AbstractVector{T}, x::AbstractVector{T}, A::DerivativeOperator{T,N,true}; overwrite = true) where {T<:Real, N}.

I check the type of each arg in the debug log below, then trigger the dispatch error:

About to run: <(DiffEqOperators.#mul!#69)(true, LinearAlgebra.mul!, ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDer...>
2|julia> typeof(x_temp) <: AbstractVector{T} where T <: Real
true

2|julia> typeof(x) <: BoundaryPaddedVector
true

2|julia> typeof(A) <: DerivativeOperator{T,N,true} where {T<:Real, N}
true
...

About to run: <(DiffEqOperators.var"#convolve_BC_left!##kw"())((overwrite = true,), DiffEqOperators.convolve_BC_left...>
1|debug> s
ERROR: MethodError: no method matching convolve_BC_left!(::Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1}, ::DiffEqOperators.BoundaryPaddedVector{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1}}, ::DerivativeOperator{Float64,1,true,Float64,StaticArrays.SArray{Tuple{4},Float64,1,4},StaticArrays.SArray{Tuple{2},StaticArrays.SArray{Tuple{4},Float64,1,4},1,2},Nothing,Nothing}; overwrite=true)
Closest candidates are:
  convolve_BC_left!(::AbstractArray{T1,1}, ::DiffEqOperators.BoundaryPaddedVector, ::DerivativeOperator{T2,N,false,T21,S1,S2,T3,F} where F where T3 where S2<:StaticArrays.SArray where S1 where T21; overwrite) where {T1, T2, N} at /home/dan/repos/DiffEqOperators.jl/src/derivative_operators/convolutions.jl:261
  convolve_BC_left!(::AbstractArray{T1,1}, ::AbstractArray{T2,1}, ::DerivativeOperator{T3,N,false,T2,S1,S2,T31,F} where F where T31 where S2<:StaticArrays.SArray where S1 where T2; overwrite) where {T1, T2, T3, N} at /home/dan/repos/DiffEqOperators.jl/src/derivative_operators/convolutions.jl:55
  convolve_BC_left!(::AbstractArray{T,1}, ::DiffEqOperators.BoundaryPaddedVector, ::DerivativeOperator{T,N,true,M,S1,S2,T3,F} where F where T3 where S2<:StaticArrays.SArray where S1; overwrite) where {T<:Real, N, M<:(AbstractArray{T,N} where N)} at /home/dan/repos/DiffEqOperators.jl/src/derivative_operators/convolutions.jl:443
  ...
Stacktrace:
 [1] mul!(::Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1}, ::DerivativeOperator{Float64,1,true,Float64,StaticArrays.SArray{Tuple{4},Float64,1,4},StaticArrays.SArray{Tuple{2},StaticArrays.SArray{Tuple{4},Float64,1,4},1,2},Nothing,Nothing}, ::DiffEqOperators.BoundaryPaddedVector{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1}}; overwrite::Bool) at /home/dan/repos/DiffEqOperators.jl/src/derivative_operators/convolutions.jl:18
 [2] mul!(::Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1}, ::DerivativeOperator{Float64,1,true,Float64,StaticArrays.SArray{Tuple{4},Float64,1,4},StaticArrays.SArray{Tuple{2},StaticArrays.SArray{Tuple{4},Float64,1,4},1,2},Nothing,Nothing}, ::DiffEqOperators.BoundaryPaddedVector{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1}}) at /home/dan/repos/DiffEqOperators.jl/src/derivative_operators/convolutions.jl:18
 [3] mul!(::Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1}, ::GhostDerivativeOperator{Float64,DerivativeOperator{Float64,1,true,Float64,StaticArrays.SArray{Tuple{4},Float64,1,4},StaticArrays.SArray{Tuple{2},StaticArrays.SArray{Tuple{4},Float64,1,4},1,2},Nothing,Nothing},RobinBC{Float64,Array{Float64,1}}}, ::Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1}) at /home/dan/repos/DiffEqOperators.jl/src/derivative_operators/ghost_derivative_operator.jl:16
 [4] *(::GhostDerivativeOperator{Float64,DerivativeOperator{Float64,1,true,Float64,StaticArrays.SArray{Tuple{4},Float64,1,4},StaticArrays.SArray{Tuple{2},StaticArrays.SArray{Tuple{4},Float64,1,4},1,2},Nothing,Nothing},RobinBC{Float64,Array{Float64,1}}}, ::Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1}) at /home/dan/repos/DiffEqOperators.jl/src/derivative_operators/ghost_derivative_operator.jl:22
 [5] *(::DerivativeOperator{Float64,1,true,Float64,StaticArrays.SArray{Tuple{4},Float64,1,4},StaticArrays.SArray{Tuple{2},StaticArrays.SArray{Tuple{4},Float64,1,4},1,2},Nothing,Nothing}, ::RobinBC{Float64,Array{Float64,1}}, ::Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1}, ::Tuple{}) at operators.jl:529
 [6] step(::Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1}, ::DiffEqBase.NullParameters, ::Float64) at /tmp/julia-tmp8xGaN7:14
 [7] (::ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing})(::Tuple{Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1},DiffEqBase.NullParameters,Float64}) at /home/dan/.julia/packages/DiffEqBase/Avuz1/src/diffeqfunction.jl:248
 [8] (::DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters})(::Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1}) at /home/dan/.julia/packages/DiffEqBase/Avuz1/src/function_wrappers.jl:30
 [9] forwarddiff_color_jacobian(::DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}, ::Array{Float64,1}, ::SparseDiffTools.ForwardColorJacCache{Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1},Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters},Float64}(),Float64,12},1},Array{Float64,1},Array{Array{NTuple{12,Bool},1},1},UnitRange{Int64},Nothing}, ::Nothing) at /home/dan/.julia/packages/SparseDiffTools/Oxsav/src/differentiation/compute_jacobian_ad.jl:100
 [10] forwarddiff_color_jacobian(::DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}, ::Array{Float64,1}; dx::Array{Float64,1}, colorvec::UnitRange{Int64}, sparsity::Nothing, jac_prototype::Nothing, chunksize::Nothing) at /home/dan/.julia/packages/SparseDiffTools/Oxsav/src/differentiation/compute_jacobian_ad.jl:73
 [11] (::SparseDiffTools.var"#forwarddiff_color_jacobian##kw")(::NamedTuple{(:colorvec, :sparsity, :jac_prototype, :chunksize),Tuple{UnitRange{Int64},Nothing,Nothing,Nothing}}, ::typeof(SparseDiffTools.forwarddiff_color_jacobian), ::DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}, ::Array{Float64,1}) at /home/dan/.julia/packages/SparseDiffTools/Oxsav/src/differentiation/compute_jacobian_ad.jl:70
 [12] jacobian_autodiff(::DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}, ::Array{Float64,1}, ::ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}, ::KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType}) at /home/dan/.julia/packages/OrdinaryDiffEq/Pvuwl/src/derivative_wrappers.jl:48
 [13] jacobian(::DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}, ::Array{Float64,1}, ::OrdinaryDiffEq.ODEIntegrator{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},false,Array{Float64,1},Nothing,Float64,DiffEqBase.NullParameters,Float64,Float64,Float64,Array{Array{Float64,1},1},ODESolution{Float64,2,Array{Array{Float64,1},1},Nothing,Nothing,Array{Float64,1},Array{Array{Array{Float64,1},1},1},ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem},KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},OrdinaryDiffEq.InterpolationData{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Array{Float64,1},1},Array{Float64,1},Array{Array{Array{Float64,1},1},1},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}}},DiffEqBase.DEStats},ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}},OrdinaryDiffEq.DEOptions{Float64,Float64,Float64,Float64,typeof(DiffEqBase.ODE_DEFAULT_NORM),typeof(LinearAlgebra.opnorm),CallbackSet{Tuple{},Tuple{}},typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN),typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE),typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK),DataStructures.BinaryHeap{Float64,DataStructures.LessThan},DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Nothing,Nothing,Int64,Array{Float64,1},Array{Float64,1},Array{Float64,1}},Array{Float64,1},Float64,Nothing,OrdinaryDiffEq.DefaultInit}) at /home/dan/.julia/packages/OrdinaryDiffEq/Pvuwl/src/derivative_wrappers.jl:73
 [14] calc_J(::OrdinaryDiffEq.ODEIntegrator{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},false,Array{Float64,1},Nothing,Float64,DiffEqBase.NullParameters,Float64,Float64,Float64,Array{Array{Float64,1},1},ODESolution{Float64,2,Array{Array{Float64,1},1},Nothing,Nothing,Array{Float64,1},Array{Array{Array{Float64,1},1},1},ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem},KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},OrdinaryDiffEq.InterpolationData{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Array{Float64,1},1},Array{Float64,1},Array{Array{Array{Float64,1},1},1},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}}},DiffEqBase.DEStats},ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}},OrdinaryDiffEq.DEOptions{Float64,Float64,Float64,Float64,typeof(DiffEqBase.ODE_DEFAULT_NORM),typeof(LinearAlgebra.opnorm),CallbackSet{Tuple{},Tuple{}},typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN),typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE),typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK),DataStructures.BinaryHeap{Float64,DataStructures.LessThan},DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Nothing,Nothing,Int64,Array{Float64,1},Array{Float64,1},Array{Float64,1}},Array{Float64,1},Float64,Nothing,OrdinaryDiffEq.DefaultInit}, ::OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}) at /home/dan/.julia/packages/OrdinaryDiffEq/Pvuwl/src/derivative_utils.jl:68
 [15] calc_W(::OrdinaryDiffEq.ODEIntegrator{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},false,Array{Float64,1},Nothing,Float64,DiffEqBase.NullParameters,Float64,Float64,Float64,Array{Array{Float64,1},1},ODESolution{Float64,2,Array{Array{Float64,1},1},Nothing,Nothing,Array{Float64,1},Array{Array{Array{Float64,1},1},1},ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem},KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},OrdinaryDiffEq.InterpolationData{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Array{Float64,1},1},Array{Float64,1},Array{Array{Array{Float64,1},1},1},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}}},DiffEqBase.DEStats},ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}},OrdinaryDiffEq.DEOptions{Float64,Float64,Float64,Float64,typeof(DiffEqBase.ODE_DEFAULT_NORM),typeof(LinearAlgebra.opnorm),CallbackSet{Tuple{},Tuple{}},typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN),typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE),typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK),DataStructures.BinaryHeap{Float64,DataStructures.LessThan},DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Nothing,Nothing,Int64,Array{Float64,1},Array{Float64,1},Array{Float64,1}},Array{Float64,1},Float64,Nothing,OrdinaryDiffEq.DefaultInit}, ::OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}, ::Float64, ::Bool, ::Bool) at /home/dan/.julia/packages/OrdinaryDiffEq/Pvuwl/src/derivative_utils.jl:465
 [16] update_W!(::OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}}, ::OrdinaryDiffEq.ODEIntegrator{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},false,Array{Float64,1},Nothing,Float64,DiffEqBase.NullParameters,Float64,Float64,Float64,Array{Array{Float64,1},1},ODESolution{Float64,2,Array{Array{Float64,1},1},Nothing,Nothing,Array{Float64,1},Array{Array{Array{Float64,1},1},1},ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem},KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},OrdinaryDiffEq.InterpolationData{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Array{Float64,1},1},Array{Float64,1},Array{Array{Array{Float64,1},1},1},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}}},DiffEqBase.DEStats},ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}},OrdinaryDiffEq.DEOptions{Float64,Float64,Float64,Float64,typeof(DiffEqBase.ODE_DEFAULT_NORM),typeof(LinearAlgebra.opnorm),CallbackSet{Tuple{},Tuple{}},typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN),typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE),typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK),DataStructures.BinaryHeap{Float64,DataStructures.LessThan},DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Nothing,Nothing,Int64,Array{Float64,1},Array{Float64,1},Array{Float64,1}},Array{Float64,1},Float64,Nothing,OrdinaryDiffEq.DefaultInit}, ::OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}}, ::Float64, ::Bool) at /home/dan/.julia/packages/OrdinaryDiffEq/Pvuwl/src/derivative_utils.jl:509
 [17] nlsolve!(::OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}}, ::OrdinaryDiffEq.ODEIntegrator{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},false,Array{Float64,1},Nothing,Float64,DiffEqBase.NullParameters,Float64,Float64,Float64,Array{Array{Float64,1},1},ODESolution{Float64,2,Array{Array{Float64,1},1},Nothing,Nothing,Array{Float64,1},Array{Array{Array{Float64,1},1},1},ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem},KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},OrdinaryDiffEq.InterpolationData{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Array{Float64,1},1},Array{Float64,1},Array{Array{Array{Float64,1},1},1},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}}},DiffEqBase.DEStats},ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}},OrdinaryDiffEq.DEOptions{Float64,Float64,Float64,Float64,typeof(DiffEqBase.ODE_DEFAULT_NORM),typeof(LinearAlgebra.opnorm),CallbackSet{Tuple{},Tuple{}},typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN),typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE),typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK),DataStructures.BinaryHeap{Float64,DataStructures.LessThan},DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Nothing,Nothing,Int64,Array{Float64,1},Array{Float64,1},Array{Float64,1}},Array{Float64,1},Float64,Nothing,OrdinaryDiffEq.DefaultInit}, ::OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}}, ::Bool) at /home/dan/.julia/packages/OrdinaryDiffEq/Pvuwl/src/nlsolve/nlsolve.jl:14
 [18] perform_step!(::OrdinaryDiffEq.ODEIntegrator{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},false,Array{Float64,1},Nothing,Float64,DiffEqBase.NullParameters,Float64,Float64,Float64,Array{Array{Float64,1},1},ODESolution{Float64,2,Array{Array{Float64,1},1},Nothing,Nothing,Array{Float64,1},Array{Array{Array{Float64,1},1},1},ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem},KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},OrdinaryDiffEq.InterpolationData{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Array{Float64,1},1},Array{Float64,1},Array{Array{Array{Float64,1},1},1},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}}},DiffEqBase.DEStats},ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}},OrdinaryDiffEq.DEOptions{Float64,Float64,Float64,Float64,typeof(DiffEqBase.ODE_DEFAULT_NORM),typeof(LinearAlgebra.opnorm),CallbackSet{Tuple{},Tuple{}},typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN),typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE),typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK),DataStructures.BinaryHeap{Float64,DataStructures.LessThan},DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Nothing,Nothing,Int64,Array{Float64,1},Array{Float64,1},Array{Float64,1}},Array{Float64,1},Float64,Nothing,OrdinaryDiffEq.DefaultInit}, ::OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}}, ::Bool) at /home/dan/.julia/packages/OrdinaryDiffEq/Pvuwl/src/perform_step/kencarp_kvaerno_perform_step.jl:831
 [19] perform_step!(::OrdinaryDiffEq.ODEIntegrator{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},false,Array{Float64,1},Nothing,Float64,DiffEqBase.NullParameters,Float64,Float64,Float64,Array{Array{Float64,1},1},ODESolution{Float64,2,Array{Array{Float64,1},1},Nothing,Nothing,Array{Float64,1},Array{Array{Array{Float64,1},1},1},ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem},KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},OrdinaryDiffEq.InterpolationData{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Array{Float64,1},1},Array{Float64,1},Array{Array{Array{Float64,1},1},1},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}}},DiffEqBase.DEStats},ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}},OrdinaryDiffEq.DEOptions{Float64,Float64,Float64,Float64,typeof(DiffEqBase.ODE_DEFAULT_NORM),typeof(LinearAlgebra.opnorm),CallbackSet{Tuple{},Tuple{}},typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN),typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE),typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK),DataStructures.BinaryHeap{Float64,DataStructures.LessThan},DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Nothing,Nothing,Int64,Array{Float64,1},Array{Float64,1},Array{Float64,1}},Array{Float64,1},Float64,Nothing,OrdinaryDiffEq.DefaultInit}, ::OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}}) at /home/dan/.julia/packages/OrdinaryDiffEq/Pvuwl/src/perform_step/kencarp_kvaerno_perform_step.jl:784
 [20] solve!(::OrdinaryDiffEq.ODEIntegrator{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},false,Array{Float64,1},Nothing,Float64,DiffEqBase.NullParameters,Float64,Float64,Float64,Array{Array{Float64,1},1},ODESolution{Float64,2,Array{Array{Float64,1},1},Nothing,Nothing,Array{Float64,1},Array{Array{Array{Float64,1},1},1},ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem},KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType},OrdinaryDiffEq.InterpolationData{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Array{Float64,1},1},Array{Float64,1},Array{Array{Array{Float64,1},1},1},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}}},DiffEqBase.DEStats},ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},OrdinaryDiffEq.KenCarp4ConstantCache{OrdinaryDiffEq.NLSolver{NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},false,Array{Float64,1},Float64,OrdinaryDiffEq.NLNewtonConstantCache{Float64,Float64,Array{Float64,2},LinearAlgebra.LU{Float64,Array{Float64,2}},DiffEqBase.UDerivativeWrapper{ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Float64,DiffEqBase.NullParameters}}},OrdinaryDiffEq.KenCarp4Tableau{Float64,Float64}},OrdinaryDiffEq.DEOptions{Float64,Float64,Float64,Float64,typeof(DiffEqBase.ODE_DEFAULT_NORM),typeof(LinearAlgebra.opnorm),CallbackSet{Tuple{},Tuple{}},typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN),typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE),typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK),DataStructures.BinaryHeap{Float64,DataStructures.LessThan},DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Nothing,Nothing,Int64,Array{Float64,1},Array{Float64,1},Array{Float64,1}},Array{Float64,1},Float64,Nothing,OrdinaryDiffEq.DefaultInit}) at /home/dan/.julia/packages/OrdinaryDiffEq/Pvuwl/src/solve.jl:424
 [21] __solve(::ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem}, ::KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType}, ::Tuple{}; kwargs::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}) at /home/dan/.julia/packages/OrdinaryDiffEq/Pvuwl/src/solve.jl:5
 [22] __solve(::ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem}, ::KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType}, ::Tuple{}) at /home/dan/.julia/packages/OrdinaryDiffEq/Pvuwl/src/solve.jl:4
 [23] (::DiffEqBase.var"#466#467"{ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem},Tuple{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType}}})() at /home/dan/.julia/packages/DiffEqBase/Avuz1/src/solve.jl:49
 [24] maybe_with_logger(::DiffEqBase.var"#466#467"{ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem},Tuple{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType}}}, ::Nothing) at /home/dan/.julia/packages/DiffEqBase/Avuz1/src/utils.jl:259
 [25] solve_call(::ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem}, ::Tuple{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType}}; merge_callbacks::Bool, kwargs::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}) at /home/dan/.julia/packages/DiffEqBase/Avuz1/src/solve.jl:48
 [26] solve_call(::ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem}, ::Tuple{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType}}) at /home/dan/.julia/packages/DiffEqBase/Avuz1/src/solve.jl:37
 [27] solve(::ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem}, ::Tuple{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType}}; kwargs::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}) at /home/dan/.julia/packages/DiffEqBase/Avuz1/src/solve.jl:66
 [28] solve(::ODEProblem{Array{Float64,1},Tuple{Float64,Float64},false,DiffEqBase.NullParameters,ODEFunction{false,typeof(step),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem}, ::Tuple{KenCarp4{0,true,DefaultLinSolve,NLNewton{Rational{Int64},Rational{Int64},Rational{Int64}},DataType}}) at /home/dan/.julia/packages/DiffEqBase/Avuz1/src/solve.jl:54
 [29] ##thunk#261() at REPL[7]:1

"""
get_coefficient(coefficients::AbstractVector, index::Int) = coeff[i]

# FIXME: I don't think this case is used anymore
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This case was covered by a subset of convolution functions, but I don't know what the use-case is.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It's for a constant coefficient, like D*A so a constant diffusion. If the coefficient is a number, then the whole setup is actually a bitstype so it'll compile to the GPU and distribute much better, and this is a common case so it's worth supporting. I think that one function is really all that's necessary to support it.

if coeff_func != nothing
compute_coeffs!(coeff_func, coefficients)
end
coefficients = init_coefficients(coeff_func, len)
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This fixes #240

@ghost
Copy link
Author

ghost commented May 6, 2020

get_coefficient(coefficients::Number, index::Int) = coefficients

# FIXME: Why use "true" here for the value 1?
get_coefficient(coefficients::Nothing, index::Int) = true
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

true is a hard 1 in Julia. true will convert to any type, so it's the recommended value if you want it to promote in further operations.

@@ -336,7 +336,7 @@ function LinearAlgebra.Array(A::DerivativeOperator{T,N,true}, len::Int=A.len) wh
end

for i in len-bpc+1:len
cur_coeff = coeff[i]
cur_coeff = get_coefficient(coeff, i)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

very clean. Beautiful.

@ChrisRackauckas
Copy link
Member

These changes are looking great. It's always nice to see a code reduction that keeps functionality, fixes bugs, and is easier to read.

@@ -13,7 +13,7 @@ function Base.copyto!(L::AbstractMatrix{T}, A::DerivativeOperator{T}, N::Int) wh

coeff = A.coefficients
get_coeff = if coeff isa AbstractVector
i -> coeff[i]
i = get_coefficient(coeff, i)
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do you reference the variable i before it's defined here? Get an error julia> Array(CenteredDifference{1}(2, 3, 0.01, 100, c_squared)) ERROR: UndefVarError: i not defined Stacktrace: [1] copyto!(::Array{Float64,2}, ::DerivativeOperator{Float64,1,false,Float64,StaticArrays.SArray{Tuple{5},Float64,1,5},StaticArrays.SArray{Tuple{1},StaticArrays.SArray{Tuple{5},Float64,1,5},1,1},Array{Float64,1},Array{Float64,1}}, ::Int64) at /home/simen/Documents/school/masteroppgave/wave_simulation/dev/DiffEqOperators/src/derivative_operators/concretization.jl:16 [2] Array(::DerivativeOperator{Float64,1,false,Float64,StaticArrays.SArray{Tuple{5},Float64,1,5},StaticArrays.SArray{Tuple{1},StaticArrays.SArray{Tuple{5},Float64,1,5},1,1},Array{Float64,1},Array{Float64,1}}, ::Int64) at /home/simen/Documents/school/masteroppgave/wave_simulation/dev/DiffEqOperators/src/derivative_operators/concretization.jl:45 (repeats 2 times) [3] top-level scope at REPL[2]:1 when using a locally merged version of this pull request.

Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Changed te = to -> in line 16 as the following lines, which works for me.

Sign up for free to subscribe to this conversation on GitHub. Already have an account? Sign in.
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants