Solving an ordinary differential equation on a fixed grid (preferably in python)
<p>I have a differential equation of the form </p>

<p><code>dy(x)/dx = f(y,x)</code></p>

<p>that I would like to solve for <code>y</code>. </p>

<p>I have an array <code>xs</code> containing all of the values of <code>x</code> for which I need <code>ys</code>. </p>

<p>For only those values of <code>x</code>, I can evaluate <code>f(y,x)</code> for any <code>y</code>. </p>

<p>How can I solve for <code>ys</code>, preferably in python?</p>

<h2>MWE</h2>

<pre><code>import numpy as np

# these are the only x values that are legal
xs = np.array([0.15, 0.383, 0.99, 1.0001])

# some made up function --- I don't actually have an analytic form like this
def f(y, x):
    if not np.any(np.isclose(x, xs)):
        return np.nan
    return np.sin(y + x**2) 

# now I want to know which array of ys satisfies dy(x)/dx = f(y,x)
</code></pre>

<p>Assuming you can use something simple like Forward Euler...</p>

<p>Numerical solutions will rely on approximate solutions at previous times.  So if you want a solution at <code>t = 1</code> it is likely you will need the approximate solution at <code>t&lt;1</code>.</p>

<p>My advice is to figure out what step size will allow you to hit the times you need, and then find the approximate solution on an interval containing those times.</p>

<pre><code>import numpy as np

#from your example, smallest step size required to hit all would be 0.0001.

a = 0 #start point
b = 1.5 #possible end point
h = 0.0001
N = float(b-a)/h

y = np.zeros(n)
t = np.linspace(a,b,n)
y[0] = 0.1 #initial condition here

for i in range(1,n):
   y[i] = y[i-1] + h*f(t[i-1],y[i-1])
</code></pre>

<p>Alternatively, you could use an adaptive step method (which I am not prepared to explain right now) to take larger steps between the times you need.</p>

<p>Or, you could find an approximate solution over an interval using a coarser mesh and interpolate the solution.</p>

<p>Any of these should work. </p>

python math mathematical-optimization numerical-methods differential-equations