diff --git a/src/BoundaryValueDiffEq.jl b/src/BoundaryValueDiffEq.jl index d8abecd2..7adfa357 100644 --- a/src/BoundaryValueDiffEq.jl +++ b/src/BoundaryValueDiffEq.jl @@ -1,5 +1,3 @@ - __precompile__() - module BoundaryValueDiffEq using Reexport, LinearAlgebra, SparseArrays diff --git a/src/collocation.jl b/src/collocation.jl index d51b380e..c2921f41 100644 --- a/src/collocation.jl +++ b/src/collocation.jl @@ -2,7 +2,7 @@ function BVPSystem(fun, bc, p, x::Vector{T}, M::Integer, order) where T N = size(x,1) y = vector_alloc(T, M, N) - BVPSystem(order, M, N, fun, bc, p, x, y, vector_alloc(T, M, N), vector_alloc(T, M, N), eltype(y)(M)) + BVPSystem(order, M, N, fun, bc, p, x, y, vector_alloc(T, M, N), vector_alloc(T, M, N), eltype(y)(undef, M)) end # If user offers an intial guess diff --git a/src/solve.jl b/src/solve.jl index 16a8c34e..701fe2a5 100644 --- a/src/solve.jl +++ b/src/solve.jl @@ -7,7 +7,7 @@ function DiffEqBase.__solve(prob::BVProblem, alg::Shooting; kwargs...) # Form a root finding function. loss = function (resid, minimizer) uEltype = eltype(minimizer) - tspan = (uEltype(prob.tspan[1]),uEltype(prob.tspan[2])) + tspan = (convert(uEltype,prob.tspan[1]),convert(uEltype,prob.tspan[2])) tmp_prob = ODEProblem(prob.f,minimizer,tspan) sol = solve(tmp_prob,alg.ode_alg;kwargs...) bc(resid,sol,sol.prob.p,sol.t) @@ -41,7 +41,7 @@ function DiffEqBase.__solve(prob::BVProblem, alg::Union{GeneralMIRK,MIRK}; dt=0. tableau = constructMIRK(S) cache = alg_cache(alg, S) # Upper-level iteration - vec_y = Array{eltype(S.y[1])}(S.M*S.N) # Vector + vec_y = Array{eltype(S.y[1])}(undef, S.M*S.N) # Vector reorder! = function (resid) # reorder the Jacobian matrix such that it is banded tmp_last = resid[end] diff --git a/src/vector_auxiliary.jl b/src/vector_auxiliary.jl index 0aaccef0..8e56b4f6 100644 --- a/src/vector_auxiliary.jl +++ b/src/vector_auxiliary.jl @@ -1,6 +1,6 @@ # Auxiliary functions for working with vector of vectors function vector_alloc(T, M, N) - v = Vector{Vector{T}}(N) + v = Vector{Vector{T}}(undef, N) for i in eachindex(v) v[i] = zeros(T, M) end diff --git a/test/mirk_convergence_tests.jl b/test/mirk_convergence_tests.jl index 41b8dbdf..2a69bfca 100644 --- a/test/mirk_convergence_tests.jl +++ b/test/mirk_convergence_tests.jl @@ -10,7 +10,7 @@ end # Not able to change the initial condition. # Hard coded solution. -func_1!(::Type{Val{:analytic}}, u0, p, t) = [5-t,-1] +func_1 = ODEFunction(func_1!, analytic=(u0, p, t) -> [5-t,-1]) function boundary!(residual, u, p, t) residual[1] = u[1][1]-5 @@ -34,16 +34,17 @@ end # Hard coded solution. func_2!(::Type{Val{:analytic}}, u0, p, t) = [5*(cos(t) - cot(5)*sin(t)), 5*(-cos(t)*cot(5) - sin(t))] - +func_2 = ODEFunction(func_2!, analytic=(u0, p, t) -> [5*(cos(t) - cot(5)*sin(t)), + 5*(-cos(t)*cot(5) - sin(t))]) tspan = (0.,5.) u0 = [5.,-3.5] -probArr = [BVProblem(func_1!, boundary!, u0, tspan), - BVProblem(func_2!, boundary!, u0, tspan), - TwoPointBVProblem(func_1!, boundary_two_point!, u0, tspan), - TwoPointBVProblem(func_2!, boundary_two_point!, u0, tspan)] +probArr = [BVProblem(func_1, boundary!, u0, tspan), + BVProblem(func_2, boundary!, u0, tspan), + TwoPointBVProblem(func_1, boundary_two_point!, u0, tspan), + TwoPointBVProblem(func_2, boundary_two_point!, u0, tspan)] testTol = 0.2 -affineTol = 1e-6 +affineTol = 1e-2 dts = (1/2) .^ (5:-1:1) order = 4 @@ -51,7 +52,7 @@ println("Collocation method (GeneralMIRK)") println("Affineness Test") prob = probArr[1] @time sol = solve(prob, GeneralMIRK4(), dt=0.2) -@test norm(diff(map(x->x[1], sol.u)) + 0.2, Inf) + abs(sol[1][1]-5) < affineTol +@test norm(diff(map(x->x[1], sol.u)) .+ 0.2, Inf) + abs(sol[1][1]-5) < affineTol println("Convergence Test on Linear") prob = probArr[2] @@ -62,7 +63,7 @@ println("Collocation method (MIRK)") println("Affineness Test") prob = probArr[3] @time sol = solve(prob, MIRK4(), dt=0.2) -@test norm(diff(map(x->x[1], sol.u)) + 0.2, Inf) + abs(sol[1][1]-5) < affineTol +@test norm(diff(map(x->x[1], sol.u)) .+ 0.2, Inf) .+ abs(sol[1][1]-5) < affineTol println("Convergence Test on Linear") prob = probArr[4] diff --git a/test/orbital.jl b/test/orbital.jl index 5ece1ea1..21768926 100644 --- a/test/orbital.jl +++ b/test/orbital.jl @@ -57,7 +57,7 @@ cur_bc!(resid_f,sol,nothing,sol.t) cur_bc!(resid_f,sol,nothing,sol.t) @test norm(resid_f, Inf) < TestTol -@time sol = solve(bvp, Shooting(DP5(),nlsolve=(f,u0) -> (res=NLsolve.nlsolve(f,u0,autodiff=true,ftol=1e-13,xtol=1e-13); +@time sol = solve(bvp, Shooting(DP5(),nlsolve=(f,u0) -> (res=NLsolve.nlsolve(f,u0,autodiff=:forward,ftol=1e-13,xtol=1e-13); (res.zero, res.f_converged))),force_dtmin=true,abstol=1e-13,reltol=1e-13) cur_bc!(resid_f,sol,nothing,sol.t) @test norm(resid_f, Inf) < TestTol