Royal Society Publishing

High-frequency homogenization for periodic media

R. V. Craster , J. Kaplunov , A. V. Pichugin


An asymptotic procedure based upon a two-scale approach is developed for wave propagation in a doubly periodic inhomogeneous medium with a characteristic length scale of microstructure far less than that of the macrostructure. In periodic media, there are frequencies for which standing waves, periodic with the period or double period of the cell, on the microscale emerge. These frequencies do not belong to the low-frequency range of validity covered by the classical homogenization theory, which motivates our use of the term ‘high-frequency homogenization’ when perturbing about these standing waves. The resulting long-wave equations are deduced only explicitly dependent upon the macroscale, with the microscale represented by integral quantities. These equations accurately reproduce the behaviour of the Bloch mode spectrum near the edges of the Brillouin zone, hence yielding an explicit way for homogenizing periodic media in the vicinity of ‘cell resonances’. The similarity of such model equations to high-frequency long wavelength asymptotics, for homogeneous acoustic and elastic waveguides, valid in the vicinities of thickness resonances is emphasized. Several illustrative examples are considered and show the efficacy of the developed techniques.

1. Introduction

Many structures are constructed from composites with periodic, or doubly periodic, variations in material parameters. The varied, and sometimes unexpected, wave propagation properties of composites (Milton 2002) have motivated a number of remarkable new applications, including, but not limited to, photonic crystals and microstructured fibres (Joannopoulos et al. 1995; Zolla et al. 2005), as well as the rapidly developing area of metamaterials (Smith et al. 2004). Thus, there is considerable interest in modelling wave propagation through media with regularly spaced defects or inhomogeneities. Direct numerical simulation using finite elements (Zolla et al. 2005) is popular, but can become intensive for large-scale media with many small inclusions even if calculated over a single cell with Floquet–Bloch conditions invoked. In the latter case, progress can also be made using Rayleigh’s multipole method (Rayleigh 1892) and modern extensions and refinements are possible (Movchan et al. 2002; Poulton et al. in press). Another popular procedure is the plane wave expansion method that amounts to expanding all relevant fields over the cell into Fourier series and solving the resulting infinite system of linear equations (Kushwaha et al. 1993; Andrianov et al. 2008).

Complementary to this literature is that of homogenization, which involves taking a medium with rapidly oscillatory material properties on a fine microscale and averaging these out, in some fashion, to obtain an equivalent homogeneous material with effective material parameters. This is naturally convenient as it buries the microstructure into coefficients, and computations (or analysis) are then just performed on the macroscale. The development of traditional techniques of asymptotic homogenization has been strongly focused on recovering the classical limiting behaviours in effective media (Sanchez-Palencia 1980; Bakhvalov & Panasenko 1989). This usually implies slow variation of relevant fields both on the microscale and on the macroscale, which limits the applicability of the associated expansions to low-frequency situations, e.g. Parnell & Abrahams (2006) or Andrianov et al. (2008). While the frequency range of such models can be extended by considering higher order correction terms (Santosa & Symes 1991; Bakhvalov & Eglit 2000; Smyshlyaev & Cherednichenko 2000), the resulting models cannot fully reproduce high-frequency dynamic behaviours characteristic of microstructured materials, such as strong dispersion, the presence of band gaps or negative refraction. The other way to describe this limitation of the traditional homogenization theory would be to say that it is only capable of describing the fundamental Bloch mode at low frequencies. When inclusions are small with respect to the period, it is possible to construct ‘wide-spacing’ approximation schemes capable of reproducing higher Bloch modes (Poulton et al. 2001; McIver 2007).

Traditional homogenization holds the frequency fixed, while the natural small parameter tends to zero, which is incompatible with the rapidly oscillating fields characteristic of higher Bloch modes. Bensoussan et al. (1978) demonstrate that in order to obtain the complete spectrum of Bloch modes, in their notation, one needs to scale frequency as the inverse of the natural small parameter, i.e. to consider the high-frequency regime. Even at leading order, fields obtained in this asymptotic limit oscillate on the microscale, hence motivating the use of the Wentzel–Kramers–Brillouin–Jeffreys (WKB) ansatz by Bensoussan et al. (1978). Intuitively, this may seem to contradict the very definition of what homogenization is usually thought to achieve. Nevertheless, for fields oscillating at the microscale, the variation of the solution from one periodicity cell to another can be very small. The high-frequency homogenization we develop is designed to exploit this situation and aims to model the modulation of the strongly oscillating field. Although very different from the physical point of view, a mathematically similar asymptotic regime is also observed in the so-called ‘double porosity limit’ of high-contrast homogenization, e.g. Arbogast et al. (1990). When the contrast in material parameters is sufficiently large, higher Bloch modes may become part of the low-frequency response; this has been recently studied in the context of wave propagation (Babych et al. 2008; Smyshlyaev 2009).

