Wave Equation Formulations
For systems described by the wave equation in 1d there are three equivalent ways to describe the system. This post explains them and covers their relation to various types of boundary conditions we typically want to apply. It is part of the theory behind the simulation in the post “Interactive Linear Shallow Water Simulation”.
Contents
Formulations in 1D
The wave equation in one dimension reads
$$\ppfrac{y}{t} - c^2 \ppfrac{y}{x} = 0$$It is a second order partial differential equation, meaning that in order to fully describe the state of a system we must specify the value of the function $y$ and its time derivative $\dot y$. Thus, any initial and boundary conditions we want to specify must also be in terms of $y$ and $\dot y$. The meaning of $y$ is clear, but to tie $\dot y$ to something more interpretable we can seek alternative descriptions of the system and then relate these back to $\dot y$.
One-way Wave Equations
The 1d wave equation can be factored as
$$\Big(\pfrac{}{t} - c\pfrac{}{x}\Big)\Big(\pfrac{}{t} + c\pfrac{}{x}\Big)y = 0$$which is the product of the left- and right-moving 1d advection operators, both with velocity $c$. A right-moving solution \(y = g(x - ct)\) is a solution to the rightward 1d advection equation:
$$\Big(\pfrac{}{t} + c\pfrac{}{x}\Big)g(x - ct) = -cg^\prime(x - ct) + cg^\prime(x-ct) = 0$$We are free to choose the order of the application of the left- and rightwards advection operators because they commute (because the partial derivatives commute.) This means the solution to the right-moving advection equation must also be a solution of the wave equation. We can see this is true because we can apply the right-moving advection operator first, and squash it to zero, leaving the left-moving with only zero to operate on.
The leftwards operator will change the right-moving solution, but it remains a right-moving solution:
$$\Big(\pfrac{}{t} - c\pfrac{}{x}\Big)g(x - ct) = -cg^\prime(x - ct) - cg^\prime(x-ct) = -2cg^\prime(x - ct)$$This must be true, since the operators commute and if the leftwards operator gave a non right-moving solution, then we wouldn’t get zero as the wave equation says we must.
Obviously, the same argument works for the left-moving solution.
Because the wave equation is linear, valid solutions are superpositions of left- and right-moving waves (in 1d). Thus a generic solution to the 1d wave equation is
$$y(x, t) = f(x + ct) + g(x - ct) \label{wave_eqn_soln}\tag{1}$$To connect this description back to $\dot y$ we can just differentiate it wrt time
$$\begin{align} \pfrac{y}{t} &= \pfrac{f(x+ct)}{t} + \pfrac{g(x - ct)}{t} \\[0.5em] &= cf^\prime(x+ct) - cg^\prime(x-ct) \end{align}$$This gives us a way to specify our system in terms of $y$ and $\dot y$ but have it display the desired left-/right-moving behaviour.
Here’s a visualization of two left- and right- moving waves and their sum, as well as the time derivative. Drag the sliders to control time, and the respective magnitudes of the left- and right-moving waves.
Density & Velocity
We can also describe the system via the function itself $y$ and an auxiliary velocity field $u$, which are governed by two coupled first order differential equations.
$$\begin{align} \pfrac{y}{t} = c\pfrac{u}{x} \label{first_order_y}\tag{2}\\[0.5em] \pfrac{u}{t} = c\pfrac{y}{x} \label{first_order_u}\tag{3} \end{align}$$The proof that this yields the wave equation is simple, we can just differentiate equation \eqref{first_order_y} wrt time, then use equation \eqref{first_order_u}:
$$\begin{align} \pnfrac{y}{t}{2} = c\pfrac{}{t}\pfrac{u}{x} = c\pfrac{}{x}\pfrac{u}{t} = c^2\pnfrac{y}{x}{2}\qquad\square \end{align}$$Since the equations are symmetric, it’s clear that $u$ is also governed by the wave equation. In this interpretation $y$ can be seen to be a density and $u$ is its flux (or current).
In the post on the shallow-water equations we saw that $y$ was equivalent to the shallow fluid’s height, and $u$ was the fluid’s average horizontal velocity. Similarly, the homogeneous linearized shallow-water equations are exactly equivalent to the equations describing the propagation of density perturbations in a quiescent compressible fluid (see). In this case $y$ is literally the fluid’s density and $u$ its velocity.
What’s the relationship between the velocity $u$ and the left- and right-moving components; $f$ and $g$ respectively? Just substitute in equation \eqref{wave_eqn_soln} to equations \eqref{first_order_y} and $\eqref{first_order_u}$.
$$\begin{align} \pfrac{u}{x} &= \frac{1}{c}\pfrac{}{t}\Big(f(x+ct) + g(x-ct)\Big) = f(x+ct) - g(x-ct) \\[0.5em] \pfrac{u}{t} &= c\pfrac{}{x}\Big(f(x+ct) + g(x-ct)\Big) = c f(x+ct) + c g(x-ct) \end{align}$$We can see from inspection that this implies
$$u = f(x+ct) - g(x-ct)$$As expected this is also a solution to the wave equation, simply with the sign of the right-moving wave flipped.
Let’s look at the velocity on the visualization:
Boundary Conditions
Usually we want to specify one of 3 types of boundary conditions:
- Reflecting (inverting)
- Reflecting (non-inverting)
- Non-reflecting
Imagine we have a boundary at the origin, \(x=0\) with a left-moving wave approaching it from the right. Let’s examine the three types of boundary conditions.
Reflecting (inverting)
The first type of boundary, reflecting (inverting), is achieved by setting $y$ to be fixed at the boundary; $y(0, t) = 0$. To understand why the reflected wave gets inverted we can imagine the following scenario. Picture a fictitious right-moving wave, which mirrors the left-moving wave’s approach to the boundary, but has opposite sign. The two waves will exactly cancel at the origin ensuring $y = 0$.
Mathematically we are saying $g(x-ct) = -f(-(x-ct))$, which gives $y = f(x+ct) - f(ct-x)$. At $x=0$ we get $y(x=0) = f(ct) - f(ct) = 0$ as desired.
Further, the part of the left-moving wave that has moved past the origin becomes part of the “fictitious” image, and similarly the right-moving part in the domain becomes “real”. Because the right-moving image had inverted sign to ensure $y = 0$, the reflected wave is inverted in sign.
If we are implementing a simulation in terms of the left- and right-moving components $f$ and $g$ respectively, then this boundary condition is equivalent to setting the right-moving component to be $g(0, t) = -f(0, t)$. If we are implementing things in terms of $y$ and $\dot y$, then it’s exactly setting $y(0, t) = 0$, which also implies $\dot y(0, t) = 0$. Finally, if we are implementing things in terms of $y$ and $u$ then this boundary condition implies we have $\partial_x u = 0$ from equation \eqref{first_order_y} and the fact $\partial_t y = 0$.
Reflecting (non-inverting)
The second type of boundary, reflecting (non-inverting), is achieved by fixing the $x$-gradient at the boundary; $\partial_x y(0, t) = 0$ and leaving $y$ free to vary. Here we can imagine a similar fictitious mirrored right-moving wave, but this time it receives the same sign as the left-moving wave, thereby ensuring zero gradient at the boundary.
In this case we are requiring $g(x-ct) = f(-(x-ct))$ so $y(x, t) = f(x+ct) + f(ct-x)$, yielding
$$\pfrac{y}{x} = f^\prime(x+ct) - f^\prime(ct-x)$$which is equal to zero at $x=0$ as required.
Again, once the left-moving wave moves past the boundary it becomes fictitious, and the right-moving image become real once it is past the boundary and inside the domain. This time, because the right-moving image received the same sign as the left-moving wave, the reflection is not inverted.
If implementing in terms of $f$ and $g$ then setting $\partial_x y(0, t) = 0$ is equivalent to requiring $g(0, t) = f(0, t)$, as we have seen. From the velocity standpoint it’s equivalent to requiring $u(0, t) = 0$ because we know $u = f - g$.
Non-reflecting
Often for our simulations we actually won’t want any reflection from the boundary. For this case we can require that the boundary obeys the appropriate one-way wave equation. In our example we have a left-moving wave approaching a boundary at the origin, so at the boundary we require
$$\pfrac{y}{t} (0, t) - c \pfrac{y}{x} (0, t) = 0 \label{bc_non_reflecting_y}\tag{4}$$which is directly specified in terms of $\dot y$. Note that the value of $y$ is also fixed by these since it has to obey $y = \int_0^t \dot y dt + y(0, 0)$. In this case there is no right-moving component so we have simply $g(0, t) = 0$.
If we are describing the system in terms of $y$ and $u$, then we can just substitute equations \eqref{first_order_y} and \(\eqref{first_order_u}\) into the boundary condition (i.e. the left-moving advection equation), yielding
$$\pfrac{u}{t}(0, t) - c\pfrac{u}{x}(0, t) = 0 \label{bc_non_reflecting_u}\tag{5}$$The velocity $u$ also just follows the left-moving advection equation. This is not surprising given the symmetry of equations \eqref{first_order_y} and \(\eqref{first_order_u}\).
Note that the boundary conditions here are in terms of $\partial_t y$ and $\partial_t u$, but using $u$ lets us describe the system as two coupled first order equations, so we cannot directly set the time derivatives in any simulations we might run. To achieve these boundary conditions in a simulation we can update the stored values of $y(0, t)$ and $u(0, t)$ using suitably discretized versions of equations \eqref{bc_non_reflecting_y} and \(\eqref{bc_non_reflecting_u}\).
Visualization
Looking at each type of boundary condition on our visualization, we can get a very intuitive sense of what the maths is describing:
Higher Dimensions
In dimensions higher than 1 we clearly can’t use the left-/right-moving description of the system, because there are infinite directions a wave front could be moving in.
However the description in terms of the velocity generalizes simply; the velocity becomes a vector matching the dimension of the system we are interested in.
The first order equations are
We can see easily that $y$ is still governed by the wave equation:
$$\begin{align} \pfrac{}{t} \pfrac{y}{t} &= \pfrac{}{t}\Big(c \nabla \cdot \mathbf u\Big) \\[0.5em] &= c \nabla \cdot \pfrac{\mathbf u}{t} \\[0.5em] &= c^2 \nabla \cdot \nabla y \\[0.5em] &= c^2 \nabla^2 y\qquad\square \end{align}$$In this case \(\mathbf u\) is governed by an equation which looks quite a lot like the wave equation, but which mixes the components via the divergence:
$$\begin{align} \pfrac{}{t} \pfrac{u}{t} &= \pfrac{}{t}\Big(c \nabla y\Big) \\[0.5em] &= c\nabla \pfrac{y}{t} \\[0.5em] &= c^2 \nabla( \nabla \cdot \mathbf u) \end{align}$$Boundary Conditions
The three types of boundary conditions discussed above are simple enough to write down in higher dimensions; we just need to account for the boundary’s normal vector. Their specification is shown in the table below, where the normal unit vector $\hat{\mathbf n}$ points into the domain.
| Boundary Condition | Second Order | Density & Velocity |
|---|---|---|
| Reflecting (inverting) | $y=0$ | $\hat{\mathbf n} \cdot \nabla u = 0$ |
| Reflecting (non-inverting) | $\hat{\mathbf n} \cdot \nabla y = 0$ | $u = 0$ |
| Non-reflecting | $\partial_t y = c \hat{\mathbf n} \cdot \nabla y$ | $\partial_t y = c \hat{\mathbf n} \cdot \nabla y$ $\partial_t u = c \hat{\mathbf n} \cdot \nabla u$ |