Wave scattering by platonic grating stacks

N. V. Movchan, R. C. McPhedran, A. B. Movchan, C. G. Poulton

Abstract

We address the problem of scattering of flexural waves obeying the biharmonic equation by a stack of a finite number of gratings. We express the solution of the scattering problem for a single grating in terms of reflection and transmission matrices, incorporating the effects of both propagating and evanescent incident waves. The plane wave expansion coefficients above and below the grating are linked to multipole coefficients within the grating using the grating sums and the Rayleigh identities. We derive the recurrence procedure giving the reflection and transmission matrices of the stack in terms of those of individual layers. Trapped waves between a pair of gratings are investigated.

1. Introduction

The topic we address here is the scattering of flexural waves by sets of periodic scatterers or gratings. This work is complementary to our recent investigations of the band structure for doubly periodic problems (Movchan et al. 2007; McPhedran et al. 2009). There are a number of reasons why it is important to study scattering problems involving a finite number of gratings. First, in any practical situation, the assemblage of scatterers will be finite, so it is essential to know how effects such as filtering scale with the size of the system. Second, the study of wave problems for finite systems of gratings permits access to the link between scattering matrices and properties of bases of Bloch functions (see, for example, Platts et al.2002, 2003; Botten et al. 2005). This link may be exploited to deduce scattering matrices for the half space filled with grating layers or for randomized arrangements of stacks in which the gratings vary from layer to layer. A third reason is that finite sets of layers may be used to demonstrate properties such as negative refraction or autocollimation, as has been done for photonic crystals (Botten et al. 2005), and may also be exploited to advantage for any other type of wave that has unusual diffractive effects.

The mathematical model we use here is based on the results of analysis of a canonical spectral problem for the biharmonic operator. Specifically, we develop a rigorous scattering model describing the interaction of plane waves with a grating composed of circular inclusions or fixed point scatterers. The solution relies on the expression of the flexural wave function in terms of both sets of plane waves and multipole expansions composed of cylindrical harmonics satisfying either Helmholtz or modified Helmholtz equations. The multipole expansions about individual cylinders are made consistent with quasi-periodicity using the Rayleigh identities (Platts et al. 2002; Botten et al. 2005), and reconstruction equations are derived to transform from the multipole coefficients back to plane wave coefficients. In this way, reflection and transmission matrices can be established for any layer in a stack. We note that the reflection and transmission matrices used here take into account evanescent input waves, as well as the effects of all propagating incident waves.

Using the reflection and transmission matrices for a single grating, we then apply a recurrence procedure to obtain the reflection and transmission matrices for a finite stack. The recurrence relations take into account coupling between gratings by both propagating and evanescent waves. The structure of the recurrence relations is well adapted to the analysis of flexural wave modes trapped between gratings. These trapped modes have been investigated previously by Evans & Porter (2007). The analysis presented here includes the detailed derivation of the conservation of energy and reciprocity relations, which may be used in the evaluation of the accuracy of the numerical computations.

The structure of the paper is as follows. The governing equations, together with the plane wave and multipole representations of solutions, are discussed in §§24. Quasi-periodic sums over multipoles and the Rayleigh identities employing them are discussed in §5. Furthermore, in §6, we derive the reconstruction equations, which link the multipole and plane wave coefficients. Section 7 presents the general form of the conservation of energy and reciprocity relationships. The particular case of pinned inclusions of zero radius and the properties of a single grating consisting of such inclusions are discussed in §§8 and 9. Recurrence relations for stacks of finite number of gratings are derived in §10. Section 11 contains a discussion of properties of a pair of gratings consisting of pinned inclusions of zero radius and, in particular, of trapped modes that can exist between the gratings.

2. Governing equations

Let W be a solution of the scattering problem for the biharmonic operator Embedded Image 2.1 where Embedded Image, h denotes the plate thickness, ρ is the mass density, D=Eh3/(12(1−ν2)) is the flexural rigidity of the plate, E is the Young’s modulus, ν is the Poisson ratio and ω is the angular frequency. The scattering system will be composed of circular cylindrical objects of radius a, although the numerical results will be confined to the limit case of a→0, with the Dirichlet (clamping) conditions being set on the circular contours: Embedded Image 2.2 The solution of (2.1) admits the representation Embedded Image 2.3 where Embedded Image 2.4 and hence W is the superposition of ‘waves’ with real and imaginary wavenumbers.