For periodic, and doubly periodic, media the interest is often in identifying the Bloch spectra and stop band structure of the solution (Zolla et al. 2005). The Bloch spectra at the edges of the irreducible Brillouin zone correspond to standing waves (Brillouin 1953). Our aim is to develop a high-frequency asymptotic procedure based upon perturbing about these standing wave solutions occurring at particular frequencies across a periodic structure. Taking advantage of the scale separation between micro- and macroscales, the standing waves can be considered upon the microscale and then ‘averaged’ to get a problem just upon the macroscale. The existence of a differential operator of this type has been recently proved rigorously (Birman 2004; Birman & Suslina 2006). The upshot of our analysis is that an effective partial differential equation is constructed on the macroscale that is valid for frequencies in the vicinity of the standing wave. The coefficients of this equation involve integrals of the standing wave solutions and an auxiliary solution along an elementary cell on the microscale. The problem is therefore homogenized in the sense that the microscale plays no explicit role in the effective model. However, in contrast to classical homogenization, the asymptotic macroscale equation is strongly dispersive and can be used to model a range of dynamic phenomena characteristic for composite media. The final form of the effective model bears remarkable similarity to high-frequency long-wave equations in asymptotic theories, for elastic and acoustic waveguides, valid in the vicinity of thickness resonances (Berdichevski 1983; Kaplunov et al. 1998; Le 1999; Gridin et al. 2005).

This paper is structured as follows. In §2, we take a two-dimensional structure and develop an asymptotic procedure where we perturb away from standing wave solutions. The exposition of the two-dimensional case is simplified by considering standing waves with period of exactly one cell size. Very similar expansions can be constructed in the vicinity of standing waves with the period of two cell sizes along one or both coordinates; motivated by the form of corresponding boundary conditions, we term such waves anti-periodic. The point is explicitly demonstrated for one-dimensional examples treated in §2a. Illustrative examples aimed at showing the accuracy of the asymptotics for Floquet–Bloch problems are then chosen. A piecewise homogeneous medium, §3a, illustrates the influence of double roots in the Bloch spectra. The piecewise homogeneous string on a Winkler foundation, §3b, has a lowest curve in the spectra that does not pass through the origin, and a low-frequency band gap, allowing us to contrast the asymptotics developed here with the classical theory. A string with periodically varying density considered in §3c leads us to a Mathieu equation for which the leading order solutions can have a degeneracy that must be overcome, and is a non-trivial example for which the asymptotics are pursued using numerical methods. Some concluding remarks are drawn together in §4.

2. General theory

For generality, we consider a regular periodic structure composed of inclusions and a matrix material. The inclusions are defined on a short scale l, that of the microstructure, whereas the overall problem is considered on a macroscale L that can be thought of as a typical wavelength dominating the dynamic response, or an overall dimension of the structure (see figure 1 for an illustration of typical structures). The ratio of these scales, ϵl/L, is assumed small and provides a natural small parameter. The microstructure is characterized by functions Embedded Image and Embedded Image that are periodic in ξ=(x1/l,x2/l) and can be smooth functions, or piecewise continuous; x=(x1,x2) are Cartesian coordinates orientated along the edges of the cell.

Figure 1.

(a) The geometry under consideration showing a typical two-dimensional structure with microstructure of scale l and macrostructure on scale L. (b) A one-dimensional model structure, a piecewise homogeneous string, considered in §2a is shown.

We consider a model wave equation for, say, sheer horizontal waves in anti-plane isotropic elasticity with periodic density Embedded Image and shear modulus Embedded Image, and with time harmonic dependence Embedded Image assumed understood, asEmbedded Image 2.1 on Embedded Image with ∇x as the gradient operator with respect to the x coordinates. The material parameters in equation (2.1) can be measured in their typical unit values, so that Embedded Image and Embedded Image, with the absence of hats indicating non-dimensionality. In this case, equation (2.1) transforms toEmbedded Image 2.2 with Embedded Image a characteristic wave speed. We adopt a multiple scales approach treating the disparate length scales X=x/L, and ξ=x/l as new independent variables with the result that equation (2.2) becomesEmbedded Image 2.3 where the separation of scales is now made explicit. The notations ∇ξ and ∇X denote (∂ξ1,∂ξ2) and (∂X1,∂X2) respectively. Classical homogenization theory is usually concerned with λ≪1, see equation (2.2); the crucial distinction of high-frequency homogenization allows for λO(1).

We develop the methodology for perturbations to standing wave solutions that are periodic on the cell. If we take the periodicity in ξ to be on a square cell with −1≤ξi≤1, for i=1,2, then periodicity conditions are to be taken in ξ along the edges of the cell, namely Embedded Image 2.4and Embedded Image 2.5 The corresponding solution u(X,ξ) will be periodic in ξ, but not necessarily in X. One can replace these periodicity conditions along one or both spatial dimensions with special ‘out-of-phase’ boundary conditions. These anti-periodic boundary conditions lead to solutions that are periodic across two cells in the appropriate direction(s) of ξ. For the sake of clarity, we specialize the two-dimensional derivations to the case periodic in both ξ1 and ξ2 for which equations (2.4) and (2.5) apply. We illustrate the application and consequences of anti-periodic boundary conditions only for the one-dimensional case presented in §2a.

