Skip to content
Closed
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -27,6 +27,7 @@ NLsolve = "2774e3e8-f4cf-5e23-947b-6d7e65073b56"
NonlinearSolve = "8913a72c-1f9b-4ce2-8d82-65094dcecaec"
Polyester = "f517fe37-dbe3-4b94-8317-1923a5111588"
PreallocationTools = "d236fae5-4411-538c-8e31-a6e3d9e00b46"
Preferences = "21216c6a-2e73-6563-6e65-726566657250"
RecursiveArrayTools = "731186ca-8d62-57ce-b412-fbd966d074cd"
Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
Expand Down Expand Up @@ -58,6 +59,7 @@ NLsolve = "4.3"
NonlinearSolve = "0.3.14"
Polyester = "0.3, 0.4, 0.5, 0.6"
PreallocationTools = "0.2, 0.3, 0.4"
Preferences = "1.3"
RecursiveArrayTools = "2.26.3"
Reexport = "0.2, 1.0"
SciMLBase = "1.50"
Expand Down
187 changes: 108 additions & 79 deletions src/OrdinaryDiffEq.jl
Original file line number Diff line number Diff line change
Expand Up @@ -83,6 +83,7 @@ using DocStringExtensions

import ArrayInterfaceStaticArrays, ArrayInterfaceGPUArrays
import FunctionWrappersWrappers
import Preferences

DEFAULT_PRECS(W,du,u,p,t,newW,Plprev,Prprev,solverdata) = nothing,nothing

Expand Down Expand Up @@ -192,90 +193,118 @@ using DocStringExtensions

import SnoopPrecompile

SnoopPrecompile.@precompile_all_calls begin
function lorenz(du,u,p,t)
du[1] = 10.0(u[2]-u[1])
du[2] = u[1]*(28.0-u[3]) - u[2]
du[3] = u[1]*u[2] - (8/3)*u[3]
end

function lorenz_oop(u,p,t)
[10.0(u[2]-u[1]),u[1]*(28.0-u[3]) - u[2],u[1]*u[2] - (8/3)*u[3]]
end

solver_list = [
BS3(), Tsit5(), Vern7(), Vern9(),

Rosenbrock23(), Rosenbrock23(autodiff=false),
Rosenbrock23(chunk_size = 1), Rosenbrock23(chunk_size = Val{1}()),

Rodas4(), Rodas4(autodiff=false),
#Rodas4(chunk_size = 1), Rodas4(chunk_size = Val{1}()),

Rodas5(), Rodas5(autodiff=false),
#Rodas5(chunk_size = 1), Rodas5(chunk_size = Val{1}()),

Rodas5P(), Rodas5P(autodiff=false),
Rodas5P(chunk_size = 1), Rodas5P(chunk_size = Val{1}()),

TRBDF2(), TRBDF2(autodiff=false),
#TRBDF2(chunk_size = 1), TRBDF2(chunk_size = Val{1}()),

KenCarp4(), KenCarp4(autodiff=false),
#KenCarp4(chunk_size = 1), KenCarp4(chunk_size = Val{1}()),

QNDF(), QNDF(autodiff=false),
#QNDF(chunk_size = 1), QNDF(chunk_size = Val{1}()),

AutoTsit5(Rosenbrock23()), AutoTsit5(Rosenbrock23(autodiff=false)),
AutoTsit5(Rosenbrock23(chunk_size = 1)),
AutoTsit5(Rosenbrock23(chunk_size = Val{1}())),

AutoTsit5(TRBDF2()), AutoTsit5(TRBDF2(autodiff=false)),
#AutoTsit5(TRBDF2(chunk_size = 1)),
#AutoTsit5(TRBDF2(chunk_size = Val{1}())),

AutoVern9(KenCarp47()), AutoVern9(KenCarp47(autodiff=false)),
#AutoVern9(KenCarp47(chunk_size = 1)),
#AutoVern9(KenCarp47(chunk_size = Val{1}())),

AutoVern9(Rodas5()), AutoVern9(Rodas5(autodiff=false)),
AutoVern9(Rodas5(chunk_size = 1)),
AutoVern9(Rodas5(chunk_size = Val{1}())),

AutoVern9(Rodas5P()), AutoVern9(Rodas5P(autodiff=false)),
AutoVern9(Rodas5P(chunk_size = 1)),
AutoVern9(Rodas5P(chunk_size = Val{1}())),

AutoVern7(Rodas4()), AutoVern7(Rodas4(autodiff=false)),
#AutoVern7(Rodas4(chunk_size = 1)),
#AutoVern7(Rodas4(chunk_size = Val{1}())),

#AutoVern7(Rodas5P()), AutoVern7(Rodas5P(autodiff=false)),
#AutoVern7(Rodas5P(chunk_size = 1)),
#AutoVern7(Rodas5P(chunk_size = Val{1}())),

AutoVern7(TRBDF2()), AutoVern7(TRBDF2(autodiff=false)),
#AutoVern7(TRBDF2(chunk_size = 1)),
#AutoVern7(TRBDF2(chunk_size = Val{1}())),
]
SnoopPrecompile.@precompile_all_calls begin
function lorenz(du, u, p, t)
du[1] = 10.0(u[2] - u[1])
du[2] = u[1] * (28.0 - u[3]) - u[2]
du[3] = u[1] * u[2] - (8 / 3) * u[3]
end

