Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Is it possible to get data within @taylorize #635

Open
AnderGray opened this issue Jun 16, 2022 · 0 comments
Open

Is it possible to get data within @taylorize #635

AnderGray opened this issue Jun 16, 2022 · 0 comments

Comments

@AnderGray
Copy link

I have a dynamics problem where I need to fetch data (e.g. using linear interpolation) within the derivative. Is such a thing possible?

using ReachabilityAnalysis, Interpolations

tspan = (0.0, 10.0)

f(x) = x^2    # real f
xs = range(tspan[1], tspan[2], length = 20);

data = f.(xs)
Lin_interp = LinearInterpolation(data, XS)  # interpolated f

@taylorize function f_real(du, u, p, t)
    du[1] = f(t)
end

@taylorize function f_interp(du, u, p, t)
    du[1] = Lin_interp(xs)
end

X0 = Interval(0, 0.1)

prob = @ivp(x' = f_real(x), x(0)  X0, dim=1)
prob_interp = @ivp(x' = f_interp(x), x(0)  X0, dim=1)

sol = solve(prob, tspan = (0.0, 10.0), alg = TMJets21a(abstol = 1e-10))
sol_interp = solve(prob_interp, tspan = (0.0, 10.0), alg = TMJets21a(abstol = 1e-10))

plot(sol,vars = (0, 1))
plot(sol_interp, vars = (0, 1))

When running solve with the interpolated version, it tries to evaluate the function with a Taylor series.

My initial thoughts on doing this is that we could take the bounding box of the Taylor model, and evaluate the interpolator with it. Perhaps using interval arithmetic, but for the case of linear interpolations I believe we can get the exact bounds, as either the endpoints or the maximum and minimum data value within the range. i.e. something like this:

using IntervalArithmetic: interval
function interpolate(Lin_interp, x::Interval)

    grid = Lin_interp.itp.knots[1]
    data = Lin_interp.itp.coefs

    data_in_x = data[grid .∈ x]
    x_lo = Lin_interp(x.lo)
    x_hi = Lin_interp(x.hi)

    points = [data_in_x; x_lo; x_hi]

    return interval(minimum(points), maximum(points))

end
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant