# Awkward ODE problem

**URL:** <https://discourse.julialang.org/t/awkward-ode-problem/88435>\
**Category:** Modelling & Simulations\
**Tags:** question, differentialequation\
**Created:** [October 8, 2022, 10:27am UTC](https://discourse.julialang.org/t/awkward-ode-problem/88435 "2022-10-08T10:27:28Z")\
**Posts on this page:** 8\
**Page:** 1

<div class="post-metadata">

**Author:** ![josuagrw](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/josuagrw/32/1015_2.png) [@josuagrw](https://discourse.julialang.org/u/josuagrw)\
**Post date:** [October 8, 2022, 10:27am UTC](https://discourse.julialang.org/t/awkward-ode-problem/88435/1 "2022-10-08T10:27:28Z")

</div>

I have a differential equation of the kind

\dot q\_1 = f(q\_1, q\_2)\\ A(q\_1,q\_2)\dot q\_2=B(q\_1,q\_2)\dot q\_1

where q\_1, q\_2 are vectors, A is a symmetric, positive-definite real matrix, and B is a real (non-square) matrix.

When I use `Tsit5` or `RK4`, the time step `dt` fluctuates a lot and sometimes the solver goes into a loop (evaluating the r.h.s. over and over without moving forward in time). The problem goes away when I set \dot q\_2=0.

So, I suspect that the solution to the second equation (the linear equation) is not good enough and thus the error estimate for adaptive methods is unreliable.

Any ideas what could be a better approach here? Am I using the wrong solvers?

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [October 8, 2022, 10:28am UTC](https://discourse.julialang.org/t/awkward-ode-problem/88435/2 "2022-10-08T10:28:54Z")

</div>

Is this setup using mass matrices? That wouldn’t work with explicit methods so I’m a bit confused there.

How do the implicit or Rosenbrock methods do?

---

<div class="post-metadata">

**Author:** ![josuagrw](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/josuagrw/32/1015_2.png) [@josuagrw](https://discourse.julialang.org/u/josuagrw)\
**Post date:** [October 8, 2022, 10:29am UTC](https://discourse.julialang.org/t/awkward-ode-problem/88435/3 "2022-10-08T10:29:49Z")

</div>

In practice I am inverting A (it should be invertible). I am not using the mass matrix solvers in OrdinaryDiffEq (yet)

EDIT: I should be more precise. I run the problem as

\dot q\_2 = A^{-1}B\dot q\_1

using `ldiv` from `LinearAlgebra`

---

<div class="post-metadata">

**Author:** ![josuagrw](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/josuagrw/32/1015_2.png) [@josuagrw](https://discourse.julialang.org/u/josuagrw)\
**Post date:** [October 8, 2022, 10:34am UTC](https://discourse.julialang.org/t/awkward-ode-problem/88435/4 "2022-10-08T10:34:15Z")

</div>

The vectors q\_1 and q\_2 have combined 4000 components. I’m not sure implicit / Rosenbrock methods are feasible here. I have no idea how I would calculate the Jacobian in a practical amount of time.

---

<div class="post-metadata">

**Author:** ![Oscar\_Smith](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/oscar_smith/32/25343_2.png) [@Oscar\_Smith](https://discourse.julialang.org/u/Oscar_Smith)\
**Post date:** [October 8, 2022, 3:26pm UTC](https://discourse.julialang.org/t/awkward-ode-problem/88435/5 "2022-10-08T15:26:31Z")

</div>

try FBDF. it’s often really good.

---

<div class="post-metadata">

**Author:** ![josuagrw](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/josuagrw/32/1015_2.png) [@josuagrw](https://discourse.julialang.org/u/josuagrw)\
**Post date:** [October 8, 2022, 3:41pm UTC](https://discourse.julialang.org/t/awkward-ode-problem/88435/6 "2022-10-08T15:41:53Z")

</div>

Thanks, I tried `ImplicitEuler` with Krylov GMRES and `autodiff=false` which works (albeit very slowly).

I’ll try FBDF next

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [October 8, 2022, 3:46pm UTC](https://discourse.julialang.org/t/awkward-ode-problem/88435/7 "2022-10-08T15:46:10Z")

</div>

> [@josuagrw](#):
>
> Thanks, I tried `ImplicitEuler` with Krylov GMRES and `autodiff=false` which works (albeit very slowly).

`ImplicitEuler` is essentially always slow. That’s kind of by definition though. `FBDF`, `QNDF`, etc. try the efficient methods if you need efficiency.

---

<div class="post-metadata">

**Author:** ![josuagrw](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/josuagrw/32/1015_2.png) [@josuagrw](https://discourse.julialang.org/u/josuagrw)\
**Post date:** [October 24, 2022, 11:14am UTC](https://discourse.julialang.org/t/awkward-ode-problem/88435/8 "2022-10-24T11:14:02Z")

</div>

So, just to close this thread out:

I changed the definition of the problem so that the integrand fluctuates less and therefore appears less stiff to the solvers. Now it runs nicely with `RK4()` (and probably `Tsit5()` but each round of testing takes a long time so I won’t let perfect be the enemy of good and leave it at that).