Next, we adopt the ansatz:Embedded Image 2.6 Each ui(X,ξ) for i=1,2,…, is periodic in ξ, cf. equations (2.4) and (2.5). This ansatz assumes variation at both the microscale and macroscale even at leading order, as opposed to the classical homogenization theory for which u0(X,ξ)≡u0(X).

Substituting equation (2.6) into equation (2.3), and equating to zero the coefficients of individual powers of ϵ, we obtain a hierarchy of equations for ui(X,ξ) and λi with associated boundary conditions from the periodicity in ξ. This hierarchy is resolved from the lowest order up. At the lowest order, we obtain an eigenvalue problem forEmbedded Image 2.7 subject to the appropriate periodicity boundary conditions. This gives rise to a discrete spectrum of eigenvalues Embedded Image for which there is no phase shift across a period of the structure and standing wave is formed. If we fix a simple eigenvalue λ0, then the corresponding eigenmode is of the formEmbedded Image 2.8 where U0(ξ,λ0) is a periodic function of ξ, and f0(X) remains to be determined; we introduce λ0 into the argument of the periodic function to emphasize that it depends upon the frequency λ0. This leading order solution is exactly periodic on the cell, hence conforming to the commonly used notion of the unit cell resonance (or microresonance). Notably, the case of coincident eigenvalues can also arise and is discussed in the simpler context for one-dimensional structures in the following section. Generally, the leading order problem would need to be solved numerically, as is the case in classical homogenization. As we desire the corrections Embedded Image, Embedded Image and the function f0(X), we continue with the hierarchy.

At next order, the equation for u1(X,ξ) isEmbedded Image 2.9 and we now invoke an orthogonality condition that involves multiplying equation (2.9) by U0 and integrating over the periodic cell in ξ. This yieldsEmbedded Image 2.10 where Embedded Image denotes integration over the cell. Using a corollary to the divergence theorem, the first integral on the right-hand side is converted to an integral along the edges of the cell, and thereafter we use periodicity of a(ξ) and U0 and it vanishes. We now subtract the integral of equation (2.7), multiplied by u1/f0, over the cell, from equation (2.10), which results inEmbedded Image 2.11 Using Green’s second identity, the left-hand side becomesEmbedded Image 2.12 where ∂S is the boundary of region S and n is the outward pointing normal. From periodicity in ξ of U0,u1 and a, this integral vanishes. Thus, the only non-zero term is that multiplying Embedded Image; therefore, λ1 must be identically zero. An explicit solution for u1(X,ξ), from equation (2.9), with ξ restricted to be in the cell S, is found asEmbedded Image 2.13 where the coefficient, f1, of the homogeneous solution does not appear, to leading order, in the final result and plays no further role here. The vector function Embedded Image satisfies the leading order equation (2.7), that is,Embedded Image 2.14 i.e. each component of V1 is the non-doubly periodic solution of equation (2.7) linearly independent of U0(ξ,λ0). Notably u1(X,ξ) must be periodic in ξ; however, both terms in the square brackets of equation (2.13) are not. Therefore, we select each individual component of V1 to be periodic along one of the ξi and then choose its boundary conditions along the other ξj, ji, in such a way that conditions (2.4) and (2.5) are enforced; the solution can then be periodically continued to the full structure. We set Embedded Image to have periodicity in ξ2 and then periodicity of u1 in ξ1 results in Embedded Image 2.15and Embedded Image 2.16 Thus, Embedded Image is a periodic function in ξ2 that must be found by solving the inhomogeneous boundary value problem given by equation (2.14) subject to boundary conditions consistent with equations (2.15) and (2.16), as well as equations (2.4)2 and (2.5)2.

Similarly, we take Embedded Image to be a solution of equation (2.14) periodic in ξ1 with boundary conditions Embedded Image 2.17and Embedded Image 2.18 The second-order equation for u2(X,ξ) isEmbedded Image 2.19 and this contains both f0(X) and the correction, Embedded Image, to the eigenvalue Embedded Image. Invoking an orthogonality condition, as at the previous order, by multiplying equation (2.19) by U0, subtracting equation (2.7) times u2/f0 from it, and integrating the result over the cell yields an eigenvalue problem for f0 and Embedded Image as the homogenized partial differential equationEmbedded Image 2.20 The components of the matrix tij are given by Embedded Image 2.21 Embedded Image 2.22and Embedded Image 2.23 Thus, given a particular structure, we solve the leading order equation (2.7) to determine λ0, U0(ξ,λ0). Then, solve equation (2.14), with boundary conditions (2.15)–(2.18) and periodicity, to find V1(ξ,λ0). Given these quantities, the differential eigenvalue problem (2.20) is formed, and thus λ2 identified.

