Stormer-Verlet integration was developed specifically to solve Newtons equation F=ma, or:
[itex]\ddot{x}=a(x(t))[/itex] for initial value problems.
Stormer-Verlet integration leads to the relation:
[itex]x_{n+1}=2x_n -x_{n-1}+a_n\Delta t^2 + O(\Delta t^4)[/itex]
That's why you usually only see it in this form.
However, it is easy to extend this to more general cases, although I've never seen people do this. For the general second order ode:
[itex]a(t)\ddot{x}+b(t)\dot{x}+c(t)x=0[/itex],
first transform to get rid of the [itex]\dot{x}[/itex] term:
[itex]\ddot{y}=ry[/itex], with [itex]r=\frac{b(t)^2+2a(t)\dot{b}(t)-2b(t)\dot{a}(t)-4a(t)c(t)}{4a(t)^2}[/itex]
then, proceed as usual to get the recursive relationship:
[itex]y_{n+1}=2y_n -y_{n-1}+r_n\Delta t^2 + O(\Delta t^4)[/itex]
Note that to get the same accuracy as the original Stormer-Verlet method, you have to construct r(t) exactly first, you should not approximate the derivatives of a(t),b(t) and c(t).
By the way, when a,b,c are constants, r reduces to [itex]r=\frac{b^2-4ac}{4a^2}[/itex]
To get the first derivative (velocity), you need to do the usual (leapfrogging or central differencing)
For other differential equations, you need to switch to other methods though.