# ODE solver save intermediate outputs (y = C\*x + D\*u)

**URL:** <https://discourse.julialang.org/t/ode-solver-save-intermediate-outputs-y-c-x-d-u/41697>\
**Category:** General Usage\
**Tags:** diffeq\
**Created:** [June 18, 2020, 8:46pm UTC](https://discourse.julialang.org/t/ode-solver-save-intermediate-outputs-y-c-x-d-u/41697 "2020-06-18T20:46:20Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![mopg](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mopg/32/3368_2.png) [@mopg](https://discourse.julialang.org/u/mopg)\
**Post date:** [June 18, 2020, 8:46pm UTC](https://discourse.julialang.org/t/ode-solver-save-intermediate-outputs-y-c-x-d-u/41697/1 "2020-06-18T20:46:20Z")

</div>

I am using DifferentialEquations.jl for solving an ODE, and it works great. However, I’d like to save intermediate outputs, and I can’t figure out how to do that.

Consider the example below:

```
using DifferentialEquations
α=1
u0=1/2
function f(t,not_used,u)

    p = u^2 # how to save p??
    
    return α*u
    
end
tspan = (0.0,1.0)
prob = ODEProblem(f,u0,tspan)

solve(prob)

```

How can I save `p` as a function of `t`?

This seems to be a pretty standard application of solving ODEs, consider the system below for instance

```julia
x_dot(t) = A * x(t) + B * u(t)
y(t) = C * x(t) + D * u(t)

```

How does one save `y(t)`?

---

<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:** [June 18, 2020, 11:10pm UTC](https://discourse.julialang.org/t/ode-solver-save-intermediate-outputs-y-c-x-d-u/41697/2 "2020-06-18T23:10:19Z")

</div>

Just calculate `y(t)` from the output. `solve` outputs a continuous function.

---

<div class="post-metadata">

**Author:** ![mopg](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mopg/32/3368_2.png) [@mopg](https://discourse.julialang.org/u/mopg)\
**Post date:** [June 18, 2020, 11:50pm UTC](https://discourse.julialang.org/t/ode-solver-save-intermediate-outputs-y-c-x-d-u/41697/3 "2020-06-18T23:50:10Z")

</div>

Sure, I can do that. However, then I have to recompute all intermediate variables a second time. My application is fairly computationally intensive, so I’d very much prefer not to compute these intermediate variables twice.

---

<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:** [June 18, 2020, 11:57pm UTC](https://discourse.julialang.org/t/ode-solver-save-intermediate-outputs-y-c-x-d-u/41697/4 "2020-06-18T23:57:57Z")

</div>

It’s at most 1/8 of the cost, if y(t) \>\> x(t) in terms of size. That’s the reason why it’s not implemented in any ODE solver. I can show some ways to do it, but if you’re looking at this for efficiency you’re looking in the wrong place.

---

<div class="post-metadata">

**Author:** ![mopg](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mopg/32/3368_2.png) [@mopg](https://discourse.julialang.org/u/mopg)\
**Post date:** [June 19, 2020, 12:09am UTC](https://discourse.julialang.org/t/ode-solver-save-intermediate-outputs-y-c-x-d-u/41697/5 "2020-06-19T00:09:42Z")

</div>

For the systems I’m interested in, typically size(y) is similar to size(x) – I’m not trying to save every single intermediate variable. Does your conclusion still apply then (or is the cost even more negligible)?

Out of curiosity, where does the 1/8 factor come from?

---

<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:** [June 19, 2020, 12:38am UTC](https://discourse.julialang.org/t/ode-solver-save-intermediate-outputs-y-c-x-d-u/41697/6 "2020-06-19T00:38:04Z")

</div>

5th order methods like `Tsit5()` do 7 `f` evaluations per step, so if you do an entire extra one to compute the `y` part once more, then it’s at most `7+1` evaluations where you’re adding 1, so 1/8. If you’re only having to do half of that work, and if you’re outputting less points then the steps, and …, then that factor keeps dropping.