The grating is aligned along the x-axis, and the structure will be built by adding more parallel gratings and analysed using the recurrence procedure. The frame of reference is shown in figure 1.

Figure 1.

(a) Elementary cell of the grating. (b) Geometry for the reciprocity theorem: return of a reflected order. (c) Geometry for the reciprocity theorem: return of a transmitted order.

The incident field is the plane wave of any of the following two types.

  1. Propagating (evanescent) solution of the Helmholtz equation: Embedded Image 2.5 where Embedded Image χ0 is real and positive for a propagating solution, and it is pure imaginary (with the positive imaginary part corresponding to a reduction of |Wi,H| with a decrease in y) for an evanescent solution. In terms of an angle of incidence θi, Embedded Image.

  2. Evanescent solution of the modified Helmholtz equation: Embedded Image 2.6 where Embedded Image

The total field is represented as Embedded Image 2.7 where Wi is the incident field of one of the two types mentioned earlier, and Ws,H and Ws,M are the scattered fields solving either the Helmholtz or the modified Helmholtz equations. The solution W also satisfies the Bloch quasi-periodicity conditions along the x-axis, as defined by the incident field.

3. Plane wave expansions

Let γ+ and γ be straight lines parallel to the x-axis in the upper and lower half-planes, respectively, placed at an arbitrary distance from the grating. On these lines, we shall use the plane wave expansions of Ws,H and Ws,M as follows.

  1. On γ+, Embedded Image 3.1 where the infinite range of p can be divided into a finite set Embedded Image of propagating waves with χp real and an infinite set Embedded Image of evanescent waves with χp imaginary with positive imaginary part. Similarly, for the evanescent waves corresponding to the modified Helmholtz equation, Embedded Image 3.2

  2. On γ, Embedded Image 3.3 and Embedded Image 3.4

4. Bessel function expansions

In the vicinity of the circular boundary C0, we use the Bessel function expansions Embedded Image 4.1 and Embedded Image 4.2 where An, Bn, En and Fn are unknown multipole coefficients. Owing to the boundary conditions (2.2), these coefficients are related as follows: Embedded Image 4.3 and Embedded Image 4.4

The following Bessel function expansions from Abramowitz & Stegun (1965) will be needed for the incident fields of the two types shown in §2 (see equations (2.5) and (2.6)):

  1. Embedded Image 4.5 where χ0 is real and positive for a propagating wave, and χ0 is imaginary with positive imaginary part for an evanescent wave.

  2. Embedded Image 4.6 where Embedded Image is imaginary with positive imaginary part.

5. Grating sums and Rayleigh identities

The above expansions around the central cylinder x=0 are related to similar expansions near other cylinders x=dn in the periodic grating by the quasi-periodic Bloch factor Embedded Image, where d is the grating period. When forming quasi-periodic Green’s functions for the grating problem, we thus encounter sums of the form: Embedded Image 5.1Embedded Image 5.2 and Embedded Image 5.3Embedded Image 5.4 In the case of l=0, we shall need the following result: Embedded Image 5.5 where the right-hand sum is a cubically convergent sum. Expressions for grating sums over Hankel functions were given by Twersky (1961), and the sums of the K-type are exponentially convergent and can be evaluated by direct summation.

The Rayleigh identities are derived as in Movchan et al. (2007) with the use of the appropriate grating lattice sums in place of the doubly periodic lattice sums (Chin et al. 1994). We have, with Embedded Image denoting the amplitude of the incident wave (2.6) and Embedded Image denoting the amplitude of the incident wave (2.5), Embedded Image 5.6 and Embedded Image 5.7

These identities are uncoupled for solutions of the Helmholtz and the modified Helmholtz equations. The necessary coupling between the two types of waves is provided by the boundary conditions (4.3) and (4.4). The result is an inhomogeneous system of algebraic equations with respect to the multipole coefficients El and Fl: Embedded Image 5.8 and Embedded Image 5.9

6. Reconstruction equations

