Solving differential equations with derivatives on both sides of the equation

Specific Domains Modelling & Simulations
,
1

I have two coupled, first order, non-linear ODEs that I’m trying to solve. However, both of the two equations contain additional derivatives on the right-hand sides of the equation (and there don’t seem to be any obvious algebraic means to rectify that). DifferentialEquations.jl requires that the equations be in the form,

\frac{du}{dt} = f(u, p, t) ,

so I’m not sure I can use DifferenetialEquations.jl to solve my ODE as it is presently formulated. I found one example problem that has derivatives on the right side: the three-body problem. Copying from the DiffEq docs, the three body problem is formulated as follows:

\frac{dy_1}{dt} = y_1 + 2\frac{dy_2}{dt} - \bar{μ}\frac{y_1+μ}{D_1} - μ\frac{y_1-\bar{μ}}{D_2}
\frac{dy_2}{dt} = y_2 - 2\frac{dy_1}{dt} - \bar{μ}\frac{y_2}{D_1} - μ\frac{y_2}{D_2}

However, the implementation that I see on Github for the right-hand side function f(u, p, t) doesn’t make much sense to me. Here is the implementation on Github:

const threebody_μ = big(0.012277471); const threebody_μ′ = 1 - threebody_μ

threebody = (du,u,p,t) -> begin
  # 1 = y₁
  # 2 = y₂
  # 3 = y₁'
  # 4 = y₂'
  D₁ = ((u[1]+threebody_μ)^2 + u[2]^2)^(3/2)
  D₂ = ((u[1]-threebody_μ′)^2 + u[2]^2)^(3/2)
  du[1] = u[3]
  du[2] = u[4]
  du[3] = u[1] + 2u[4] - threebody_μ′*(u[1]+threebody_μ)/D₁ - threebody_μ*(u[1]-threebody_μ′)/D₂
  du[4] = u[2] - 2u[3] - threebody_μ′*u[2]/D₁ - threebody_μ*u[2]/D₂
end

It appears to me that the lines defining du[3] and du[4] are actually defining the second derivatives of y_1 and y_2, rather than the first derivatives that appear in the ODE.

So, I have a few questions:

  • Is this implementation correct?
  • If it is correct, can anyone explain to me why it is correct and how it works?
  • If it’s incorrect, is there a correct way that one can do this sort of thing with DifferentialEquations.jl?

Thanks in advance.

2

That looks like a typo to me: the equations should probably have second derivatives on the left-hand sides. (It’s basically “acceleration = Force/mass”, right?). That is, the code is correct but the equations are wrong.

Can you write down the actual equations you are solving?

1 Like
3

In the most general case, this is a differential-algebraic equation (DAE), for which DifferentialEquations.jl provides several solvers.

(For example, the incorrect “three-body equations” above correspond to a mass-matrix DAE.)

2 Likes
4

Thanks! I’ll take a look at the DAE stuff.

Here’s the equations that I have:

x' = \frac{1}{\sin \theta} \left[ \frac{1}{r} ( c_x(t) - x) \sqrt{(x')^2 + (y')^2} + (\cos \theta) y' \right]

y' = \frac{1}{\sin \theta} \left[ \frac{1}{r} ( c_y(t) - y) \sqrt{(x')^2 + (y')^2} + (\cos \theta) x' \right]

(where the derivatives are with respect to t).

5

Yes, looks at first glance like an implicit DAE.

1 Like