We now specialize to one-dimensional or quasi-one-dimensional problems for which this scheme is more readily applied and for which some explicit details are immediately apparent.

(a) One-dimensional periodic media

Assuming a(ξ) is constant, then for one-dimensional structures equation (2.2) simplifies toEmbedded Image 2.24 where we define c2(ξ)=a/ρ(ξ). The assumption that a(ξ) is constant in equation (2.24) is not essential as any homogeneous second-order differential equation can be transformed into this form by an appropriate change of variables.

The independent variables are now X=x/L,ξ=x/l, and the solution u(X,ξ) will be periodic, or anti-periodic, in ξ, but not necessarily in X. The two-scales approach yieldsEmbedded Image 2.25 This equation is solved subject to either ξ-periodicity conditions u(X,1)=u(X,−1) and uξ(X,1)=uξ(X,−1) or, what we have named, anti-periodicity conditions u(X,1)=−u(X,−1) and uξ(X,1)=−uξ(X,−1), which actually result in the solutions being ξ-periodic across two cells. As we shall see shortly, these two cases correspond to standing waves that are either in-phase or completely out-of-phase at the end of each cell and, in Floquet–Bloch theory, to wavenumbers situated at the opposite ends of the Brillouin zone.

The separation of scales is loosely analogous to the long-wave high-frequency asymptotics outlined for infinite, deformed, waveguides (Berdichevski 1983; Kaplunov et al. 1998; Le 1999) and this analogy can be pursued to obtain a physical interpretation of the results and this is discussed in §§2b and 4. As in the two-dimensional theory, we form a hierarchy of equations to be solved order-by-order. To leading orderEmbedded Image 2.26 and if λ0 is a simple eigenvalue of the problem with an appropriate set of boundary conditions, thenEmbedded Image 2.27 In the periodic case, U0(1,λ0)=U0(−1,λ0), U0ξ(1,λ0)=U0ξ(−1,λ0) and in the anti-periodic case, U0(1,λ0)=−U0(−1,λ0), U0ξ(1,λ0)=−U0ξ(−1,λ0). From Floquet theory, all band gaps present in the Bloch spectra of equation (2.26) are bound by the eigenvalues corresponding to periodic or anti-periodic solutions. Hence, the high-frequency asymptotics developed here are expected to provide accurate description of the solution’s band gap structure in a wide range of problems.

To next orderEmbedded Image 2.28 with the compatibility condition giving λ1=0. The solution isEmbedded Image 2.29 where W1(ξ,λ0) is a non-periodic solution of the leading order equation andEmbedded Image 2.30 with the upper (lower) sign being for the periodic (anti-periodic) case and the constant A is chosen to enforce that behaviour. At next order,Embedded Image 2.31 with the compatibility condition givingEmbedded Image 2.32 where the coefficient T isEmbedded Image 2.33 The differential eigenvalue problem (2.32) then encapsulates the essential physics on the macroscale and the correction to the frequency squared, Embedded Image, is its eigenvalue.

If the leading order solution has U0(±1,λ0)=0, then equation (2.29) is not valid; however, this is easy to overcome with equation (2.33) replaced byEmbedded Image 2.34 where A=2U0ξ(1,λ0)/[W1ξ(1,λ0)∓W1ξ(−1,λ0)] and the non-periodic function W1 satisfies the Dirichlet conditions, W1(±1,λ0)=1 (or any non-zero constant).

A degeneracy occurs whenever λ0 is not a simple eigenvalue, as the leading order solution (2.26) must consist of two linearly independent periodic solutionsEmbedded Image 2.35 The compatibility condition for the first-order term then gives two coupled ordinary differential equations (ODEs) for Embedded Image from which the eigenvalue correction Embedded Image is obtained; the correction to the eigenvalue is now linear rather than quadratic. The coupled equations areEmbedded Image 2.36 for i,j=1,2 and ji.

The asymptotic model (2.32) is not uniform in the sense that when material parameters are continuously varied, such that two eigenvalues approach one another, one has then to switch to the model equations (2.36). Examples from the waveguide theory suggest that it is possible to construct uniformly valid composite expansions (Moukhomodiarov et al. in press) to overcome this.

(b) Analogy with homogeneous waveguides

Placing equation (2.32) in dimensional formEmbedded Image 2.37 then this equation governs the one-dimensional long-wave motion of a periodic medium near the resonance frequencies of a cell. It may be subject to macroscopic boundary conditions along the edges and may also involve inhomogeneous terms on the right-hand side corresponding to edge and surface loadings, see Kaplunov et al. (1998) and references therein.

The long-wave equation (2.37) may be interpreted as describing a homogeneous string attached to an elastic foundation. Although the microscale problem is governed by a Helmholtz equation with periodic coefficients, equation (2.24), the proposed model cannot be formally identified as a Helmholtz equation with an averaged wave speed, which contrasts with traditional homogenization theory. At the same time, the model does not contain microscale variables operating only with long-scale phenomena. It is worth noting that non-local equations are also known to occur as homogenized limits for certain high-contrast composite media (Cherednichenko et al. 2006).