The system of equations (5.8) and (5.9) can be truncated and solved for the multipole coefficients El and Fl. However, many quantities of physical interest require the knowledge of the plane wave coefficients above and below the grating or the stack of gratings. These can be readily evaluated once the multipole coefficients are known using the reconstruction equations, which we will now derive.

The method for constructing the reconstruction equations is based on the use of Green’s theorem applied to a test function and the part of the total field corresponding to either the Helmholtz or the modified Helmholtz equation.

For the Helmholtz part of the field, the test function is Embedded Image 6.1 where Embedded Image We note that Vm,H has the opposite quasi-periodicity property to the field WH from equation (2.3). Thus, if we apply Green’s theorem to a rectangular cell, with width d equal to the period of the grating and an arbitrary finite extent along the y-axis, the vertical sides of the rectangle together do not contribute to the line integral Embedded Image Let us next consider the integral along γ+, the segment of length d on the line y=const.>a (figure 1a): Embedded Image 6.2 The result is independent of y in the two cases μm=−χm and μm=χm.

Similar remarks apply to the integral over the lower contour γ, the segment on the line y=const.<−a.

The remaining contribution to the integral comes from the boundary C0 of the circular cylinder of radius a. While evaluating this integral, we use the expansion (4.1) together with the following expansion of the test function: Embedded Image We shall also use the Wronskian formula (Abramowitz & Stegun 1965) Embedded Image 6.3

With the calculations done, the final results for the reflection and transmission coefficients, Rm and Tm, are Embedded Image 6.4 and Embedded Image 6.5

For the modified Helmholtz part of the field, the test function is Embedded Image 6.6 where Embedded Image Its expansion in terms of the modified Bessel functions is Embedded Image We proceed as before, with two possible options of Embedded Image, and the required Wronskian formula (Abramowitz & Stegun 1965) is Embedded Image 6.7 The reconstruction equations for the reflection and transmission coefficients, Embedded Image and Embedded Image, of the modified Helmholtz equation take the form Embedded Image 6.8 and Embedded Image 6.9

7. Conservation of energy and reciprocity relations

We now consider conservation of energy and reciprocity relations, which govern the scattering by a single grating. We stress that, while our interest here is in the case of gratings having cylindrical scatterers, the relations we derive will hold irrespective of the shape of the scatterer. We derive these relationships for both propagating and evanescent incident waves and also for both types of reflected and transmitted waves.

We use the following integral identity pertaining to two solutions u and v of a problem of time-harmonic flexural vibrations of a thin plate: Embedded Image 7.1 We apply this relationship to a unit cell of the grating, with ∂Ω comprising two horizontal segments (γ+ and γ), two vertical segments (γl and γr), as well as the boundary C0 of a scatterer. In the calculations that follow, u and v will be chosen to have opposite quasi-periodicity factors, so that the contributions from γl and γr will cancel. We also assume that u and v satisfy the boundary conditions (2.2) and thus the boundary integral along C0 is equal to zero. For such fields, the above integral relation reduces to Embedded Image 7.2 where Helmholtz and modified Helmholtz terms appear separately.

We first apply equation (7.2) to W=WH+WM and its complex conjugate and use the plane wave expansions (2.5), (2.6) and (3.1)–(3.4) to derive the following conservation of energy relation: Embedded Image 7.3 Here δp denotes the amplitude of the Helmholtz type incident wave, whereas Embedded Image is the amplitude of the modified Helmholtz type incident wave, both in the direction specified by αp=α0+2πp/d. The set ΩH contains Helmholtz waves of propagating type, whereas Embedded Image denotes Helmholtz waves of evanescent type. The set Embedded Image comprises all integers p, since all modified Helmholtz waves are evanescent. This energy relationship agrees with that given by Evans & Porter (2007) (see their formula (4.18)) in the case of a propagating incident wave.

We now consider the reciprocity relations applying to the case of the return of the reflected order p*, as shown in figure 1b. These relations are derived using equation (7.2), with u being the original field corresponding to α0, and v is the returned field corresponding to α0′=−α0−2πp*/d. We use the plane wave expansions as mentioned earlier, which now incorporate the following expressions for direction sines and cosines Embedded Image 7.4 We thus find that the integral along γ is zero and the integral along γ+ gives a set of reciprocity constraints for the reflection coefficients, depending on whether the incident wave in the original problem is a propagating or evanescent Helmholtz wave or (evanescent) modified Helmholtz wave: Embedded Image 7.5Embedded Image 7.6Embedded Image 7.7

