Hi all,

I’m trying to use the julia DifferentialEquations package to integrate an ODE system over gridded data. The system is something like the following:

\frac{dU_1}{dx} = f(\underline{U}, \vec{\underline{Y}}, \underline{p}, x)

\frac{dU_2}{dx} = g(\underline{U}, \vec{\underline{Y}}, \underline{p}, x)

where \underline{U} = \{U_1, U_2\} is the state vector, \underline{p} is the vector containing the ODE parameters, and \vec{\underline{Y}} = \{\underline{Y}_1, \underline{Y}_2, ..., \underline{Y}_n\} is a vector containing the data (say, pressure, temperature, velocity or something like that) known at each grid point.

My first attempt involved constructing cubic interpolations using Interpolations.jl through the data to transform each \vec{\underline{Y}} into \underline{Y}(x), so that all of my known data variables become functions of position only. However, the system is blowing up and still seems fairly unstable (only advances a few iterations). I think there may be errors in my setup causing this, but as I address those, I wanted to make sure that this approach is a good one, and if not, what might I do instead?