# Integrating Hamilton's equations
**URL:** https://discourse.julialang.org/t/integrating-hamiltons-equations/96147
**Category:** Numerics
**Tags:** question, differentialequation
**Created:** [March 15, 2023, 7:52pm UTC](https://discourse.julialang.org/t/integrating-hamiltons-equations/96147 "2023-03-15T19:52:37Z")
**Posts on this page:** 5
**Page:** 1
### Author:  [@roi.holtzman](https://discourse.julialang.org/u/roi.holtzman)
#### Post date: [March 15, 2023, 7:52pm UTC](https://discourse.julialang.org/t/integrating-hamiltons-equations/96147/1 "2023-03-15T19:52:37Z")
I want to integrate the equations of motion of a classical physical system.
My Hamiltonian is changing in time, and actually, I am interested in considering the general case where it is **not** necessarily of the form H(q, p) = \frac{p^2}{2m} + U(q).
My goal is to get high accuracy in the solution, so I guess that using a symplectic integrator is the way-to-go. (Or maybe not? Note that since the Hamiltonian is changing in time, energy is not conserved, but still the symplectic structure holds).
From the documentation of [DifferentialEquations.jl](https://docs.sciml.ai/DiffEqDocs/stable/examples/classical_physics/#Second-Order-Linear-ODE) I could not understand what is a good approach.
It seems to me like the only options to use symplectic integrators in my case (time changing Hamiltonian, and non-standard form of q, p) are `SecondOrderODEProblem` or `DynamicalODEFunction`. However, I could not find an example that showcases `DynamicalODEFunction`. Is there such an example?
I tried to use `SecondOrderODEProblem` from the [example here](https://docs.sciml.ai/DiffEqDocs/stable/examples/classical_physics/#Symplectic-Integration), and to modify from the two coordinates x,y to just one coordinate x. However, I get that the velocity does not change in time. I am surely doing something wrong. I would love to get some help.
Here is my modification of the example:
```julia
function HH_acceleration!(dv,v,u,p,t)
x = u
dx = dv
dv = -cos.(x)
end
initial_positions = [0.1]
initial_velocities = [0.5]
prob = SecondOrderODEProblem(HH_acceleration!,initial_velocities,initial_positions,(0.0, 10.0))
sol2 = solve(prob, KahanLi8(), dt=1/10)
```
And what I get back is:
```julia
retcode: Success
Interpolation: 3rd order Hermite
t: 101-element Vector{Float64}:
0.0
0.1
0.2
0.30000000000000004
0.4
0.5
0.6
0.7
0.7999999999999999
0.8999999999999999
0.9999999999999999
1.0999999999999999
1.2
⋮
8.899999999999984
8.999999999999984
9.099999999999984
9.199999999999983
9.299999999999983
9.399999999999983
9.499999999999982
9.599999999999982
9.699999999999982
9.799999999999981
9.89999999999998
10.0
u: 101-element Vector{ArrayPartition{Float64, Tuple{Vector{Float64}, Vector{Float64}}}}:
([0.5], [0.1])
([0.5], [0.15])
([0.5], [0.20000000000000004])
([0.5], [0.25000000000000006])
([0.5], [0.3000000000000001])
([0.5], [0.35000000000000014])
([0.5], [0.4000000000000002])
([0.5], [0.45000000000000023])
([0.5], [0.5000000000000002])
([0.5], [0.5500000000000003])
([0.5], [0.6000000000000003])
([0.5], [0.6500000000000004])
([0.5], [0.7000000000000004])
⋮
([0.5], [4.550000000000033])
([0.5], [4.600000000000033])
([0.5], [4.650000000000034])
([0.5], [4.700000000000035])
([0.5], [4.7500000000000355])
([0.5], [4.800000000000036])
([0.5], [4.850000000000037])
([0.5], [4.900000000000038])
([0.5], [4.950000000000038])
([0.5], [5.000000000000039])
([0.5], [5.05000000000004])
([0.5], [5.1000000000000405])
```
---
### Author:  [@SteffenPL](https://discourse.julialang.org/u/SteffenPL)
#### Post date: [March 16, 2023, 6:48am UTC](https://discourse.julialang.org/t/integrating-hamiltons-equations/96147/2 "2023-03-16T06:48:13Z")
You need to write
```julia
dx .= dv
dv .= -cos.(x)
```
to overwrite the input vectors instead of creating new ones.
---
### Author:  [@empet](https://discourse.julialang.org/u/empet)
#### Post date: [March 16, 2023, 7:01am UTC](https://discourse.julialang.org/t/integrating-hamiltons-equations/96147/3 "2023-03-16T07:01:01Z")
The second order ODE equation is defined as follows:
```julia
function eqn!(ddu, du, u, p, t)
ddu[1]=-cos(u[1])
end
function getsol()
#du0=[0.5]
#u0=[0.1]
#tspan=(0, 10.0)
probl=SecondOrderODEProblem(eqn!, [0.5], [0.1], (0.0,10.0))
solve(probl, KahanLi8(), dt=0.05)
end
s=getsol()
```
---
### Author:  [@roi.holtzman](https://discourse.julialang.org/u/roi.holtzman)
#### Post date: [March 16, 2023, 1:17pm UTC](https://discourse.julialang.org/t/integrating-hamiltons-equations/96147/4 "2023-03-16T13:17:03Z")
> [@SteffenPL](#):
>
> to overwrite the input vectors instead of creating new ones.
Thanks! Actually, this is what is working for me
```julia
function HH_acceleration!(dv,v,u,p,t)
x = u
dx = dv
dv .= -cos.(x)
end
```
---
### Author:  [@SteffenPL](https://discourse.julialang.org/u/SteffenPL)
#### Post date: [March 16, 2023, 1:24pm UTC](https://discourse.julialang.org/t/integrating-hamiltons-equations/96147/5 "2023-03-16T13:24:13Z")
Oh, yes, you are right! You can actually remove `dx = dv` since it doesn’t do anything (for second order problems) 😉