# How to efficiently calculate the Jacobian of the solution of a differential equation with respect to the initial conditions?

**URL:** <https://discourse.julialang.org/t/how-to-efficiently-calculate-the-jacobian-of-the-solution-of-a-differential-equation-with-respect-to-the-initial-conditions/76650>\
**Category:** Performance\
**Tags:** diffeq, forwarddiff\
**Created:** [February 17, 2022, 8:02pm UTC](https://discourse.julialang.org/t/how-to-efficiently-calculate-the-jacobian-of-the-solution-of-a-differential-equation-with-respect-to-the-initial-conditions/76650 "2022-02-17T20:02:15Z")\
**Posts on this page:** 2\
**Page:** 1

<div class="post-metadata">

**Author:** ![marcsgil](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/marcsgil/32/33908_2.png) [@marcsgil](https://discourse.julialang.org/u/marcsgil)\
**Post date:** [February 17, 2022, 8:02pm UTC](https://discourse.julialang.org/t/how-to-efficiently-calculate-the-jacobian-of-the-solution-of-a-differential-equation-with-respect-to-the-initial-conditions/76650/1 "2022-02-17T20:02:16Z")

</div>

I need to numerically perform an integral of the form

\int\_{0}^{t\_{max}} F \left( \mathbf{u}(t), \frac{\partial \mathbf{u}(t)}{ \partial \mathbf{u}(0)} \right) dt

where F is a certain function, \mathbf{u}(t) is the solution of a simple ODE (it is a 1D Hamiltonian system) and \frac{\partial \mathbf{u}(t)}{ \partial \mathbf{u}(0)} is the Jacobian of \mathbf{u}(t) with respect to the initial condition \mathbf{u}(0). For that reason, I need a fast way of evaluating \mathbf{u}(t) and \frac{\partial \mathbf{u}(t)}{ \partial \mathbf{u}(0)} for many values of t.

The following code allows me to calculate these values, but it isn’t efficient:

```julia
using DifferentialEquations, ForwardDiff

function f!(du,u,p,t)
    du[1] = exp(-u[2])*(1-exp(-u[2]))
    du[2] = 2*u[1]
end

function get_trajectory(u0,t_max)
    problem = ODEProblem(f!,u0,(0.0,t_max))
    solve(problem)
end

function get_jacobian(u0,t_max)
    problem(u) = ODEProblem(f!,u,(0.0,t_max))
    t->ForwardDiff.jacobian(u->solve(problem(u))(t),u0)
end

u = get_trajectory([1.0,0.0], 1.0)
u(0.5)

jac = get_jacobian([1.0,0.0], 1.0)
jac(0.5)

```

It works fine for calculating \mathbf{u}, as, when I call `get_trajectory`, the differential equation is solved once and then each call of `u(t)` is fast. On the other hand, each time I call `jac(t)` the differential equation is solved once again, which makes the call slow. I was wondering if there is a way to solve the differential equation only once and get all the values I need.

On the [DifferentialEquations.jl documentation](https://diffeq.sciml.ai/stable/analysis/sensitivity/) there is a page on Sensitivity Analysis that appears to treat a similar question, but the Jacobian is calculated with respect to parameters. I also found [this topic](https://discourse.julialang.org/t/sensitivities-with-respect-to-initial-conditions-in-differentialequations-jl/25555) which may have a solution to my problem, but the links containing the references to the solution appear broke to me.

---

<div class="post-metadata">

**Author:** ![goerch](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/goerch/32/29122_2.png) [@goerch](https://discourse.julialang.org/u/goerch)\
**Post date:** [February 21, 2022, 9:44pm UTC](https://discourse.julialang.org/t/how-to-efficiently-calculate-the-jacobian-of-the-solution-of-a-differential-equation-with-respect-to-the-initial-conditions/76650/2 "2022-02-21T21:44:24Z")

</div>

> [@marcsgil](#):
>
> On the [DifferentialEquations.jl documentation](https://diffeq.sciml.ai/stable/analysis/sensitivity/) there is a page on Sensitivity Analysis that appears to treat a similar question, but the Jacobian is calculated with respect to parameters. I also found [this topic](https://discourse.julialang.org/t/sensitivities-with-respect-to-initial-conditions-in-differentialequations-jl/25555) which may have a solution to my problem, but the links containing the references to the solution appear broke to me.

Hm, the problem of URL stability. Here is what I [found](https://diffeq.sciml.ai/stable/analysis/sensitivity/).