We now consider the reciprocity relations for the return of the transmitted order p*, as in figure 1c. Applying equation (7.2) to the original and returned fields (corresponding to α0 and Embedded Image respectively), then gives Embedded Image 7.8Embedded Image 7.9 and Embedded Image 7.10

8. Zero radius limit

In the zero radius limit of the algebraic system (5.8) and (5.9), it may be shown that all multipole coefficients, except for n=0, are equal to zero. For small a, this statement is valid to leading order approximation, in which the leading multipole coefficients E0 and F0 satisfy two equations Embedded Image 8.1 and Embedded Image 8.2 The coefficient near (βa)−1 in equation (8.2) must be zero, which implies leading order Embedded Image 8.3 Substituting equation (8.3) into (8.1), we get Embedded Image 8.4 Taking into account that the term in square brackets in equation (8.4) is 1+O(β2a2), we arrive at the result Embedded Image 8.5

According to equations (6.4) and (6.5), the reflection and transmission amplitudes for waves of the Helmholtz type are then Embedded Image 8.6 The amplitudes for waves of the modified Helmholtz type are Embedded Image 8.7 according to equations (6.8), (6.9) and (8.3).

For the discussion of grating stacks, which follows, we will need to characterize the reflection and transmission properties of a single grating for the additional case of a wave incident from below rather than from above. If the value of α0 is the same for the wave coming from below, we can show from the Rayleigh identities (8.1) and (8.2) that the values of E0 and F0 are the same for the two problems. We then need to adapt the argument of §6 for the case when the incident wave is in the region y<0. The results for the reflection and transmission coefficients, using the superscript ‘’ to denote the incidence from below, are: Embedded Image 8.8Embedded Image 8.9

9. Properties of a single grating

In figure 2, we show the normalized energy reflected |R0|2 as a function of β for a single grating of pinned cylinders in the zero radius limit. Results are given for both normal incidence and an angle of incidence of 30°. The normal incidence curve shows a single Wood anomaly (McPhedran & Waterworth 1972), which is associated with the Rayleigh value β=2π (this is the value at which the orders p=±1 become propagating). Before this is the reflection peak at which |R0|=1, employed by Evans & Porter (2007) in their construction of trapped waves for a pair of gratings. The curve for an angle of incidence of 30° has two Wood anomalies (McPhedran & Waterworth 1972), associated with β=4π/3d (p=−1) and β=8π/3d (p=−2). Only the first of these is preceded by a point of total reflectance. Note that the Rayleigh values are associated with zero reflectance points on the curves.

Figure 2.

Normalized reflected energy as a function of β for a single grating of zero radius cylinders, with an incident wave of propagating Helmholtz type, for angles of incidence of 0° (solid line) and 30° (dashed line).

These properties can be understood on the basis of the lattice sum formula (5.5) and the formula for R0 in the zero radius limit, assuming a propagating input wave: Embedded Image 9.1 As β increases, we reach successive singularities of the lattice sums at each Rayleigh value, which occur first in the imaginary part and then in the real part (figure 3). These singularities in the denominator of equation (9.1) ensure that R0 is zero at Rayleigh values. Zeros of transmittance T0 require that R0=−1, which in turn leads to the requirement that as β approaches the lowest Rayleigh value, the following identity holds: Embedded Image 9.2 Now, for small β, the sum on the left-hand side of equation (9.2) is much larger than the sum on the right-hand side owing to the element Embedded Image, while the latter diverges as β approaches the lowest Rayleigh value, so the two must satisfy equation (9.2) in this β range.

Figure 3.

The real (left) and imaginary (right) parts of the denominator in equation (9.1) as a function of β for an incident wave of propagating Helmholtz type, for angles of incidence of 0° (solid lines) and 30° (dashed lines).

Consider next the long wavelength limit, in which β tends to zero, with α0 also tending to zero in such a way as to keep the incident wave of propagating Helmholtz type. Then in the denominator of equation (9.1), the dominant terms are 1/χ0 and Embedded Image, leading to the estimate Embedded Image 9.3 This gives R0≃(−1+i)/2, for normal incidence and Embedded Image, for an angle of incidence of 30°.

