Source: wave-equation.pdf

To describe the implementation of the wave equation in Octave we first have to derive an algorithm to numerically solve the differential equation. This is the general one-dimensional case of the wave equation:

\[ \frac{ \partial^2 u(x, t) }{ \partial x^2 } = \frac{1}{c^2} \frac{ \partial^2 u(x, t)}{\partial t^2}\]

Where \(u(x,t)\) stands for the deflection at the point \(x\) and the time \(t\) and \(c\) is the wave propagation speed in the medium. The second derivation can be approximated by a central difference:

\[ \frac{d^2 f(x)}{d x^2}\approx \frac{f(x+\Delta x) - 2f(x) + f(x-\Delta x)}{\Delta x^2}\]

which can be substituted in the first equation:

\[ \frac{u(x+\Delta x,t) - 2u(x,t) + u(x-\Delta x,t)}{\Delta x^{2}} = \frac{1}{c^{2}} \frac{u(x,t+\Delta t) - 2u(x,t) + u(x,t-\Delta t)}{\Delta t^{2}}\]

To maintain a clear view we will rename some functions as follows:

\[ c(x) = u(x,t) p(x) = u(x, t-\Delta t) f(x) = u(x, t+\Delta t)\]

So the equation will simplify to:

\[ \frac{c(x+\Delta x) - 2c(x) + c(x-\Delta x)}{\Delta x^{2}} = \frac{1}{c^{2}} \frac{f(x) - 2c(x) + p(x)}{\Delta t^{2}}\]

Now we have to transform the equation to express f(x). The reason for this step is to get a recursive form, which can be easily implemented into a program.

\[ step 1: \\ c^{2}\Delta t^{2} \frac{c(x-\Delta x) - 2c(x) + c(x+\Delta x)}{\Delta x^{2}} = f(x) - 2c(x) + p(x) \\ step 2: \\ f(x) = 2c(x) - p(x) + r^{2} \lbrack c(x-\Delta x) - 2c(x) + c(x+\Delta x) \rbrack\]

Where \(r=\frac{c \Delta t}{\Delta x}\). This is the final equation which can now be used as an algorithm for a simulation.