For a fluid dynamics problem, I need to integrate a first-order differential equation for a fixed time for many points in space: $$ \dot {\mathbf x} = {\mathbf f}({\mathbf x}) $$ Suppose that $\mathbf x$ is two-dimensional. Solving for one trajectory for a fixed amount of time is straightfoward with NDSolve[]. However, I want to compute quantities that require millions of adjacent initial conditions to be integrated for the same amount of time.

One way to do this is to use ParallelTable[]. However, I am curious if there is a way to put fewer calls into NDSolve[] by vectorizing the equations. For example, when using a fixed-step Euler integrator I could treat the entire space of N initial points as a vector of length $2 N$, and then at each timestep update the values in the entire vector. Is there a straightforward way to do this with NDSolve[], and would it yield a performance benefit? If not, are there another tricks for NDSolve[] that may help in this use case?

  • 3
    $\begingroup$ Yes, vectorization is good and you should look at ParametricNDSolve, for get ParallelTable until you have a very efficient NDSolve call. But since you do not provide any code, nothing definite can be said. $\endgroup$ – user21 Sep 15 '16 at 20:08
  • $\begingroup$ For a simple fixed step integration you might get better performance if you forgo NDSolve and implement it yourself. $\endgroup$ – george2079 Sep 15 '16 at 21:30
  • $\begingroup$ @user21 Thank you for the suggestion, I did not know about ParametricNDSolve. I've tried it, and based just on looking at the timing, it looks like ParametricFunction objects just delay evaluation of all the integrals until a value of the ParametricFunction is requested for specific parameter values, and so I get similar timing to running NDSolve within a Table. $f(x)$ is complicated (but analytic) expression in terms of polynomials and common functions. I have tried compiling it as well. How else would parallelization be accomplished using NDSolve? $\endgroup$ – wil3 Sep 16 '16 at 11:18