Skip to content

Commit

Permalink
test a non-parameterized function
Browse files Browse the repository at this point in the history
  • Loading branch information
ChrisRackauckas committed Mar 22, 2020
1 parent 6478139 commit 2307296
Show file tree
Hide file tree
Showing 5 changed files with 41 additions and 28 deletions.
1 change: 1 addition & 0 deletions .gitignore
Original file line number Diff line number Diff line change
@@ -1,3 +1,4 @@
*.jl.cov
*.jl.*.cov
*.jl.mem
Manifest.toml
4 changes: 3 additions & 1 deletion Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -5,12 +5,14 @@ version = "0.1.0"
[deps]
DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e"
MATLAB = "10e44e05-a98a-55b3-a45b-ba969058deb6"
ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78"
Reexport = "189a3867-3050-52da-a836-e630ba90ab69"

[compat]
Reexport = "0.2"
DiffEqBase = "6.5"
MATLAB = "0.7"
ModelingToolkit = "1.4.2"
Reexport = "0.2"
julia = "1"

[extras]
Expand Down
33 changes: 21 additions & 12 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -23,20 +23,19 @@ Pkg.add("MATLABDiffEq")

## Using MATLABDiffEq.jl

MATLABDiffEq.jl is simply a solver on the DiffEq common interface, so for details see the [DifferentialEquations.jl documentation](https://juliadiffeq.github.io/DiffEqDocs.jl/dev/). However, there are three things to know:
MATLABDiffEq.jl is simply a solver on the DiffEq common interface, so for details see the [DifferentialEquations.jl documentation](https://juliadiffeq.github.io/DiffEqDocs.jl/dev/).
However, the only options implemented are those for error calculations
(`timeseries_error`), `saveat` and tolerances.

1. The only options implemented are those for error calculations (`timeseries_error`), `saveat` and tolerances.
2. The input function must be defined by a `ParameterizedFunction`
3. The input function must not use parameters

Note that the algorithms are defined to have the same name as the MATLAB algorithms, but are not exported. Thus to use `ode45`, you would specify the algorithm as `MATLABDiffEq.ode45()`.
Note that the algorithms are defined to have the same name as the MATLAB algorithms,
but are not exported. Thus to use `ode45`, you would specify the algorithm as
`MATLABDiffEq.ode45()`.

## Example

```julia
using MATLABDiffEq, ParameterizedFunctions


f = @ode_def LotkaVolterra begin
dx = 1.5x - x*y
dy = -3y + x*y
Expand All @@ -46,6 +45,16 @@ tspan = (0.0,10.0)
u0 = [1.0,1.0]
prob = ODEProblem(f,u0,tspan)
sol = solve(prob,MATLABDiffEq.ode45())

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
u0 = [1.0;0.0;0.0]
tspan = (0.0,100.0)
prob = ODEProblem(lorenz,u0,tspan)
sol = solve(prob,MATLABDiffEq.ode45())
```

## Measuring Overhead
Expand All @@ -56,7 +65,7 @@ is done. Thus you can simply call the same ODE function and time it
directly. This is done by:

```julia
@time MATLABDiffEq.eval_string("[t,u] = $(algstr)(f,tspan,u0,options);")
@time MATLABDiffEq.eval_string("[t,u] = $(algstr)(diffeqf,tspan,u0,options);")
```

To be even more pedantic, you can play around in the actual MATLAB
Expand Down Expand Up @@ -97,16 +106,16 @@ julia> @time sol = solve(prob,alg);
julia> @time sol = solve(prob,alg);
0.065460 seconds (38.84 k allocations: 1.556 MB)

julia> @time MATLABDiffEq.eval_string("[t,u] = $(algstr)(f,tspan,u0,options);")
julia> @time MATLABDiffEq.eval_string("[t,u] = $(algstr)(diffeqf,tspan,u0,options);")
0.058249 seconds (11 allocations: 528 bytes)

julia> @time MATLABDiffEq.eval_string("[t,u] = $(algstr)(f,tspan,u0,options);")
julia> @time MATLABDiffEq.eval_string("[t,u] = $(algstr)(diffeqf,tspan,u0,options);")
0.060367 seconds (11 allocations: 528 bytes)

julia> @time MATLABDiffEq.eval_string("[t,u] = $(algstr)(f,tspan,u0,options);")
julia> @time MATLABDiffEq.eval_string("[t,u] = $(algstr)(diffeqf,tspan,u0,options);")
0.060171 seconds (11 allocations: 528 bytes)

julia> @time MATLABDiffEq.eval_string("[t,u] = $(algstr)(f,tspan,u0,options);")
julia> @time MATLABDiffEq.eval_string("[t,u] = $(algstr)(diffeqf,tspan,u0,options);")
0.058928 seconds (11 allocations: 528 bytes)
```

Expand Down
22 changes: 7 additions & 15 deletions src/MATLABDiffEq.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@ module MATLABDiffEq

using Reexport
@reexport using DiffEqBase
using MATLAB
using MATLAB, ModelingToolkit

abstract type MATLABAlgorithm <: DiffEqBase.AbstractODEAlgorithm end
struct ode23 <: MATLABAlgorithm end
Expand All @@ -22,10 +22,6 @@ function DiffEqBase.__solve(

tType = eltype(tupType)

if !(typeof(prob.f) <: DiffEqBase.AbstractParameterizedFunction)
error("Functions must be defined via ParameterizedFunctions.jl to work with this package.")
end

if prob.tspan[end]-prob.tspan[1]<tType(0)
error("final time must be greater than starting time. Aborting.")
end
Expand All @@ -43,15 +39,11 @@ function DiffEqBase.__solve(
u0 = prob.u0
end

strs = [string(f.fex.args[2i].args[2]) for i in 1:length(f.syms)]
matstr = ""
for i in 1:length(strs)
matstr *= strs[i]
i < length(strs) && (matstr *= "; ")
end
matstr = replace(matstr,"["=>"(")
matstr = replace(matstr,"]"=>")")
matstr = "f = @(t,internal_var___u) ["*matstr*"];"
sys = first(modelingtoolkitize(prob))

matstr = ModelingToolkit.build_function(sys.eqs,sys.dvs,
sys.ps,sys.iv,
target = ModelingToolkit.MATLABTarget())

# Send the variables
put_variable(get_default_msession(),:tspan,tspan)
Expand All @@ -66,7 +58,7 @@ function DiffEqBase.__solve(
eval_string("options = odeset('RelTol',reltol,'AbsTol',abstol);")
algstr = string(typeof(alg).name.name)
#algstr = replace(string(typeof(alg)),"MATLABDiffEq.","")
eval_string("[t,u] = $(algstr)(f,tspan,u0,options);")
eval_string("[t,u] = $(algstr)(diffeqf,tspan,u0,options);")
ts = jvector(get_mvariable(:t))
timeseries_tmp = jarray(get_mvariable(:u))

Expand Down
9 changes: 9 additions & 0 deletions test/runtests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,16 @@ p = [1.5,1,3,1]
tspan = (0.0,10.0)
u0 = [1.0,1.0]
prob = ODEProblem(f,u0,tspan,p)
sol = solve(prob,MATLABDiffEq.ode45())

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
u0 = [1.0;0.0;0.0]
tspan = (0.0,100.0)
prob = ODEProblem(lorenz,u0,tspan)
sol = solve(prob,MATLABDiffEq.ode45())

algs = [MATLABDiffEq.ode23
Expand Down

0 comments on commit 2307296

Please sign in to comment.