prob_list = [
ODEProblem(lorenz,[1.0;0.0;0.0],(0.0,1.0))
ODEProblem{true,false}(lorenz,[1.0;0.0;0.0],(0.0,1.0))
ODEProblem{true,false}(lorenz,[1.0;0.0;0.0],(0.0,1.0),Float64[])
ODEProblem(lorenz_oop,[1.0;0.0;0.0],(0.0,1.0))
#ODEProblem{false,false}(lorenz_oop,[1.0;0.0;0.0],(0.0,1.0))
#ODEProblem{false,false}(lorenz_oop,[1.0;0.0;0.0],(0.0,1.0),Float64[])
]
function lorenz_oop(u, p, t)
[10.0(u[2] - u[1]), u[1] * (28.0 - u[3]) - u[2], u[1] * u[2] - (8 / 3) * u[3]]
end
nonstiff_solver_options =
["BS3" => BS3(), "Tsit5" => Tsit5(), "Vern7" => Vern7(), "Vern9" => Vern9()]
stiff_solver_options = [
("Rosenbrock23", "true;true;true;true", Rosenbrock23),
("Rodas4", "true;true;false;false", Rodas4),
("Rodas5", "true;true;false;false", Rodas5),
("Rodas5P", "true;true;true;true", Rodas5P),
("TRBDF2", "true;true;false;false", TRBDF2),
("KenCarp4", "true;true;false;false", KenCarp4),
("QNDF", "true;true;false;false", QNDF),
]
autoswitch_solver_options = [
("AutoTsit5Rosenbrock23", "true;true;true;true", AutoTsit5, Rosenbrock23),
("AutoTsit5TRBDF2", "true;true;false;false", AutoTsit5, Rosenbrock23),
("AutoTsit5Rodas5P", "false;false;false;false", AutoTsit5, Rodas5P),
("AutoVern9KenCarp47", "true;true;false;false", AutoVern9, KenCarp47),
("AutoVern9Rodas5", "true;true;true;true", AutoVern9, Rodas5),
("AutoVern9Rodas5P", "true;true;true;true", AutoVern9, Rodas5P),
("AutoVern7Rodas4", "true;true;false;false", AutoVern7, Rodas4),
("AutoVern7Rodas5P", "false;false;false;false", AutoVern7, Rodas5P),
("AutoVern7TRBDF2", "true;true;false;false", AutoVern7, TRBDF2),
]


solver_list = []
for (solvername, solver) in nonstiff_solver_options
if Preferences.@load_preference(solvername, "true") == "true"
Copy link
Member

Choose a reason for hiding this comment

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

A bit simpler would be

Suggested change
if Preferences.@load_preference(solvername, "true") == "true"
if Preferences.@load_preference(solvername, true)

or

Suggested change
if Preferences.@load_preference(solvername, "true") == "true"
if Preferences.@load_preference(solvername, true) === true

Copy link
Contributor Author

@chriselrod chriselrod Aug 22, 2022

Choose a reason for hiding this comment

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

I haven't used or experimented with preferences much, but the README gave me the impression it could only load/save strings, while default could be anything.

Copy link
Member

Choose a reason for hiding this comment

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

ForwardDiff e.g. uses Bool and Int: https://github.com/JuliaDiff/ForwardDiff.jl/blob/61e4dd486a67cd0bdc6f61e58e452afcea388b7c/src/prelude.jl#L2-L3 I think everything that is supported by TOML might work (e.g., tuples didn't work but arrays did).

Copy link
Contributor Author

