So I have an ODE that needs to be solved a few thousand times on MATLAB and am wondering what the most efficient method to use would be. I am changing a constant term each time. My ODE is of the form $$x''= f(x) + g(t) + x'$$ With initial conditions of $x(0)=x0$ and $x'(0)=0$. $$f(x)=\frac{a^2}{4b}-a(\frac{d^2}{dx^2}(\frac{x}{\sqrt{c+x^2}}))$$ $$g(t)=d \cos{\omega t}$$ with $a$ being the constant that is being varied, $b$, $c$, and $d$ are constant for every iteration.
Since $f(x)$ is a nasty term involving derivatives W.R.T. $x$, I am solving it for a bunch of $x$ values before integration and using the value nearest to the previous ODE solution (accuracy here doesn't seem to be a problem), pulling the corresponding $f(x)$ value from an array. The constant term that is varied for each solver iteration is included in $f(x)$. $f(x)$ is solved with a bunch of matrix operations so it's relatively quick but since my $f(x)$ solutions are in an array corresponding to the constant term and x value, I have to solve using ode45 in a for loop which takes a long time (ode45 accepts arrays as systems of ODEs and doesn't have a 'dim' parameter option so a loop must be used as far as I know.) I know the range that the solutions will be in and it is non-stiff. Can anyone recommend a faster way to do this? Thanks.