As we have mentioned earlier, the model represents an analogue of the high-frequency long-wave theories for elastic and acoustic homogeneous waveguides (e.g. Kaplunov et al. 1998; Gridin et al. 2005). In fact, the problem on a cell (2.26) corresponds to one on the transverse cross section of a thin rod, plate or a shell. In this case, the cell size l in equation (2.24) corresponds to the half thickness and the cell eigenvalue λ0 corresponds to a thickness resonance frequency. A similar analogy also occurs in conventional homogenization theory, which is, in a sense, a counterpart of the classical low-frequency theories for rods, plates and shells associated with the names of Euler, Bernoulli, Kirchhoff and Love (e.g. Love 1944; Landau & Lifshitz 1970; Graff 1975). The discussion above is also applicable to two-dimensional periodic media.

3. Illustrative examples

To illustrate the efficacy of this approach, we turn our attention to some examples. The first example, a periodic piecewise homogeneous material, has the advantage of being explicitly solvable and serves to illustrate typical features of Bloch waves and pass-stop bands in periodic media. For some parameters, double roots occur allowing this feature to also be explored. This example is extended by introducing a parameter leading to the formation of a low-frequency stop band, hence allowing us to compare and contrast the high-frequency theory with a traditional low-frequency homogenization. Our last example is of a string with a periodic density that leads to Mathieu’s equation, allowing us to present a non-trivial application of our asymptotics.

(a) Periodic piecewise homogeneous media

For a material with piecewise constant variation in c(ξ),Embedded Image 3.1 for r a positive constant, an exact solution is readily obtained. Floquet–Bloch conditions (Brillouin 1953; Kittel 1996) are set at ξ=±1; so Embedded Image and Embedded Image, the solution and its derivatives are continuous at ξ=0, and κ is the Bloch parameter. The dispersion relation relating frequency, λ, to Bloch parameter, κ, isEmbedded Image 3.2 This dispersion relation was seemingly first obtained by Kronig & Penney (1931) for electrons in crystal lattices; it also naturally appears in many guises in one-dimensional photonic and phononic crystals (with layers of infinite height), e.g. Movchan et al. (2002) and Adams et al. (2008). There are an infinite number of discrete eigenfrequencies, λ(κ), to the dispersion relation that we denote as λ(n) with n=0,1,2… ordered from the lowest curve upward; typical Bloch spectra are shown in figure 2. The asymptotic technique we describe extracts the behaviour of the dispersion curves near Embedded Image, which is the frequency at the edge of the Brillouin zone where κ=0. The standing waves that occur when κ=0 correspond to solutions periodic across each elementary cell of width 2ϵ. The other end of the Brillouin zone is at 2ϵκ=π and corresponds to standing waves out-of-phase across the cell, the anti-periodic solutions that are periodic on a double cell and the frequencies at which these occur are denoted by Embedded Image. The asymptotics we develop are valid near each Embedded Image and all Embedded Image except the low-frequency fundamental mode passing through Embedded Image. This mode is the one described by classical homogenization, whose characterization is straightforward, Embedded Image, cf. equation (3.14), at β=0, which corresponds to the substitutionEmbedded Image 3.3 in equation (2.24). We compare the classical philosophy with the high-frequency theory at the end of §3b.

We now turn to the asymptotic procedure. The leading order solution is determined asEmbedded Image 3.4 for Embedded Image and Embedded Image with the upper, lower signs for θ=0,π, respectively. We also need a linearly independent solution, which does not satisfy periodicity or anti-periodicity conditions at ξ=±1, Embedded Image. We take it asEmbedded Image 3.5

Note that W1 can be taken to be any solution of the leading order equation with fixed λ0, which is neither periodic nor anti-periodic. Once W1 is found, it can be re-used in the asymptotics for either end of the Brillouin zone.

Substituting these into equation (2.33) gives T0 and Tπ, that is, the T values associated with the periodic and anti-periodic solutions. These two are written using a single formula asEmbedded Image 3.6

where we append a subscript to T to denote its dependence upon Embedded Image and Embedded Image; (Embedded Image) is the value with the upper (lower) sign and Embedded Image.

We now turn our attention to the ODE for f0, namely equation (2.32), and note that the Floquet–Bloch conditions lead to Embedded Image. For the periodic case in ξ, this forces Embedded Image. Thus, we connect the Bloch parameter with the frequency, found from equation (2.32), via Embedded Image and deduce thatEmbedded Image 3.7 Similarly, in the anti-periodic case, we deduce Embedded Image andEmbedded Image 3.8 The same results are found directly by expanding the exact dispersion relation (3.2), thereby validating the asymptotic scheme for an example that can be explicitly analysed; the advantage of the asymptotic procedure is that it is applicable even when the dispersion relation is unwieldy or unavailable. As the derived expansion makes no specific assumptions about the form of the leading order solution, it is valid for arbitrary volume fractions, although the resulting numerical accuracy at fixed wavenumbers may vary.

