# Two-dimensional wave propagation in layered periodic media

###### Abstract

We study two-dimensional wave propagation in materials whose properties vary periodically in one direction only. High order homogenization is carried out to derive a dispersive effective medium approximation. One-dimensional materials with constant impedance exhibit no effective dispersion. We show that a new kind of effective dispersion may arise in two dimensions, even in materials with constant impedance. This dispersion is a macroscopic effect of microscopic diffraction caused by spatial variation in the sound speed. We analyze this dispersive effect by using high-order homogenization to derive an anisotropic, dispersive effective medium. We generalize to two dimensions a homogenization approach that has been used previously for one-dimensional problems. Pseudospectral solutions of the effective medium equations agree to high accuracy with finite volume direct numerical simulations of the variable-coefficient equations.

## 1 Introduction

Consider the propagation of acoustic waves in a two-dimensional medium whose properties vary in one coordinate direction (say, ). Such waves are described by the PDE

(1) |

Here is the pressure, is the bulk modulus, and is the material density. We focus on the initial value problem in an unbounded spatial domain. We are interested in materials whose spatial variation is periodic:

Here denotes the period. In all numerical experiments and plots, we set .

A simple example of such a medium is shown in Figure 1. We refer to these as layered materials, though the coefficients need not be piecewise-constant. In subsequent sections, we frequently use the terms normal propagation and transverse propagation to refer to propagation normal to or parallel to the axis of homogeneity, respectively (see Figure 1).

We consider the propagation of waves with characteristic wavelength over a distance in a periodic medium with period where

Because the wavelength is larger than the material period , the waves “see” the medium as nearly homogeneous and travel at an effective velocity related to averages of the material properties. The study of wave propagation in this regime has been the subject of much study; see [2, 5] and references therein. Many works focus exclusively on the lowest-order terms in the homogenized equations. The potential for dispersive higher-order terms due to material periodicity was derived using Bloch expansions in [12], and computed explicitly for the case of a one-dimensional layered medium. Later works have also studied effective dispersion in one-dimensional periodic media, and further developed the relevant high-order homogenization techniques for time-dependent problems [4, 1, 3, 14].

The motivation for the present work comes from the discovery of a new kind of effective dispersion, which is described and demonstrated briefly in the rest of the introduction. In Section 2 we derive an effective medium approximation to the variable-coefficient wave equation (1). The technique we use is an extension of those appearing in [1, 14, 10]. This seems to be the first explicit application of such high-order homogenization to multidimensional materials. In Section 3 we explore the dispersion relation implied by the effective medium equations, showing that the effective medium is anisotropic and dispersive. In Section 4, we examine the complementary roles played by variation in the impedance and the sound speed. We also validate our homogenized model by comparing pseudospectral solutions of the effective medium equations with finite volume direct simulations of the variable-coefficient wave equation (1).

All code used for computations in this work, along with Mathematica worksheets used to derive the homogenized equations, are available at http://github.com/ketch/effective_dispersion_RR.

### 1.1 Effective dispersion in a layered medium

The qualitative behavior of waves propagating in a periodic medium depends on whether the sound speed and the impedance vary in space. In general both will vary, but special choices can be made so that either is constant. Figure 2 shows the typical behavior in each of the four possible types of media: homogeneous (top left); constant and variable (top right); constant and variable (bottom left); variable and (bottom right). Each medium consists of alternating horizontal layers:

(2) |

In each case, the solution shown corresponds to one quadrant of the evolution of a an initially Gaussian perturbation:

(3) |

with . The line plots show traces of the solution along the lines (normal) and (transverse). The solutions are computed using highly-resolved finite volume simulations.

Two important effects are evident. First, the speed of wave propagation is anisotropic in the heterogeneous media; typically, normally-incident waves travel more slowly. This is not surprising, given that such waves undergo partial reflection at each material interface. This effect is explained well by the lowest order homogenization theory. The second effect is that, depending on the nature of the medium, different components of the wave may develop a dispersive tail. These dispersive effects cannot be described by the lowest-order homogenization.

