PDF version of this book

Discretizing the wave equation

A discrete model can approximate a continuous one to any desired degree of accuracy. Developing such approximations is an important field in applied mathematics. These approximations are widely used in quantum mechanics. It is not possible to model a continuous equation on a digital computer. Thus discrete approximations provide the only practical approach to a great many problems.

Discrete wave models can be studied to understand their intrinsic properties and not just how they can be used to model a continuous reality that may not exist. There is little work in this area. No doubt the complexity of such models is one reason. To understand how such models might look we start with the wave equation which might be called the universal equation of physics because it occurs in so many contexts. We can generate a discrete wave model starting with the continuous wave equation.

$$\frac{\partial^2f}{\partial t^2} = c^2\nabla^2f$$

This notation defines a partial differential equation^7.3. In spite of its foreboding appearance it says something simple. Think of $f$ as the level of water in at a single point in a lake. This equation describes how the level changes based on conditions in its immediate neighborhood. The equation applies to every point on the lake so we can use it to model the dynamic behavior of a wave.

The term on the left of the equal sign is the rate at which the level is accelerating^7.4up or down at this point. The term on the right hand sums acceleration across each spatial dimension at the same point^7.5. For the surface of a lake there are two spatial dimensions. $c$ is the velocity of the wave. The equation says the rate at which the level is accelerating in time at a give point is proportional to the sum of the rates at which the level is accelerating across each dimension in space at that point.

The wave equation is the universal equation of physics. It works for light, sound, waves on the surface of water and a great deal more. It is the relativistic Schrödinger equation that describes the quantum mechanical evolution of the wave function of a single particle with zero rest mass^7.6.

There are many ways to discretize the wave equation. One of the simplest is to define a grid or array of points as shown in Figure 7.1. Instead of defining the value of a function everywhere we consider only selected points. The more closely these points are spaced the more accurate an approximation to the continuous case and the more time consuming the computation. To keep things simple we will consider two dimensions in space and one in time. It is straightforward to move to three spatial dimensions.

Indices are used to locate points in the grid. $f_{x,y,t}$ is the location $x,y$ in space at time $t$. This position has four immediate neighbors in space: $f_{x+1,y,t}$, $f_{x-1,y,t}$, $f_{x,y+1,t}$ and $f_{x,y-1,t}$. See Figure 7.1. Similarly it has two immediate neighbors in time: $f_{x,y,t+1}$ and $f_{x,y,t-1}$.

Continuous differential equations are defined by taking the limit^7.7of finitely spaced locations as the distance between points goes to zero. The first order difference is computed by subtracting neighboring values along the relevant dimension (time, x position or y position). The first order difference in time is either $f_{x,y,t} - f_{x,y,t-1}$ or $f_{x,y,t+1} - f_{x,y,t}$. The wave equation does not use the first order difference or rate of change. It uses the second order difference or rate of acceleration. To get the second order difference we compute a difference of differences.^7.8

In generating the difference equation from the differential equation we must take into account the time and distance scale of the points on the grid. For this illustration we combine these constants with the velocity of the wave $c$ generating a new constant $c_d$.

The second order difference is computed by subtracting one first order difference from the other. This is $f_{x,y,t+1} + f_{x,y,t-1} - 2 f_{x,y,t}$. There are two second order spatial differences for the $x$ and $y$ dimensions. They are computed in a similar way and added together. That is what the notation $\nabla^2$ implies. The result is the following finite difference equation.

$$f_{x,y,t+1} + f_{x,y,t-1} - 2 f_{x,y,t} = c_d^2 (f_{x+1,y,t} + f_{x-1,y,t} + f_{x,y+1,t} + f_{x-1,y,t} - 4 f_{x,y,t})$$

Converting to a form that allows us to compute the next point in time ($f_{x,y,t+1}$) from previous values yields the following equation.

$$f_{x,y,t+1} = c_d^2 (f_{x+1,y,t} + f_{x-1,y,t} + f_{x,y+1,t} + f_{x-1,y,t} - 4 f_{x,y,t}) - f_{x,y,t-1} + 2 f_{x,y,t}$$

To fully digitize the equation we must restrict $f$ to a discrete set of values like the integers. The expression $c_d^2$ must be less than one if the difference equation is to approximate the continuous wave equation. Thus we must add an element to the equation to eliminate the possibility of non integer values. For simplicity we truncate toward $0$. This means $1.8$ is truncated to $1$ and $-1.8$ is truncated to $-1$. We use $T$ to represent this truncation function.

Table 7.1 shows the fully discretized finite difference equation. This figure also shows how an initial state evolves for a few time steps using this equation with $c_d^2$ set to $1/4$.

Figure 7.1: Discretizing a continuous equation

Table 7.1: Discretized wave equation example

Time -1
0	0	0
0	100	0
0	0	0

Time 0
0	0	0
0	100	0
0	0	0

Time 1
0	25	0
25	0	25
0	25	0

Time 2
0	0	6	0	0
0	12	25	12	0
6	25	-75	25	6
0	12	25	12	0
0	0	6	0	0

Time 3
0	0	0	1	0	0	0
0	0	4	12	4	0	0
0	4	24	-11	24	4	0
1	12	-11	-50	-11	12	1
0	4	24	-11	24	4	0
0	0	4	12	4	0	0
0	0	0	1	0	0	0

Time 4
0	0	1	4	1	0	0
0	2	13	6	13	2	0
1	13	9	-34	9	13	1
4	6	-34	14	-34	6	4
1	13	9	-34	9	13	1
0	2	13	6	13	2	0
0	0	1	4	1	0	0

Time 5
0	0	0	0	1	0	0	0	0
0	0	0	5	5	5	0	0	0
0	0	8	14	-7	14	8	0	0
0	5	14	-25	-14	-25	14	5	0
1	5	-7	-14	30	-14	-7	5	1
0	5	14	-25	-14	-25	14	5	0
0	0	8	14	-7	14	8	0	0
0	0	0	5	5	5	0	0	0
0	0	0	0	1	0	0	0	0

Time 6
0	0	0	1	2	1	0	0	0
0	0	3	9	2	9	3	0	0
0	3	13	-3	-9	-3	13	3	0
1	9	-3	-34	13	-34	-3	9	1
2	2	-9	13	2	13	-9	2	2
1	9	-3	-34	13	-34	-3	9	1
0	3	13	-3	-9	-3	13	3	0
0	0	3	9	2	9	3	0	0
0	0	0	1	2	1	0	0	0

Time 7
0	0	1	3	2	3	1	0	0
0	1	8	5	-1	5	8	1	0
1	8	5	-22	0	-22	5	8	1
3	5	-22	-4	9	-4	-22	5	3
2	-1	0	9	-15	9	0	-1	2
3	5	-22	-4	9	-4	-22	5	3
1	8	5	-22	0	-22	5	8	1
0	1	8	5	-1	5	8	1	0
0	0	1	3	2	3	1	0	0

$$f_{x,y,t+1} = 2 f_{x,y,t} - f_{x,y,t-1} + T(\frac{f_{x+1,y,t} + f_{x-1,y,t} + f_{x,y+1,t} + f_{x-1,y,t} - 4 f_{x,y,t}}{4})$$

$T$ is truncation toward 0. For example $T(1.7) = 1$ and $T(-12.9) = -12$.

PDF version of this book