We order the eigenvalues for Embedded Image from lowest to highest in magnitude and note that the asymptotics do not apply to the lowest eigenvalue, Embedded Image, which passes through the origin and which have been approached separately, but is accurate for high frequencies. Figure 2 shows the band gap structure and associated asymptotics that emanate from the frequencies at the edges of the Brillouin zone. The quadratic corrections are highly accurate even for values of κ far away from the edges of the zone and this is shown in figure 2c. A double root occurs at Embedded Image and the asymptotics shown in figure 2b are the linear corrections from equation (2.36). An essential ingredient of the asymptotic procedure is the identification of the leading order solutions Embedded Image for the standing waves across the structure. Typical curves are shown in figure 3 for the U0 at each edge of the Brillouin zone, and as the frequency increases the solutions gain more spatial structure.

Figure 2.

(a) The dispersion curves from equation (3.2), for r=1/4, shown across the Brillouin zone. The solid line is from the numerical solution, the asymptotics for Embedded Image are represented by the dashed curve and that for Embedded Image are by the dotted curve. Asymptotics for the double root are shown as dot–dashed lines. (b) The detail for the Bloch spectra at the double root at Embedded Image with the asymptotics shown as the dot–dashed lines. (c) Plots Embedded Image on log–log axes to illustrate the accuracy: numerics are given by the solid line and the asymptotics by the dashed line.

Figure 3.

The leading order solutions U0 for λ(1), λ(2), λ(3) from figure 2 shown in (a) at the κ=0 edge of the Brillouin zone (periodic) and in (b) at the 2ϵκ=π edge (anti-periodic). (a) Solid line, Embedded Image; dashed line, Embedded Image; dotted line, Embedded Image. (b) Solid line, Embedded Image; dashed line, Embedded Image; dotted line, Embedded Image.

The degenerate case of a double root is again illustrated in figure 4. In figure 4b, we note that Embedded Image diverges for some values of r and further analysis shows that these values correspond to the cases of a double root for Embedded Image, which occur whenever r=1/m,m (for integer m). For instance, for r=1/3, as in the figure, Embedded Image, and there is no band gap (as shown in figure 4a). The asymptotics are extracted using the coupled equations (2.36) asEmbedded Image 3.9 and as noted earlier, the correction is linear in κ rather than the quadratic behaviour found for distinct roots. A similar double root behaviour occurs at the other end of the Brillouin zone, and is shown in figure 2b, together with the asymptotics from equation (2.36).

Figure 4.

The dispersion curves from equation (3.2), for r=1/3, shown across the Brillouin zone in (a); the asymptotic solution is denoted by the dashed curve while the numerical solution of the exact dispersion relation by the solid curve. This figure shows the presence of a double root for λ0=3π∼9.42 and the asymptotics for this degenerate case come from equation (3.9). (b)–(d) show Embedded Image versus r for n=1…6. Even (odd) n are positive (negative), and dotted (solid) lines that are almost indistinguishable in each pairing on the scale are shown.

(b) A model structure with non-uniform low-frequency behaviour

The piecewise homogeneous structure has its lowest Bloch spectra curve λ(0)(κ) passing through Embedded Image and as it passes exactly through zero, this low-frequency mode is not captured within our asymptotics. The classical homogenization theory describes the behaviour of this low-frequency curve both when it passes through zero and when an extra (small) parameter exists, as in the example to be considered in this section, that moves the curve upward to expose a low-frequency stop band. We now consider an adaption of the piecewise homogeneous model that allows us to contrast the classical theory with our high-frequency asymptotic approach for this lowest curve of the Bloch spectra.

We previously mentioned that the model equations obtained during the high-frequency homogenization of one-dimensional periodic structures describe vibration of an effective homogeneous string on a Winkler foundation. In this context, it is instructive to study homogenization of a periodically piecewise homogeneous string on a Winkler foundation. This setup also describes a striped waveguide (Adams et al. 2008, 2009) of finite thickness with Dirichlet, Neumann or impedance (Robin) boundary conditions along the guide walls. In both cases, the governing equation is equivalent toEmbedded Image 3.10 The additional term −β2u corresponds to the constant elastic restoring parameter in the Winkler model, or to the separation constant that is related to the transverse mode number for the striped waveguide, and its presence creates a low-frequency stop band. The dispersion relation is a minor adjustment of equation (3.2) toEmbedded Image 3.11 with Embedded Image, Embedded Image. For non-zero β, there is no fundamental mode passing through the origin in the Bloch spectra for this model; hence, high-frequency asymptotics are capable of describing every Floquet–Bloch mode. The two-scales analysis goes through precisely as before with Embedded Image deduced from equation (3.11) with κ=0. The asymptotic correction Embedded Image is deduced from equation (2.32) by integrating the leading order solution U0 and the associated W1 to find T(n) from equation (2.33). Results for typical values are shown in figure 5. Importantly, there is no curve passing through the origin and the lowest curve cuts the frequency axis at Embedded Image. The variation of Embedded Image versus change in wave speed, r, is shown in figure 5b for various β. When β=0, Embedded Image always and the curve then passes through the origin and is a conventional low-frequency mode. As β increases, the curves move away from the r axis and there is always a low-frequency band gap bounded above by Embedded Image. If β=0, then the leading order solution Embedded Image is just a flat line; that is, it has no spatial structure and is the state that one would use in classical homogenization, and as β increases, the mode shape varies and begins to take on more structure; the curves in figure 5c are normalized to have max(|U0|)=1. Further details of the Bloch spectra for this example are in Adams et al. (2008, 2009) together with the details of numerical schemes and other asymptotic techniques that can be applied.

