-
-
Notifications
You must be signed in to change notification settings - Fork 78
/
Sundials.jl
131 lines (105 loc) · 3.43 KB
/
Sundials.jl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
# __precompile__()
module Sundials
import Reexport
Reexport.@reexport using DiffEqBase
using SciMLBase: AbstractSciMLOperator
import DataStructures
import Logging
import DiffEqBase
import SparseArrays
import LinearAlgebra
import Libdl
using CEnum
const warnkeywords = (:isoutofdomain,
:unstable_check,
:calck,
:internalnorm,
:gamma,
:beta1,
:beta2,
:qmax,
:qmin,
:qoldinit)
warnlist = Set(warnkeywords)
warnida = union(warnlist, Set((:dtmin,)))
using Sundials_jll
export solve,
SundialsODEAlgorithm, SundialsDAEAlgorithm, ARKODE, CVODE_BDF, CVODE_Adams, IDA,
KINSOL
# some definitions from the system C headers wrapped into the types_and_consts.jl
const DBL_MAX = prevfloat(Inf)
const DBL_MIN = nextfloat(-Inf)
const DBL_EPSILON = eps(Cdouble)
const int64_t = Int64
const MPI_DOUBLE = nothing
const MPI_INT64_T = nothing
macro SUNDIALS_F77_FUNC(name, NAME)
Symbol(string(name) * "_64_")
end
# Issue 328: these symbols exist in other libraries
const libsundials_sundials = libsundials_cvode
const libsundials_sunlinsol = libsundials_cvode
const libsundials_sunnonlinsol = libsundials_cvode
# sunmatrix has been renamed to sunmatrix[dense/sparse/band]
# const libsundials_sunmatrix = libsundials_cvode
const SPARSE_SOLVERS = (:KLU,)
include("../lib/libsundials_common.jl")
include("types_and_consts_additions.jl")
include("handle.jl")
include("nvector_wrapper.jl")
include("../lib/libsundials_api.jl")
for ff in names(@__MODULE__; all = true)
fname = string(ff)
if occursin("SetLinearSolver", fname) &&
!occursin("#", fname) && # filter out compiler generated names
!occursin("Dls", fname) && !occursin("Spils", fname) # filter out old names
@eval function $ff(mem, LS::SUNLinearSolver, A::Ptr, args...)
$ff(mem, LS,
convert(SUNMatrix, A),
args...)
end
end
end
include("simple.jl")
include("common_interface/function_types.jl")
include("common_interface/verbosity.jl")
include("common_interface/algorithms.jl")
include("common_interface/integrator_types.jl")
include("common_interface/integrator_utils.jl")
include("common_interface/solve.jl")
import PrecompileTools
PrecompileTools.@compile_workload 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 = [
ARKODE(), CVODE_Adams(), CVODE_BDF(),
]
prob_list = [ODEProblem(lorenz, [1.0; 0.0; 0.0], (0.0, 1.0)),
ODEProblem(lorenz, [1.0; 0.0; 0.0], (0.0, 1.0), Float64[])]
for prob in prob_list, solver in solver_list
solve(prob, solver)(0.5)
end
prob_list = nothing
function f(out, du, u, p, t)
out[1] = -0.04u[1] + 1e4 * u[2] * u[3] - du[1]
out[2] = +0.04u[1] - 3e7 * u[2]^2 - 1e4 * u[2] * u[3] - du[2]
out[3] = u[1] + u[2] + u[3] - 1.0
end
u₀ = [1.0, 0, 0]
du₀ = [-0.04, 0.04, 0.0]
tspan = (0.0, 100000.0)
differential_vars = [true, true, false]
prob = DAEProblem(f, du₀, u₀, tspan, differential_vars = differential_vars)
sol = solve(prob, IDA())
prob = nothing
end
##################################################################
# Deprecations
##################################################################
end # module