ModelingToolkit with DataInterpolation in time-varying external force

General Usage
,
1

Hi there, I’m trying to get to work the simple example of an ODE with the external varying force with a DataInterpolation function, but I get an error building the ODESystem. The example is as follow:

using ModelingToolkit
using DifferentialEquations: solve
using DataInterpolations: LinearInterpolation

@variables t v(t) f(t) 
@parameters C      
D = Differential(t) 

a = float.(collect(0:10))
b = 2.0 .* a .+ randn() 
f_fun(t) = LinearInterpolation(b, a)(t)

@named im = ODESystem([f ~ f_fun(t), D(v) ~ f/C])

The last line errors with:

MethodError: isless(::Float64, ::Num) is ambiguous. Candidates:
  isless(a::Real, b::Num) in Symbolics at C:\Users\flvisa\.julia\packages\Symbolics\h8kPL\src\num.jl:156
  isless(x::AbstractFloat, y::Real) in Base at operators.jl:169
Possible fix, define
  isless(::AbstractFloat, ::Num)
lt(o::Base.Order.ForwardOrdering, a::Float64, b::Num) at ordering.jl:109
searchsortedfirst at sort.jl:184 [inlined]
searchsortedfirst at sort.jl:295 [inlined]
#searchsortedfirst#4 at sort.jl:297 [inlined]
searchsortedfirst at sort.jl:297 [inlined]
_interpolate(A::LinearInterpolation{Vector{Float64}, Vector{Float64}, true, Float64}, t::Num) at interpolation_methods.jl:3
(::LinearInterpolation{Vector{Float64}, Vector{Float64}, true, Float64})(t::Num) at DataInterpolations.jl:32
f_fun(t::Num) at untitled-35d453aeb3d55920c72726264c4c44d7:23
top-level scope at untitled-35d453aeb3d55920c72726264c4c44d7:25
eval at boot.jl:360 [inlined]

I made the interpolation to work for this example.

Thanks.

2

Did you try @register f_fun(t)?

3

Oh dear, yeah, now the example does work.

I still don’t get a good solution for the piecewise f_fun(t) I have, as shown below:

profile

When I run the solver inside a for loop I get a fair solution, not good enough for the model, but closer than that from ModelingToolkit. Below the data y2 I am trying to reproduce, and in blue the solution with the example above.

results

Not sure if I am approaching correctly with the MTK though.

Thanks for your help Chris. Cheers!

4

Lower the tolerance and use tstops?

1 Like
5

The tstops did the trick.
Thanks for the help again and for the great tool you guys developed.
Cheers.