Higher Order PDE using Finite Difference How to approximate higher-order partial differential equation using finite difference method?
$$\frac{\partial^{2} y}{\partial t^{2}}+\frac{\partial^{4} y}{\partial x^{4}}=0$$
 A: This is an equation of a dynamic Euler–Bernoulli beam. It can be simply approximated explicitly using standard second and fourth derivative approximation.
$$
\frac{\partial^2 u}{\partial t^2} \sim \frac{u^{n+1} - 2u^{n} + u^{n-1}}{\Delta t^2}\\
\frac{\partial^4 u}{\partial x^4} \sim \frac{u_{m-2} - 4u_{m-1} + 6u_m -4u_{m+1}+u_{m+2}}{\Delta x^4}
$$
Both formulas have truncation error of second order $O(\Delta t^2)$ and $O(\Delta x^2)$ due to symmetry.
There are plenty of methods to derive the approximation for the fourth derivative. The simplest one would be applying second derivative formula twice
$$
\frac{\partial^4 u}{\partial x^4} = \frac{\partial^2}{\partial x^2}\frac{\partial^2 u}{\partial x^2} \sim
\frac{1}{\Delta x^2}\left[
\left(\frac{u_{m+2} - 2u_{m+1}+u_{m}}{\Delta x^2}\right)
-2\left(\frac{u_{m+1} - 2u_m+u_{m+1}}{\Delta x^2}\right)
+\left(\frac{u_{m} - 2u_{m-1}+u_{m-2}}{\Delta x^2}\right)
\right]
$$
The coefficient are the same as in binomial theorem for $(a-b)^4$, since 
$$
\frac{u_{m-2} - 4u_{m-1} + 6u_m -4u_{m+1}+u_{m+2}}{\Delta x^4} = \left(\frac{T_{+1/2} - T_{-1/2}}{\Delta x}\right)^4 u_m
$$
where $T_k$ are shift operators: $T_k u_m \equiv e^{k\Delta x \frac{\partial}{\partial x}}u_m = u_{m+k}$.
So now we can construct a $O(\Delta t^2, \Delta x^2)$ explicit scheme of form
$$
\frac{u^{n+1}_m - 2u^n_m + u^n_{m+1}}{\Delta t^2}
+ K\frac{u^n_{m-2} - 4u^n_{m-1} + 6u^n_m -4u^n_{m+1}+u^n_{m+2}}{\Delta x^4} = 0.
$$
Here $K = \frac{EI}{\rho} = 1$ for your equation.
Let's investigate whether it is stable using von Neumann method. Since the equation is a difference equation with constant coefficients we can seek a solution in form $u^n_m = \lambda^n e^{i\alpha m}$. Plugging that into the equation and canceling $\lambda^n e^{i\alpha m}$ gives
$$
\lambda -2 + \frac{1}{\lambda} + \sigma \left(e^{i\alpha/2} - e^{-i\alpha/2}\right)^4 = 0
$$
where $\sigma$ is the Courant's number for your problem, specifically $\sigma = \frac{K \Delta t^2}{\Delta x^4}$. Changing $e^{i\alpha/2} - e^{-i\alpha/2} = 2\sin \frac{\alpha}{2}$ gives a quadratic equation with real coefficients for $\lambda$
$$
\lambda^2 -2\lambda + 1 + 16\sigma\lambda \sin^4\frac{\alpha}{2} = 0.
\tag{*}
$$
For a scheme to be stable every $\lambda$ for any $\alpha$ should not exceed 1 by its absolute value. If the $(*)$ equation has two real roots they could not both satisfy $|\lambda_{1,2}| < 1$ since their product is $\lambda_1 \lambda_2 = 1$ due to Vieta's formulas. So there left two possibilities. Either $\lambda_{1,2}$ is a complex-conjugated pair or the roots are real and $\lambda_1 = \lambda_2 = \pm 1$. So the discriminant of the quadratic should be negative or equal to zero.
$$
\frac{D}{4} = \left(1 - 8\sigma \sin^2 \frac{\alpha}{2}\right)^2 - 1 = 
-16\left(\sigma \sin^2 \frac{\alpha}{2}\right)\left(1 - 4\sigma \sin^2 \frac{\alpha}{2}\right) \leqslant 0.
$$
So for the scheme to be stable the following condition must hold for every $\alpha$
$$
4\sigma \sin^2 \frac{\alpha}{2} \leqslant 1
$$
which would hold iff
$$
\frac{K\Delta t^2}{\Delta x^4} \equiv \sigma \leqslant \frac{1}{4}\\
\Delta t \leqslant \frac{\sqrt{K}h^2}{2}.
$$