For waves propagating in the normal direction (), we observe that dispersion occurs when the material impedance varies (bottom plots) and not when the impedance is constant (top plots). Indeed, the propagation of a plane wave along the -axis (normal propagation) reduces to a one-dimensional problem that is well-studied. Over long distances, periodic variation in the material impedance leads to a dispersive effect – higher frequencies travel more slowly, due to reflection [12]. We refer to this as reflective dispersion. On the other hand, when the impedance does not vary, there is no reflection and no effective dispersion [12, 10]. Instead, all wavelengths travel at the harmonic average of the sound speed.

Next, let us examine the propagation of waves in the transverse direction (). From the left plots, we see that such waves undergo no dispersion when the material sound speed is constant. Remarkably, the right plots show that transverse waves are dispersed when the sound speed varies. In this work we show that diffraction can play a role similar to that of reflection in periodic media, leading again to a dispersive effect in which higher frequencies travel more slowly. Thus an effective dispersion arises even in materials with constant impedance. This diffractive dispersion is an inherently multidimensional effect, with no one-dimensional analog. Whereas reflective dispersion depends on variation in the material impedance, diffractive dispersion depends on variation in the material sound speed.

In Figure 3, we plot streamlines of the velocity field superimposed on a color plot of the pressure, the dashed lines represent the material interfaces. The presence of diffraction is evident in the velocity streamlines. It should be noted that the streamlines do not represent particle trajectories; they are merely a helpful tool for visualizing the vertical velocity components created by diffraction.

The dotted black lines in the slice plots are based on pseudospectral solutions of the effective medium equations (32), derived in Section 2 (no corresponding line is shown in the top left quadrant, since the “effective medium” is exact in the homogeneous case). At the resolution shown in the figure, they are indistinguishable from the “exact” finite volume solutions. The agreement between the effective and variable-coefficient equations is explored in Section 4.

## 2 Homogenization

Homogenization theory can be used to derive an effective PDE for waves in a periodic medium when the wavelength is larger than the period of the medium . The effective PDE is derived through a perturbation expansion, using as a small parameter. We will see that the homogenized PDE depends only on and not on . Since the effective PDE has constant coefficients, it can be used to determine an effective dispersion relation for plane waves in the periodic medium.

The lowest-order homogenized equation for (1) (containing only terms of ) is well understood already, but it contains no dispersive terms and so cannot describe even qualitatively the results shown in the introduction. In this section we derive homogenized equations including terms up to , which include dispersive terms. Additional terms up to are derived for a plane wave propagating in the -direction. Our approach is based on the technique used in [3] for one-dimensional wave propagation. We will see that some care is required to extend the technique used there to the case of two-dimensional media of the type shown in Figure 1. We remark that further difficulties arise when considering periodic media in which the coefficients depend on both and ; we do not pursue high-order homogenization for such materials here.

It will be convenient to deal with the wave equation (1) in first order form:

(4a) | ||||

(4b) | ||||

(4c) |

Here are the velocities in the and coordinate directions respectively.

We start by introducing a fast scale , and adopting the formalism that and are independent. The scale defined by is the scale on which the material properties vary, so we formally replace the -periodic functions with -periodic functions and , which are independent of the slow scale . The dependent variables are assumed to vary on both the fast and slow scales, and are assumed to be periodic in with period . The idea now is to average over the fast scale to obtain constant-coefficient equations involving only . These equations will not capture the details of the solution on the fast scale, but will include (to some degree) the influence of the fast scale on the slow scale.

Using the chain rule we find that . Therefore, system (4) becomes:

(5a) | ||||

(5b) | ||||

(5c) |

For simplicity, from here on, we omit the hats over the coefficient function and , with the understanding that all material coefficients are -periodic functions of . Next we assume that , and can be formally written as power series in ; e.g. . Plugging these expansions into (5) yields

(6a) | ||||

(6b) | ||||

(6c) |

where denotes differentiation of with respect to . Next we equate terms of the same order in ; at each order we apply the averaging operator

(7) |

to obtain the homogenized leading order system and corrections to it. Thus the homogenized equations don’t depend on the fast scale . Note that the averaging operator averages over one period in ; i.e., from to ; therefore, in it averages from to ; i.e., since , we have

(8) |

For brevity in longer equations, we will sometimes use a bar instead of brackets to denote this average (i.e., ).

We present the derivation of the first two corrections in detail. Since the derivation of higher-order terms is similar (but increasingly tedious), we give the higher-order results without detailed derivations. Most of the process is mechanical, but for each system we must make an intelligent ansatz to obtain an expression for the non-homogenized solution of the corresponding system.