It is interesting to compare the reflectance curves in figure 2 with those in figure 4, for the case of an incident wave of modified Helmholtz type. The latter case is essentially a rescaled (diminished) version of the former, with the denominator Embedded Image (larger for off-axis incidence) replacing Embedded Image. Evanescent incident waves cannot be generated in an unbounded region, but can occur in a finite region between interacting gratings, as we now proceed to show.

Figure 4.

Normalized reflected energy as a function of β for a single grating of zero radius cylinders, with an incident wave of modified Helmholtz type, for angles of incidence of 0° (solid line) and 30° (dashed line).

10. Recurrence relations for grating stacks

We consider a stack of Ng gratings with the common period d along the x-axis. The recurrence procedure adopted here for the evaluation of the transmission and reflection matrices of the stack of gratings corresponds formally to that used in Platts et al. (2002). In this earlier investigation, reflection and transmission matrices comprised partitions corresponding to shear and dilatational waves and their cross coupling, while here the division is into Helmholtz and modified Helmholtz waves and their coupling.

For a particular grating, we define the matrices R± and T±, characterizing reflection and transmission of an incident plane wave coming from above or below the grating, respectively. It is necessary to specify the reflection and transmission coefficients Rp and Tp for Helmholtz type plane waves and Embedded Image and Embedded Image for the modified Helmholtz type waves. This must be done for the whole range of incidence channels enumerated by an integer q.

If N stands for the order of truncation of the multipole expansions, then the matrices R± and T± consist of four (2N+1)×(2N+1) blocks arranged as follows: Embedded Image 10.1 Here, the (p,q)th element of Embedded Image stands for the reflection coefficient Rp when the incident wave comes from above the grating, and its only non-zero amplitude component is Embedded Image for channel q. The (p,q)th element of Embedded Image contains Embedded Image for the same incident wave data. Correspondingly, the (p,q)th element of Embedded Image contains Rp when the incident wave comes from above the grating, and its only non-zero amplitude component is Embedded Image for channel q. Finally, the (p,q)th element of Embedded Image stands for the reflection coefficient Embedded Image with the same incident wave as for Embedded Image. The matrix R is filled in the same way, but with the incident wave now coming from below the grating, and the corresponding definitions are valid for the four blocks of the transmission matrices T±.

Next, the reflection and transmission matrices we have introduced have a phase origin at the centre of the grating. When we use the grating as part of a stack, it is convenient to have such matrices with the phase origin half way between the successive elements. This is done using a diagonal propagation matrix Embedded Image. Let η be the distance between the gratings. If the phase origins for reflected and transmitted fields are placed at y=±(1/2)η, respectively, then the entries of the propagation matrix are Embedded Image 10.2 where Embedded Image if u corresponds to a plane wave of Helmholtz type and Embedded Image if u corresponds to a plane wave of modified Helmholtz type. The reflection and transmission matrices with the phase origins at y=±(1/2)η are then Embedded Image 10.3 which can be used for reflection and transmission matrices corresponding to incidence from above or below the grating.

Referring to fig. 3 of Platts et al. (2002), we suppose we know two matrices Embedded Image and Embedded Image characterizing reflection and transmission by a stack of s grating layers. We add a single layer characterized by reflection and transmission matrices Embedded Image and Embedded Image on top of the stack to give a new stack of s+1 layers. The reflection and transmission matrices of this new stack are given by the recurrence relations Embedded Image 10.4 and Embedded Image 10.5 The recurrence procedure is to fill Embedded Image and Embedded Image using the formulations of §§6 or 8 and deduce the reflection and transmission matrices for incidence from below (Embedded Image and Embedded Image) using the symmetry of the grating. These matrices are converted to phase-shifted form using equation (10.3) and used to start the recurrence. After Ng−1 applications of the relations (10.4) and (10.5), we have formed the desired matrices Embedded Image and Embedded Image characterizing the stack composed of Ng gratings. These are rephased to have their phase origins at the centres of the top and bottom layers of the stack using Embedded Image 10.6

11. Modes trapped in a pair of gratings

