David Widmann

I tried to figure out what in DelayDiffEq causes this problem (and wanted to transfer the issue to the DelayDiffEq repo) but it seems it is caused by https://github.com/SciML/OrdinaryDiffEq.jl/blob/2c9f42d7481fe0cc3718257b034b6fe4dc926711/src/initdt.jl#L22-L23 in...

Unfortunately, not all solvers support general number types (such as units) yet, so probably at the moment it's best to specify a suitable solver explicitly (such as in `solve(prob, Tsit5())`)....

> I meant that only the lines with #error behind them actually give an error. So the second line actually works as you also indicate. Sorry for not being clear...

> I ran it and got an error. I am not sure it's fixed? Sorry for the confusion, what I meant was: I get an error but ONLY for the...

It's not obvious to me that `U_smooth` would always return non-negative values since `tanh` (and hence `h` and `h_inv`) could potentially return negative values. However, if `U` could become negative,...

If it's problematic to solve the ODE in the original space and you've already tried `isoutofdomain` and `PositiveDomain` without success, my initial guess would be to perform the computations in...

BTW since `U_smooth` tends to zero and is quite small according to your comment above, it might be beneficial if you evaluate it in log space as well and use...

I just saw that actually your initial value for `u[2]` is zero, so you might want to transform only `u[1]` (initially).

BTW it seems you can even completely decouple the ODEs by modelling `v[2] = u[2] / u[1]` directly (instead of `u[2]`) since `dv[2] = pa * (c * U -...

I just quickly checked and only modelling the ratio `u[1] / u[2]` directly by `dv[1] = pa * (c * U - v[2])` yields the following plots for the final...