### 2.1 Homogenized system

Taking only terms in (6) gives

(9a) | ||||

(9b) |

which implies and . Thus the leading order pressure and vertical velocity are independent of the fast scale . Note that we can’t assume that is independent of ; we will see in the following sections that it is not. This is in contrast to the homogenization of similar systems in 1D where, to leading order, all dependent variables are independent of the fast scale [3].

Taking only the terms in (6) gives

(10a) | ||||

(10b) | ||||

(10c) |

Next we apply the averaging operator to (10). This eliminates the terms and , which are periodic with mean zero. We have no way to determine the average of because both and depend on . We therefore divide (10b) by and then apply , yielding

(11a) | ||||

(11b) | ||||

(11c) |

where Here and elsewhere the subscripts and denote the arithmetic and harmonic average, respectively:

We see already that the effective medium is anisotropic, as indicated by the appearance of these different averages of in (11b) and (11c). In particular, we see that plane waves propagating parallel to the -axis travel with speed while plane waves propagating parallel to the -axis travel with speed . We discuss this in more detail in section 3.1.

Combining (10b) and (11b) yields

(12) |

where is a time independent constant. We choose so that . This confirms that varies on the fast scale . More importantly, this indicates that propagation in is affected by the heterogeneity in even at the macroscopic scale.

Next we obtain expressions for , and in (10). To do so, we use the following ansatz:

(13a) | ||||

(13b) |

This ansatz is chosen in order to reduce system (10) to a system of ODEs. Substituting the ansatz (13), the relation (12), and the homogenized leading order system (11) into the the system (10) we get:

(14) | ||||

(15) |

Based on this, it is convenient to choose , and to satisfy

(16a) | ||||

(16b) | ||||

(16c) |

Equations (16) represent boundary value ODEs with the normalization conditions that . Note that , which implies that , and are -periodic. To solve these boundary value problems we must specify the material functions . In Appendix A we show the fast-variable functions , and for a layered medium. It is convenient to introduce the following linear operators (see [10, p. 1554]):

(17a) | ||||||

(17b) |

Then we have

(18a) | ||||

(18b) | ||||

(18c) |

### 2.2 Derivation of system

Taking only the terms in (6) gives

(19a) | ||||

(19b) | ||||

(19c) |

Substituting the ansatz for and from (13) into (19) and averaging gives

(20a) | ||||

(20b) | ||||

(20c) |

where and similarly for and . For any piecewise-constant functions it can be shown that [10]. Thus, for the piecewise-constant materials that we consider in Section 4, we have . These averages also vanish for the sinusoidal materials we will consider. For more general materials, these terms may be non-zero, but in the following we assume they vanish. Then we obtain:

(21a) | ||||

(21b) | ||||

(21c) |

Since the boundary conditions are imposed in the leading order homogenized system (11), system (21) should be solved with homogeneous Dirichlet boundary conditions; therefore, its solution vanishes:

(22) |

Taking , we make the following ansatz for the solutions and :

(23) |

which is chosen in order to reduce system (19) to a system of ODEs. From (19b) we get . The ansatz for from (13b) gives and using the homogenized leading order equation (11b) we get . Finally we get an expression for the non-homogenized solution :

(24) |

Next we substitute the ansatz for and from (13), the ansatz for and from (2.2) and the non-homogenized solution from (24) into system (19). Then we substitute the leading order homogenized system (11) and the coefficients (18a). Finally, in the resulting expression, we make the fast variable coefficients to vanish to obtain

(25a) | ||||

(25b) | ||||

(25c) | ||||

(25d) |

### 2.3 Derivation of system

### 2.4 Higher order corrections

Following similar, but more involved steps, we find the and corrections.

#### 2.4.1 homogenized correction

The third homogenized correction is:

(28a) | ||||

(28b) | ||||

(28c) | ||||

where and similarly for and . The fast-variable functions , , , , , are solutions of the BVPs

(29a) | ||||

(29b) | ||||

(29c) | ||||

(29d) | ||||

(29e) | ||||

(29f) |

For the two types of media considered in this work all coefficients on the right hand side of (28) vanish. Since the boundary conditions are fulfilled by the leading order homogenized system (11), the third homogenized correction vanishes; i.e.,