We now consider the reflection and transmission of a pair of identical gratings, each composed of clamped inclusions of zero radius. Our treatment generalizes that of Evans & Porter (2007), in that we include the effects of multiple plane waves of both Helmholtz and modified Helmholtz types in the interaction between the gratings, rather than just those of zero order.

Using equations (8.6)–(8.9), we can fill in the reflection and transmission matrices of the identical gratings, and thus the matrices padded to take into account the phase changes on propagation across half the distance η between the gratings: Embedded Image 11.1 The reflection and transmission matrices of the pair of gratings are then given by equations (10.4) and (10.5): Embedded Image 11.2

We show in figure 5 the energy transmitted as a function of β through a pair of gratings with unit separation for an angle of incidence θi=30° (compare with figure 2 for a single grating). The Rayleigh anomaly (McPhedran & Waterworth 1972) associated with the onset of the order p=1 is now a very sharp spike in transmittance, followed by a relatively broad transmission resonance. Of course, the most intriguing feature of figure 5 is the very sharp transmission resonance at β=βt=3.58221, where the transmittance attains 100%. The neighbouring half-power points occur at β=3.58187 and β=3.58253, giving a quality factor for the resonance of around 5400. This is extraordinarily large, given that the structure consists of only two gratings, but may be understood with reference to figure 2. Indeed, the resonant value of β corresponds in figure 2 to a point of reflectance quite close to 1 (in fact, 0.99955). Energy entering the region between the gratings is trapped owing to the high grating reflectance and adds up owing to constructive interference, with the plane wave amplitudes between the gratings exceeding the amplitudes outside the gratings by a large factor. The matrices giving plane wave amplitudes for the downward and upward directions between the gratings are in fact given by Embedded Image 11.3

Figure 5.

(a) Normalized transmitted energy as a function of β for a pair of gratings (solid line) and a single grating (dashed line) of zero radius cylinders, separated by a distance equal to unity, with an incident wave of Helmholtz type, for an angle of incidence of 30°. (b) Detail in the neighbourhood of the transmission resonance (five plane waves of both Helmholtz and modified Helmholtz types used in calculations).

Note that the determinant of the matrix Embedded Image attains its minimum value (around 0.0015) for β corresponding exactly to the value (3.58221) for which the transmittance of the grating pair is unity. The fact that this determinant value is small but non-zero shows that the value of β for which a true resonant solution occurs is in fact complex, with real part close to 3.58221, but with a small imaginary part.

As β increases above βt, the transmittance falls rapidly and in fact has a zero at β0=3.59194. This sort of behaviour of the reflectance or transmittance passing through the extreme values of zero and unity in a narrow frequency range is typical of Fano resonances (Fano & Cooper 1968), where a system has a rapid variation with frequency superimposed on a slower variation. In this case, the slower variation is due to the individual grating, as shown in figures 2, 4 and 5, while the rapid variation is due to the phase term incorporated in the propagation matrix Embedded Image of equation (10.2).

As the angle of incidence is changed from θi=30°, the coincidence between the values of β at which the single grating reflectance is unity and the pair of gratings has the transmission resonance becomes less striking. For example, for θi=27°, the energy transmitted for the grating pair and the single grating is shown in figure 6. The energy transmitted for the grating pair resonates at β=3.59999, while the maximum of the single grating reflectance now occurs at β=3.696. This less perfect alignment of these values results in a somewhat smaller quality factor of the resonance (now around 1700).

Figure 6.

(a) Normalized transmitted energy as a function of β for a pair of gratings (solid line) and a single grating (dashed line) of zero radius cylinders, separated by a distance equal to unity, with an incident wave of Helmholtz type, for an angle of incidence of 27°. (b) Detail in the neighbourhood of the transmission resonance (five plane waves of both Helmholtz and modified Helmholtz types used in calculations.)

Acknowledgements

R.C.M. acknowledges support from the Australian Research Council under its Discovery Grants Program and partial travel support from the Research Centre in Mathematics and Modelling of the University of Liverpool. C.G.P. acknowledges support from a start-up grant from the University of Technology, Sydney.

Footnotes

    • Received June 5, 2009.
    • Accepted July 16, 2009.

References

View Abstract