Hi everyone,

I would like to use Symplectic Integrators such as `VelocityVerlet`

or `VerletLeapfrog`

from DifferentialEquations.jl to ensure the time-reversibility of the equation of motion.

Nonetheless, my problem can not be set as a `DynamicalODEProblem`

.

My equation of motion is just:

`\frac{du}{dt} = f2(u)`

where v = f2(u) is a model for the velocity which depends only of the position u )

Since my function f2 is not a function of the velocity v but only of the position u, I can not use these symplectic integrators.

Is there a way to get around this issue?

Thank you for your time,

where did you get this idea? I think you flipped the equations. The docs say

http://docs.juliadiffeq.org/latest/types/dynamical_types.html#Mathematical-Specification-of-a-Dynamical-ODE-Problem-1

This is a Partitioned ODE partitioned into two groups, so the functions should be specified as `f1(dv,v,u,p,t)`

and `f2(du,v,u,p,t)`

(in the inplace form), where `f1`

is independent of `v`

(unless specified by the solver), and `f2`

is independent of `u`

and `t`

.

which matches your problem. Position depends on velocity, velocity depends on position and possibly time.

Thank you for your answer,

Sorry I should have been clearer, I mean that I donâ€™t have 2 ODEs to solve (velocity v and position u ) but only one ODE to solve for the position u. Indeed, the velocity is given by a kinematic model and not by an ODE, i.e. I know in advance the velocity v at each position u : `v=f2(u)`

,

Therefore I donâ€™t know what equation to set for \frac{dv}{dt}

I want to do something like this:

```
function f1(dv,v,u, p ,t)
dv[1:2] = [0.0, 0.0]
end
function f2(du,v,u,p,t)
#Velocity(u) known model for velocity
du[1:2] = Velocity(u)
end
```

Maybe I donâ€™t need to use `DynamicalODEProblem`

, I just want to ensure the time-reversibility of my solution.

Oh, then itâ€™s not a good problem for these integrators. Youâ€™ll want to use a reversible method for first order ODEs, like `Trapezoid()`

. You can set that to fixed point iteration if the problem isnâ€™t stiff.

Thanks for the advice,

My model is an interpolated function using `CubicSplineInterpolation`

of Interpolations.jl, therefore it doesnâ€™t support ForwardDiff used by default in the implicit time scheme `Trapezoid`

.

Do you know a derivative-free method I can use in `Trapezoid()`

?

Else I am thinking that I can probably provide the Jacobian of the model but I think we can only get Gradients with Interpolations.jl and not Jacobians.

Yes, try `Trapezoid(nlsolve = NLFunctional())`

or `Trapezoid(nlsolve = NLAnderson())`

Thank you @ChrisRackauckas for your help!

1 Like