Figure 5.

The dispersion curves for the piecewise homogeneous string on a Winkler foundation. (a) is for r=1/4, β=1 showing the absence of the fundamental mode and the curves from the full numerics (solid) versus the asymptotics (dashed). (b) The variation of the lowest frequency cut-off at κ=0, namely Embedded Image, versus the change in wave speed r for various values of β. (c) The variation in the leading order solution, Embedded Image as β increases for fixed r: r=1.5 in (c). (b,c) Solid line, β=0.25; dashed line, β=0.5; dot-dashed line, β=0.75; dotted line, β=1.

The lowest curve of the Bloch spectra can be approximated using quasi-static distributions along the cell corresponding to classical homogenization. This requires low frequencies and Embedded Image for which the appropriate asymptotic ordering is thatEmbedded Image 3.12 The leading order equation then gives u0(X,ξ)=f0(X), which is uniform along the cell; that is, it does not depend upon the small-scale ξ at all and this difference, and the scaling of λ, is a major difference between the high- and low-frequency models. The function f0(X), from the second-order equation, satisfies the differential eigenvalue problemEmbedded Image 3.13 In this classical model, the inverse-squared wave speed c is replaced in equation (3.10) by an averaged quantity 〈1/c2〉, defined in equation (3.3), and homogenization replaces the variable wave speed by this averaged quantity. The high-frequency asymptotics that we employ are fundamentally different; they are not limited by low-frequency quasi-static variations along the cell and do not simply replace periodic inverse wave speed squared by a constant. The asymptotics take the solution of the standing wave and construct an effective parameter, T, as an integral over ξ that involves the wave speed and the standing wave solution.

Returning to the problem at hand, the Bloch conditions yield Embedded Image and thus this lowest curve is given asymptotically in the low-frequency limit asEmbedded Image 3.14 for β≪1 and small κ. We compare this result with the high-frequency asymptotics in figure 6; from (b) and (c), we see that at low frequencies, equivalently small β, equation (3.14) performs well in predicting the dispersion curves. The high-frequency theory is applicable for small κ but diverges from the numerics sooner than the low-frequency results. For large β, figure 6a, the low-frequency theory drifts away from the solid curve, and becomes inaccurate, and the high-frequency results remain accurate and are so over a longer domain in κ.

Figure 6.

A comparison of the low-frequency homogenization formula (3.14) (dotted lines) with the numerical solution (solid) and the high-frequency asymptotics (dashed). In all panels, r=1/4 and in (a), we show λ(0) for β=2 and β=1. Likewise in (b) but for β=0.5 and β=0.25. (c) The variation in Embedded Image versus β predicted by the low-frequency asymptotics (dotted) and from the numerics (solid).

(c) A continuous periodic variation

We now consider a string with periodic variation in density that leads to the wave speed Embedded Image, where α,Θ are positive constants; this variation gives the classical Mathieu equation (Abramowitz & Stegun 1964; McLachlan 1964) and to more readily connect with the standard theory for that equation, we choose the period to be πϵ rather than 2ϵ. To obtain the Bloch spectrum, the Mathieu equationEmbedded Image 3.15 is solved with Floquet–Bloch conditions at ξ=0,π, namely Embedded Image and Embedded Image. The resulting differential eigenvalue problem for λ is solved numerically using a spectral collocation scheme with Chebyshev basis functions (Weideman & Reddy 2000; Adams et al. 2008), and a typical Bloch spectrum is shown in figure 7. In the language associated with Mathieu’s equation, κ is the characteristic exponent and λ2α,λ2Θ are usually, denoted as a,q; usually the characteristic exponent is found in terms of fixed a,q, whereas in the current application, this is reversed with κ known (and α,Θ known) but with λ2 to be determined.

Figure 7.

The dispersion curves for Mathieu’s equation for α=1,Θ=1/2, shown across the Brillouin zone in (a); the asymptotic solution from Embedded Image at the κ=0 edge of the Brillouin zone is represented by the dashed curve, the asymptotics from κ=π by the dotted curve, while the numerical solution is shown as the solid curve. (b) and (c) show the leading order solutions for Embedded Image and Embedded Image, respectively, for n=2,3. (b) Solid line, Embedded Image; dashed line, Embedded Image. (c) Solid line, Embedded Image; dashed line, Embedded Image.