@chriselrod chriselrod Aug 22, 2022

Choose a reason for hiding this comment

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

Okay, how can one actually set these preferences to try without deving ForwardDiff?
My understanding is that @set_preferences! sets preferences with respect to the module in which it is located, meaning there has to be some function within ForwardDiff that includes @set_preferences!. As there are not, it is impossible to actually set?

Copy link
Member

Choose a reason for hiding this comment

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

@set_preferences! uses the current module but users would want to use set_preferences!: https://juliadiff.org/ForwardDiff.jl/dev/user/advanced/#Fixing-NaN/Inf-Issues There the first argument is the module or the UUID of the project for which you want to set preferences.

Copy link
Contributor Author

@chriselrod chriselrod Aug 22, 2022

Choose a reason for hiding this comment

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

Importantly (for the diffeq ecosystem), it also accepts the uuid, thus we can set_preferences! without first loading the module.

push!(solver_list, solver)
end
end
for (solvername, default, solvertype) in stiff_solver_options
options = split(Preferences.@load_preference(solvername, default), ';')
Copy link
Member

Choose a reason for hiding this comment

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

A bit simpler here would also be

Suggested change
options = split(Preferences.@load_preference(solvername, default), ';')
options = Preferences.@load_preference(solvername, default)

and use defaults of e.g. [true, false, false, true].

options[1] == "true" && push!(solver_list, solvertype())
options[2] == "true" && push!(solver_list, solvertype(autodiff = false))
options[3] == "true" && push!(solver_list, solvertype(chunk_size = 1))
options[4] == "true" && push!(solver_list, solvertype(chunk_size = Val{1}()))
end
for (solvername, default, nonstifftype, stifftype) in autoswitch_solver_options
options = split(Preferences.@load_preference(solvername, default), ';')
options[1] == "true" && push!(solver_list, nonstifftype(stifftype()))
options[2] == "true" &&
push!(solver_list, nonstifftype(stifftype(autodiff = false)))
options[3] == "true" && push!(solver_list, nonstifftype(stifftype(chunk_size = 1)))
options[4] == "true" &&
push!(solver_list, nonstifftype(stifftype(chunk_size = Val{1}())))
end

for prob in prob_list, solver in solver_list
solve(prob,solver)(5.0)
end
eltypes = []
Preferences.@load_preference("Float64", "true") == "true" && push!(eltypes, Float64)
Preferences.@load_preference("Float32", "false") == "true" && push!(eltypes, Float64)
Comment on lines +254 to +255
Copy link
Member

Choose a reason for hiding this comment

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

Suggested change
Preferences.@load_preference("Float64", "true") == "true" && push!(eltypes, Float64)
Preferences.@load_preference("Float32", "false") == "true" && push!(eltypes, Float64)
Preferences.@load_preference("Float64", true) === true && push!(eltypes, Float64)
Preferences.@load_preference("Float32", false) === true && push!(eltypes, Float64)

prob_list = []
inplace_problems =
split(Preferences.@load_preference("iip_problems", "true;true;true"), ';') .==
"true"
out_of_place_problems =
split(Preferences.@load_preference("oop_problems", "true;false;false"), ';') .==
"true"
for T in eltypes
if inplace_problems[1]
push!(prob_list, ODEProblem(lorenz, T[1.0; 0.0; 0.0], (zero(T), one(T))))
end
if inplace_problems[2]
push!(
prob_list,
ODEProblem{true,false}(lorenz, T[1.0; 0.0; 0.0], (zero(T), one(T))),
)
end
if inplace_problems[3]
push!(
prob_list,
ODEProblem{true,false}(lorenz, T[1.0; 0.0; 0.0], (zero(T), one(T)), T[]),
)
end
if out_of_place_problems[1]
push!(prob_list, ODEProblem(lorenz_oop, T[1.0; 0.0; 0.0], (zero(T), one(T))))
end
if out_of_place_problems[2]
push!(
prob_list,
ODEProblem{false,false}(lorenz_oop, T[1.0; 0.0; 0.0], (zero(T), one(T))),
)
end
if out_of_place_problems[3]
push!(
prob_list,
ODEProblem{false,false}(
lorenz_oop,
T[1.0; 0.0; 0.0],
(zero(T), one(T)),
T[],
),
)
end
end

prob_list = nothing
for prob in prob_list, solver in solver_list
solve(prob, solver)(5.0)
end

prob_list = nothing
end

#General Functions
export solve, solve!, init, step!

Expand Down