# Improving DAE performance

**URL:** <https://discourse.julialang.org/t/improving-dae-performance/124775>\
**Category:** Modelling & Simulations\
**Created:** [January 14, 2025, 8:44pm UTC](https://discourse.julialang.org/t/improving-dae-performance/124775 "2025-01-14T20:44:27Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![clm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/clm/32/45923_2.png) [@clm](https://discourse.julialang.org/u/clm)\
**Post date:** [January 14, 2025, 8:44pm UTC](https://discourse.julialang.org/t/improving-dae-performance/124775/1 "2025-01-14T20:44:27Z")

</div>

I’m looking to solve two implicitly-defined nonlinear ODEs. The system diverges toward the final integration z=H, but in a way that I expect, T\sim1/(H-z) and u\sim\log(H-z), plus constants.

What I’m hoping to do is integrate the system to something close to 1 and subtract off this diverging behaviour from the solution to find the non-singular constants that remain.

The first equation depends only on dT/dz and z, while the second equation defines du/dz as a function of dTdz and z.

This is a MWE:

```julia
using DifferentialEquations, Sundials, ComponentArrays

function ddz!(res,dy,y,p,z)
    Q = p.Q
    T_z,U_z = dy
    
    res[1] = Q*abs(Q) + z^2*(p.H-z) * T_z^2 * abs(1 - Q^2/((p.H-z)^2*T_z))^p.α
    res[2] = U_z + (p.H-z)*T_z / Q

end

p = ComponentArray(Q=-0.5, α = 1.8, H = 1.)

dy0 = [3., 3.]
y0 = [1., 1.]
zspan = [0.01, p.H - 0.001]
prob1 = DAEProblem(ddz!, dy0, y0, zspan, p, differential_vars = [true,true])

solve(prob1,IDA(linear_solver=:GMRES),abstol=1e-9,reltol=1e-9)

```

I’m having two main problems:

1. Perhaps most problematic, my solution changes when I make (seemingly) innocuous changes to `ddz!`. For example, the above example gives u(1 - 0.001)=8.37. However, if I multiply the last line of the function by a negative:  
`res[2] = - U_z - (p.H-z)*T_z / Q`,  
I get u(1 - 0.001)=0.094. The results continue to change if I try to multiply the second equation by `Q` or manipulate the first equation in any way.  
This seems to stop happening for different linear\_solver options, but I don’t know how to pick one.

2. I think I need quite low tolerances. Increasing them to, say, `1e-4` changes the final `u`. However, for different parameters in the full problem, I sometimes get convergence errors even for relatively high tolerances.

I was thinking this would all be improved with a manually defined jacobian, but I can’t get it to run.

```julia
function ddz_jac(J,du,u,p,gamma,z)
    Q = p.Q
    T_z,U_z = du
    
    J[1,1] = gamma*z^2*(p.H-z) * (2*T_z + Q^2/(p.H-z)^2 * 2 * T_z / (T_z - Q^2/(p.H-z)^2)) * abs(1-Q^2/((p.H-z)^2*T_z))^2
    J[1,2] = 0.
    J[2,1] = gamma * (p.H-z) / Q
    J[2,2] = gamma
    nothing
end

DAE_F = DAEFunction(ddz!,jac = ddz_jac)
prob2 = DAEProblem(DAE_F, dy0, y0, zspan, p, differential_vars = [true,true])

```

Running `solve(prob2,IDA())` gives me an error:  
`[IDAS ERROR] IDACalcIC Newton/Linesearch algorithm failed to converge.`

And running `solve(prob2,DFBDF())` gives:  
`type Nothing has no field α`, which I don’t understand in this context. I tried manually setting my exponent in the function and jacobian to a fixed number, and I still get this error.

Any help would be greatly appreciated!

---

<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:** [January 14, 2025, 9:27pm UTC](https://discourse.julialang.org/t/improving-dae-performance/124775/2 "2025-01-14T21:27:49Z")

</div>

Have you tried defining this with ModelingToolkit and seeing what structural\_simplify and a mass matrix method would do? Generally that does a lot better than the implicit DAEProblem.

---

<div class="post-metadata">

**Author:** ![clm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/clm/32/45923_2.png) [@clm](https://discourse.julialang.org/u/clm)\
**Post date:** [January 15, 2025, 5:51pm UTC](https://discourse.julialang.org/t/improving-dae-performance/124775/3 "2025-01-15T17:51:46Z")

</div>

I wasn’t aware of structural\_simplify, but this greatly helps convergence. It also gives me many more solver options to play with to find the one that works best for my problem.

Thank you!

---

<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:** [January 16, 2025, 9:20am UTC](https://discourse.julialang.org/t/improving-dae-performance/124775/4 "2025-01-16T09:20:37Z")

</div>

No problem.

And yes the benchmarks very consistently show that the implicit DAE form of IDA is just slower than the transformed version.

[https://docs.sciml.ai/SciMLBenchmarksOutput/stable/DAE/TransistorAmplifier/](https://docs.sciml.ai/SciMLBenchmarksOutput/stable/DAE/TransistorAmplifier/)

[https://docs.sciml.ai/SciMLBenchmarksOutput/stable/DAE/ChemicalAkzoNobel/](https://docs.sciml.ai/SciMLBenchmarksOutput/stable/DAE/ChemicalAkzoNobel/)

[https://docs.sciml.ai/SciMLBenchmarksOutput/stable/DAE/OregoDAE/](https://docs.sciml.ai/SciMLBenchmarksOutput/stable/DAE/OregoDAE/)

[https://docs.sciml.ai/SciMLBenchmarksOutput/stable/DAE/ROBERDAE/](https://docs.sciml.ai/SciMLBenchmarksOutput/stable/DAE/ROBERDAE/)

There are some clear numerical/mathematical arguments as to why that’s the case as well, but we’re starting to have enough benchmarks to make it clear without requiring any background in how the methods work. So yeah… we’re going to be generally trying to push more things down this path over time, and this is why we haven’t put as much work into DFBDF. IDA is fine, but there’s a fundamental limit to its performance by keeping the `f` definition implicit.

---

<div class="post-metadata">

**Author:** ![clm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/clm/32/45923_2.png) [@clm](https://discourse.julialang.org/u/clm)\
**Post date:** [January 16, 2025, 7:51pm UTC](https://discourse.julialang.org/t/improving-dae-performance/124775/5 "2025-01-16T19:51:37Z")

</div>

Totally makes sense.

I had even looked at those benchmarks, but just thought I was stuck with IDA.

It’s totally possible I’m alone in my misunderstanding as I lack a rigorous background in DEs. But if you think others might fall victim to the same problem, a sentence in the DiffEqDocs faq might be helpful (not to suggest more work for you. Or perhaps later as you continue to push down this path). It’s always my go-to when I get stuck.

---

<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:** [January 17, 2025, 9:26am UTC](https://discourse.julialang.org/t/improving-dae-performance/124775/6 "2025-01-17T09:26:54Z")

</div>

If you missed it than others would. I’ll add something to the FAQ.