I have this system of differential equation.
$my''+bx'=f$ ... (A)
$bx'=k(y-x)$ ... (B)
From (B), (A) is
$my''+k(y-x)=f$ ... (C)
For simplicity, $f=0$ .
To solve this numercially, I did this.
By (C), $y$'s value at next time step can calculated by
current time step value of $y$ and $x$ .
Similarly, by (B), $x$'s value at next time step can calculated by
current value of $y$ and $x$ .
Then I did leap frog and calculated time response of $x$ and $y$ .
However, the result diverged.
How can I fix it?
The best course of action is to reduce to a system of first order differential equations and use a standard ode solver on that system.
$$ \begin{cases} y' = z\\ z' = -\frac{k}{m} (y-x) +\frac fm \\ x' = \frac kb (y-x) \end{cases} $$
Edit: With the set of parameters given in the comments below, Wolfram returns a numerical solution with a maximum error in the order of $10^{-8}$. The error is measured with respect to the analytical solution.