(30) |

#### 2.4.2 homogenized correction

The fourth correction is given by:

(31a) | ||||

(31b) | ||||

(31c) |

where and similarly for and . Expressions for the coefficients are given in appendix B

### 2.5 Combined homogenized equations

Once we have the homogenized leading order system and the homogenized corrections we can combine them into a single system. This is done by taking and similarly for and . Combining homogenized systems (11), (21), (27), (28) and (31) we obtain:

(32a) | ||||

(32b) | ||||

(32c) |

where expressions for the coefficients are given in appendix B. Unlike the lowest-order homogenized equation (10), this system is dispersive and the dispersion depends on the direction of propagation, see section 3.

In general, each coefficient of a term in (32) contains a matching factor (contained within the and ; see appendix B.1 for an example of the coefficients of the first order solution for a layered medium). As a result, each term on the right hand side of (32) is proportional to (and independent of ). This explains the observation in [10] that the homogenized equations are valid for any choice of the material period . Nevertheless, (32) is only valid for small ; i.e., for relatively long wavelengths . In all numerical simulations in this work, we take .

## 3 Effective dispersion relations

Equation (32) is a linear system of PDEs with constant coefficients. Hence its solutions can be completely described by the dispersion relation, which relates frequency and wavenumber for a plane wave:

(33) |

Here is the amplitude, is the angular frequency and is the wave vector. Let with , where is the direction of propagation (see Figure 1). Because (32) describes a medium that is anisotropic and dispersive, the speed of propagation of a plane wave depends on both the angle and the wavenumber magnitude .

We can combine (32) into a single second-order equation by differentiating (32) with respect to , differentiating (32b) and (32c) with respect to and respectively, and equating mixed partial derivatives. By substituting (33) into the result, we obtain the effective dispersion relation, up to :

(34) |

It is important to keep in mind that although (34) is accurate to , it is obtained by homogenization and is not expected to be valid for wavelengths shorter than , the medium period. This is because is assumed to be small; in particular, .

### 3.1 Effective sound speed

Taking only terms in (34) we obtain the effective sound speed:

(35) |

The effective sound speed, which indicates the speed of very long wavelength perturbations, depends on the direction of propagation. For normally incident waves (), we have , which is the effective sound speed in a 1D layered medium [12]. For transverse waves (), we have . Since the harmonic average is less than or equal to the arithmetic average, long-wavelength normal waves never travel faster than their transverse wave counterparts. This is intuitively reasonable since transverse propagating waves undergo no reflection.

In Figure 4 we plot as a function of for material parameters , and different values of . When , the sound speed is the same in all directions so we obtain the blue line in Figure 4. As increases (corresponding to strong impedance variation and thus more reflection), the effective speed in decreases and we obtain the red, cyan and black lines in Figure 4. In Figure 5 we take the initial condition

(36) |

with and a piecewise-constant medium (2) with , and which gives and , corresponding to the black line in Figure 4. We show the homogenized leading order pressure (left) and the finite volume pressure (right) at , which correspond to the solution of (11) and (1) respectively. The predicted anisotropic behavior is observed in both solutions.

### 3.2 Normally and transversely incident waves

The dispersion relation (34) has some special properties for waves that are aligned with the coordinate axes. Recall the definitions of normal and transverse wave propagation illustrated in Figure 1.

A normally-incident plane wave corresponds to initial data that is constant in . For such waves, (32) reduces to a one-dimensional equation:

(37a) | ||||

(37b) |

This system was obtained previously in [12, 3, 10]. For a piecewise-constant medium, the coefficients on the right hand side are all proportional to the difference of squared impedance; for instance

(38) |

where and .

In general all of the dispersive terms in (37) vanish when the impedance is constant (see Appendix B). This is because no reflection occurs in such media. This explains why the normally-propagating waves in Figure 2 are dispersed in the lower two plots (with variable ) but not in the upper two plots (with constant ).

Transversely-incident plane waves correspond to initial data that is constant in . For such waves, (32) simplifies to

(39a) | ||||

(39b) |

Here we have included additional 6th-order corrections, because this case will be of particular interest in what follows. For a piecewise-constant medium, all coefficients on the right hand side are proportional to the difference of squared sound speeds; for instance

(40) |

where and .