The periodic and anti-periodic leading order solutions Embedded Image and Embedded Image are extracted numerically for κ=0,π, respectively, and W1 is determined as a non-periodic solution of equation (3.15). In the periodic case, the leading order solutions fall into the two linearly independent solutions described in §20.3 of Abramowitz & Stegun (1964): one being characterized by Embedded Image, Embedded Image, Embedded Image (these lead to the n odd cases) and the other by Embedded Image, Embedded Image, Embedded Image (these lead to the n even cases). The n even cases have Embedded Image and the formula for T from equation (2.34) is used. The non-periodic solution W1 is taken to have boundary conditions Embedded Image for the odd case and Embedded Image for the even case. The anti-periodic case also naturally falls into two linear independent solutions and these too are found numerically.

Given the leading order (U0) and non-periodic (W1) solutions, the integrals for T are evaluated numerically. The change of length of the periodic domain leads to minor changes owing to a rescaling of the domain, and after performing this change the corrections Embedded Image to Embedded Image and Embedded Image are readily found and the resulting asymptotics are shown versus the complete numerics in figure 7. The asymptotics are highly accurate near the edges of the Brillouin zone, as illustrated in figure 7a. The curve for λ(0) passes through the origin and the asymptotics for this curve near zero are found using the classical low-frequency approach, as in equation (3.13), and Embedded Image.

The leading order solutions for U0 are shown in figure 7b with the n=2 solution clearly zero at the ends of the domain, and the solutions chosen for W1 are shown in figure 7c. Figure 8a,b shows that as Embedded Image, physically the material variation decreasing, the values of Embedded Image shown increase dramatically and the difference between the consecutive Embedded Image decreases until the band gap disappears and one obtains a double root; the figure shows the results for n=2,3 and similar results hold for higher n and at the other end of the Brillouin zone.

Figure 8.

(a) and (b) show the variation of Embedded Image and Embedded Image versus Θ, for α=1, for n=2,3. (a) Solid line, Embedded Image; dashed line, Embedded Image. (b) Solid line, Embedded Image; dashed line, Embedded Image.

4. Concluding remarks

The high-frequency asymptotic theory that we present extends classical homogenization, breaking free of the static or low-frequency limitation on the solution variation along the cell. The examples chosen show that by perturbing away from the standing wave solutions, the Bloch spectra are identified through a simple differential eigenvalue problem (2.20) in two-dimension and (2.32) in one-dimension. This differential eigenvalue problem is characterized by a constant parameter whose definition involves the integrations over the short scale of the periodic cell and this short scale plays no further role in the problem; the methodology differs from conventional homogenization theory in several critical ways, and the main one is that the basic state has spatial dependence and so the integrated quantities are not simply averaged wave speeds or simple averaged quantities.

Remarkably, the final ODE equation (2.32) is exactly that which arises in the high-frequency long-wave asymptotics in, say, a straight acoustic waveguide (Gridin et al. 2004). Near the thickness resonance frequencies for the waveguide, i.e. near the eigenvalues of the transverse resonance problem for a homogeneous waveguide, a wave bounces across the guide width, forming a near-standing wave that barely propagates along the guide. Therefore, despite being at high frequency, the wavelength is long. In the periodic situation considered in this article, the vision we have of the wave is that it bounces within a periodic cell of the structure with no phase change, or complete phase change, across the period, to leading order, again forming a standing wave barely propagating along the structure. In both situations, the transverse resonance and standing waves, the wave about which we perturb has a large wavelength, but can occur at a high frequency and so they are both amenable to similar asymptotic methods. A major benefit of having an asymptotic theory at hand is that it uncovers the physics and also complements numerical schemes. For instance, by evaluating the parameter T directly from the standing wave solutions, it is possible to both determine the sign and estimate the value of group velocity of Bloch modes near the edges of the Brillouin zone.

In the long-wave high-frequency theory of waveguides, the ODE is often augmented by terms that account for curvature or geometrical variations along the guide (Gridin et al. 2005; Kaplunov et al. 2005) and which may lead to trapping or localization phenomena. Similar adjustments to the periodic theory can be undertaken and this approach leads to a theory identifying localized modes in media with weakly varying periodic behaviour. The two-scales approach outlined in this article provides a general methodology for treating doubly continuous periodic media, and can be extended to discrete periodic models consisting of point masses and springs. Such models are commonplace in solid state physics and they also exhibit band gap phenomena (Brillouin 1953; Kittel 1996). This, and other, extensions of the theory are underway.


The authors thank NSERC (Canada) and the EPSRC (UK) for support through the Discovery Grant Scheme and grant EP/H021302, respectively. J.K. thanks the University of Alberta for its hospitality and support under the Visiting Scholar programme. The authors gratefully acknowledge useful conversations with Sebastien Guenneau, Evgeniya Nolde and Valery Smyshlyaev.


    • Received November 17, 2009.
    • Accepted February 4, 2010.


    View Abstract