A homogenization route towards square cylindrical acoustic cloaks

We analyse the acoustic properties of a square cylindrical heterogeneous orthotropic cloak deduced from a geometric transform. We propose to achieve some of the properties of this cloak with an elastic square coating divided into 256 cracks preserving its four-fold symmetry. The cloaking mechanism renders freely vibrating and clamped inclusions neutral to some extent and is inherently broadband: provided that the working wavelength is at least three times larger than the outermost sectors of the coating, its homogenized elastic parameters display some transverse anisotropy. Finally, we compare such a cloak to an absorbing coating consisting of elliptical inclusions.


Transformed optics, perfect lenses and electromagnetic cloaks
In 1967, the Russian physicist Victor Veselago wrote a visionary paper in which materials with simultaneously negative permittivity (ε) and magnetic permeability (µ) were shown to display a negative refractive index [1]. It was demonstrated by a ray analysis that a slab of such a negative refractive index material can act as a flat lens that images a source on one side to a point on the other. Thirty years afterwards, Sir John Pendry et al in UK [2,3] proposed designs of structured materials which would have effectively negative ε and µ: the meta-materials. The subsequent experimental demonstration of a negative refractive index at GHz frequencies by an American team led by David Smith [4] opened new vistas in optics, once the early stage controversy following these discoveries died out [5,6]. To set the scene, one should note that negative refractive index materials cannot be found in nature, but are structured at subwavelength lengthscales (typically one-tenth of the wavelength). It is therefore possible to regard them as homogeneous and describe their response by dispersive effective medium parameters (see [7,8] for a comprehensive survey on homogenization theory). Potential applications of negative refraction came when Pendry demonstrated in 2000 that the Veselago slab lens not only involves the propagating waves but also the evanescent near-field components of a source in the image formation [9]. While the resolution of such a perfect lens will be necessarily limited by its constituents [10], Pendry was able to show that it is possible to overcome the Rayleigh diffraction limit for a given optical wavelength with a material at hand: a thin layer of silver (the poor man's lens [9]) was experimentally shown to display a resolution of less than a sixth of wavelength by Zhang's group at Berkeley in 2005 [11]. Pendry and Ramakrishna further explored the possibility to create perfect lenses through mapping from Cartesian to generalized coordinate systems, such as polar or spherical ones [12]. Among 3 others, they theoretically investigated isotropic corner lenses, and this was in essence the birth of transformation optics (its genesis can be attributed to [13]). It was later shown that two corners of anisotropic heterogeneous negative index material (NIM) combined in a checkerboard fashion can act as a unique resonator displaying an infinite local density of states (LDOS) in the limit of zero absorption [14]. Such plasmonic resonances associated with infinitely degenerated eigenmodes also occur in two-dimensional (2D) and 3D periodic checkerboards fulfilling certain skew-symmetric properties [15]. Such checkerboards, when finite in one direction, display some subwavelength imaging [16]. Interestingly, one can trap light to certain extent in photonic crystal checkerboards [17,18], but the LDOS is now bound to be finite [19].
In 2006, Milton and Nicorovici proposed to cloak a countable set of line sources using anomalous resonance when it lies in the close neighbourhood (yet outside) of a cylindrical coating filled with a negative refractive index material: a cylindrical perfect lens [20,21]. A geometric route to cloaking was independently proposed the same year by Pendry, et al [22] and Leonhardt and Philbin [23,24] in the ray optics limit. A numerical verification of cloaking was subsequently provided with full electromagnetic wave simulations in [25,26]. The transformation-based solutions to the Maxwell equations in curvilinear coordinate systems first reported in [22,23] enable one to bend electromagnetic waves around arbitrarily sized and shaped solids. In [27], a cylindrical electromagnetic cloak was constructed using specially designed concentric arrays of split ring resonators which enable to meet among others the prerequisite artificial magnetism property. This locally resonant micro-structured cloak was shown to conceal a copper cylinder around 8.5 GHz, as predicted by numerical simulations [25]. Interestingly, there are two different solutions to Maxwell's equations for an object with and without the cloak that have in theory the exact same field distributions outside the cloak, but this is not in contradiction to the uniqueness theorem which only applies to isotropic media [28], see also [29,30] in the electrostatic context. The effectiveness of the transformation-based cloak was demonstrated theoretically and numerically by solving the Schrödinger equation which is valid in the geometric optic limit [23]. Evidence of cloaking at any frequency was analysed theoretically in [31]- [37] and solving Maxwell's equations with finite elements for an incident plane wave (far-field limit) in [25] and a line current source (near-field limit) in [26]. Square and elliptical electromagnetic cloaks were also investigated numerically in [38]- [41]. In [42], a reduced set of material parameters was introduced to relax the constraint on the permeability, necessarily leading to an impedance mismatch with vacuum which was shown to preserve the cloak effectiveness to a good extent. Other routes to invisibility include reduction of backscatter [43], and applications of transformation optics are also in computationally effective ways of modelling domains of infinite extent via reflectionless interfaces [44,45].

Acoustic lenses and inside-out cloaking
In 2000, Liu et al [46] provided the first numerical and experimental evidence of localized resonant structures for elastic waves in 3D arrays of thin coated spheres. This work by the team of Sheng paved the way towards acoustic analogues of electromagnetic meta-materials, such as fluid-solid composites [47,48]. Movchan and Guenneau subsequently proposed to use arrays of cylinders with a split ring cross section as building blocks for 2D localized resonant acoustic and elastic structures [49]- [51]. Li and Chan [52] independently proposed a similar type of negative acoustic meta-material. In a recent work, Fang et al [53] experimentally demonstrated a dynamic effective negative stiffness for a chain double C resonator (DCR) for ultrasonic waves. 4 Milton et al provided a thorough mathematical frame for such effects including cloaking for certain types of elastodynamic waves in structural mechanics [54]. Milton and Norris independently furthered the theory of acoustic cloaking analysing the underlying continuum elastodynamic governing equations [55,56]. Interestingly, similar elastic localization effects have also been observed experimentally by Russell et al at MHz frequencies in sonic band gap crystal fibre preforms with defects, and were shown to enhance Brillouin scattering [57].
Such neutral inclusions have been also studied in the elastostatic context using asymptotic and computational methods in the case of anti-plane shear and in-plane coupled pressure and shear polarizations [58]. Cummer and Schurig demonstrated that acoustic waves in a fluid also undergo the same geometric transform for a 2D geometry [59], which has been since then generalized to 3D acoustic cloaks for pressure waves [60,61]. Such cloaks require an anisotropic mass density.
In this paper, we propose an alternative route towards acoustic cloaks through a new design that combines some of the advantages of cloaking through geometric transforms with some broadband features. In section 2, we review the elastic material properties required for an ideal cylindrical cloak of square cross section, in the case of anti-plane shear waves. In section 3, we propose a design of a multiply perforated (disintegrating) cloak which approximates to a certain extent the ideal orthotropic heterogeneous cloak, and we further explain its mechanism via homogenization theory. In section 4, we study the reflectionless properties of an orthotropic absorbing cloak consisting of an array of elongated inclusions with certain orientation. In section 5, we draw some conclusions on the pros and cons of the two types of cloaks.

Formulation of the problem in the transformed medium
The design of neutral inclusions in the general 3D elasticity case was theoretically investigated by Milton et al [54]. These authors found that all three components of elastodynamic waves were fully coupled in the transformed medium, and in general do not fulfil a Navier equation unless some restrictions are made on the problem [54]. In the present paper, we focus on the case of out-of-plane shear elastic waves propagating within an elastic material of cylindrical geometry: in-plane pressure and shear waves can be studied independently.

A singular geometric transform in Cartesian coordinates
Let us consider a cylindrical elastic object we want to cloak located within a square of side length 2a. We can use the same route to cloaking as that proposed in electromagnetism by Rham et al [38]. For this, we consider a geometric transformation which maps the field within the square where x 1 , x 2 and x 3 are the contracted Cartesian coordinates. Moreover, this transformation maps the field for x ∞ b onto itself by the identity transformation.
This change of coordinates is characterized by the transformation of the differentials through the Jacobian [38]: Note that R(θ) is the rotation matrix through an angle θ = 0, π/2, π, 3π/2 in the x 1 x 2 -plane with the well-known properties: We thus require four different Jacobian matrices to design our cloak. Note that θ = 2 arctan(x 2 /(x 1 + x 2 1 + x 2 2 )). In our cylindrical case, the elastic properties of the transformed medium are described by the matrix which is a representation of the metric tensor. The transformed coordinates replace the materials (often homogeneous and isotropic) by equivalent ones that are inhomogeneous and orthotropic, whose elastic parameters are determined by The inverse T −1 of the matrix T is given explicitly by [38]: where R(θ) is again the rotation matrix through an angle θ = π/2, π, 3π/2 in the x 1 x 2 -plane. The matrix T of (5) was computed for transverse electric (TE) waves in Cartesian coordinates in [38]. This representation still holds in the present situation. Note that α appears only in the last two diagonal components of T −1 , which further displays a fully symmetric block diagonal part (unlike J x x ). Indeed, the considered out-of-plane shear waves decouple from in-plane pressure and shear waves in cylindrical media, so that they satisfy a generalized Helmholtz equation (with weak derivatives) which can be expressed as where µ T is a 2 × 2 positive definite symmetric matrix consisting of the entries {µ i j } i, j=1,2 of µ , and ρ 33 is the third diagonal element of ρ = ρT −1 . Also, I s denotes the imposed out-of-plane displacement at the source location r s . Except for the first simulation, where we consider a point source located within the square cloak and compare its radiation pattern with the same source located within perfectly matched absorptive layers, see figure 1, we will consider a sheet located on one side of the computational domain, since we want to generate a plane wave.  (4) and (5). (b) Same line source radiating within a square absorptive cloak with orthotropic homogeneous elastic parameters µ and ρ defined by (4) and (16). (c) Real part Re(u) of the displacement field along the x 2 -axis for panel (a) (blue curve), for panel (b) (red curve) and in the case of a homogeneous elastic bulk (black curve). The origin of the x 2 -axis is taken at point (2d, 4d) and it ends at point (2d, 0). The normalized frequency of the acoustic source is ωd/c = 3.5.
The result of figure 1(c), even if it is not central to the paper, indicates that the geometric transform for the square cloak leads to singular material parameters inducing a phase shift between a wave radiated in a homogeneous elastic medium, and within the cloak, which was less apparent for a circular cloak [26]. But it is obvious from figure 1(c) that if we now enclose 7 the source within absorptive perfectly matched layers (PMLs) the magnitude of the radiated wave is no longer comparable with that in a homogeneous medium, while the phase shift is much reduced. It might thus be preferable to use PMLs over the cloak, depending upon the application we have in mind, invisibility or stealth technology, the latter creating a shadowed region in the landscape (hence achieving invisibility only on a dark background).
We note that there is no change in the shear-wave speed v in the transformed medium, since the density ρ and the shear modulus µ undergo the same transformation: Here, I denotes the 3 × 3 identity matrix, and ρ and µ are 3 × 3 matrices.

Setup of the normalized finite element model
As already mentioned, thanks to the cylindrical geometry of the problem, elastodynamic waves (governed by the vector Navier equations), decouple into anti-plane shear waves and in-plane coupled pressure and shear waves. In this paper, we focus on the analysis of anti-plane shear waves. The implementation within the finite element package is fairly straightforward. We first multiply equation (6) by a smooth test function v and using Green's formula, we obtain the weak form of the acoustic equation (6) and discretize it using test functions v taking values on nodes of a triangular mesh, and enforcing its value at the location of the source. We note that setting traction free boundary conditions on an inclusion (or crack) amounts to assuming Neumann (natural) homogeneous data, and Dirichlet data are in order for a clamped inclusion. We shall analyse cloaking for both clamped and stress free square inclusions. The finite element formulation was implemented in FEMLAB (COMSOL Multiphysics package), which requires some changes in the readily available boundary conditions. In this paper, we consider some incoming plane wave, which we model as an imposed displacement I s = 1 at the straight interface r s = {(x 1 , 2) : |x 1 | 2} between the isotropic homogeneous bulk elastic material (in our case silica) surrounding the obstacle and the cloak, which is located in a square of side length 2, and the annulus of orthotropic heterogeneous material surrounding it. Importantly, all numerical results presented in this paper were obtained for a bulk material defined by normalized elastic parameters ρ = 1 (mass density) and µ = 1 (shear modulus or the first Lamé coefficient). Let us further assume that the other Lamé coefficient for this material is λ = 2.3 (pressure modulus), even though it does not appear within the governing equation (6) for shear waves. In this case, the elastic material parameters used are those of fused silica, ρ = 2.2 × 10 3 kg m −3 , µ = 16.05 × 10 9 Pa and λ = 31.15 × 10 9 Pa.

A square disintegrating cloak
In this section, we propose a possible route towards broadband acoustic cloaking based on multiple perforations of an elastic bulk of material. The cracks need to preserve the four-fold symmetry of the square cloak introduced in [38] and also become larger as we walk away from the centre, see figure 2. Nevertheless, to demonstrate phenomenologically the way it works (through effective transverse anisotropy), we will simplify its design and divide the cloak evenly into a square array of rectangular voids: our multiple scale analysis requires these simplifications to be carried out. . We note the phase shift behind the cloak, which can be attributed to the acoustical paths followed by the acoustic wave as shown by the red curves on panel (d), which are larger than that of a plane wave propagating in a homogeneous isotropic elastic bulk.

Homogenization model for a multiply perforated cloak at fixed frequency
When the acoustic wave penetrates the structured cloak c , it undergoes fast periodic oscillations: c (a |x i | b, i = 1, 2), is evenly divided into a large number N (η) ∼ η −2 of small sectors ηY , where η is a small positive real parameter. The smaller η, the larger the number of small sectors ηY . The homogenization technique amounts to looking at the limit when η goes to zero, whereas the frequency ω in equation (6) remains fixed. To filter these oscillations, we consider an asymptotic expansion of the displacement field solution of the equation (6) in terms of a macroscopic (slow) variable x = (x 1 , x 2 ) and a microscopic (fast) variable x/η where u (i) : c × Y * −→ C is a smooth function of four variables, independent of η, such that Rescaling the differential operator in equation (6) accordingly as ∇ = ∇ x + 1 η ∇ y , and collecting terms of same powers of η, we obtain the following homogenized problem in the limit when η tends to zero (see also [8]): From equation (9) (where derivatives are understood in weak sense, i.e. they include discontinuities), we deduce the effective transmission conditions for the values u (−) hom and u (+) hom of the homogenized displacement u hom on the inner and outer boundaries ∂ − c and ∂ + c of the cloak c and for its normal derivative (the flux). Our result shows that the displacement field is solution of equation (9) with a shear anisotropic matrix given by Here, A(Y * ) denotes the area of the region Y * surrounding a void in an elementary cell Y of the periodic array, and ψ i j represent corrective terms where n is the unit outward normal to the boundary ∂ S of the freely vibrating inclusion in the cell Y . Furthermore, j , j ∈ {1, 2}, are periodic potentials which are unique solutions (up to an additive constant) of the following two Laplace equations (L j ): which are supplied with the effective boundary condition ∂ j ∂n = −n · e j on the boundary ∂ S of the inclusion. Here, e 1 and e 2 denote the vectors of the basis in Cartesian coordinates.

Numerical illustration: a meta-material for reduction of forward scattering
Modelling of the invisibility cloak has been carried out using the commercial finite elements package COMSOL. Figure 2 shows the wave picture for both the structured (c) and homogenized (d) cloaks in the case of 256 sectors. The ideal homogenized coating and the structured one show that even if not perfect the structure actually works as expected. The displacement field within the homogenized cloak is characterized by the shear matrix where θ = 0, π/2, π and 3π/2. In figure 2(b), we plot the periodic potential 2 used in the numerical computation of the second entry of the matrix [µ hom ] by solving equations (13)- (22) with i = j = 2. It transpires from figure 2(b) that this potential takes larger values along the direction tangent to the square sides, in agreement with the resulting transverse anisotropy. We then compute the diffraction of a plane wave by a square obstacle when it is clamped or freely vibrating, cloaked or not, and rotated or not, see figures 3 and 4. We note that there is a phase shift between the plane waves propagating in homogeneous silica and forward scattered by the structured cloak. This is not surprising since the acoustic path followed by the plane wave is longer in the second case. But unlike for a circular structured cloak [27], the profile of the wave scattered by the obstacle is somehow comparable to what we achieve when we surround it by the square cloak for a clamped obstacle rotated by an angle π/2, and a freely vibrating square obstacle which is rotated or not. The improvement in forward scattering is only convincing for the case of a clamped (not rotated) obstacle, see figure 4(a). We attribute these poor results to the corners of the cloak, which certainly worsen the diffraction, and therefore compete with the effective anisotropy which bends the wave around the square obstacle.

PMLs versus square acoustic cloak
PMLs provide a reflectionless interface between the region of interest and absorptive layers at all incident angles. PMLs were originally introduced by Berenger in electromagnetism, in which case the regions of PML consist of anisotropic absorptive media [44]. In our case, we consider some orthotropic absorptive media defined by where T is a representation of the metric (rank two) tensor with entries given by Here, s 1 , s 2 and s 3 represent complex stretched coordinates along the axes x 1 , x 2 and x 3 defined by in the case of a wave propagating in the x 1 -direction. The real parameters a and b are chosen ad hoc so that the amplitude of the outgoing acoustic wave vanishes smoothly within the PML region which was meshed with a mesh size about a tenth of the wavelength. In our numerical simulations reported in figure 5, we took the normalized elastic parameters for silica ρ = µ = 1 in (16), which combined with λ = 2.3 corresponds to realistic normalized elastic constant (that of fused silica, as mentioned earlier). We further chose a = b = 1 for the working normalized frequency ωd √ ρ/µ = 3.5. From figure 5, we conclude that PMLs reduce considerably the backscattering of a clamped square obstacle irrespective of its orientation. However, absorption inherently leads to a shadow region behind the PMLs, thereby revealing the presence of the clamped obstacle: PMLs can be used in stealth technology, but not for transparency. We now wish to introduce some microstructured metamaterial enabling us to mimic the acoustic property of the ideal PMLs.

A homogenization model for a cloak with elongated absorptive inclusions
Our result shows that the displacement field is now solution of equation (9) with a shear anisotropic matrix given by Here, the brackets denote averaging over the elementary cell Y of the periodic array, and ψ i j represent corrective terms given by which in the case of a two-phase phononic crystal can be simplified to [62,63] ∀i, j ∈ {1, 2}, which are supplied with the effective boundary conditions on the boundary ∂ S of the inclusion. We note that if we now let µ tend to zero inside the inclusion S, the expression of [µ hom ] in equation (19) reduces to with ψ i j given by equation (13). This expression differs from that in equation (19) by a factor A(Y * ), which is known in the literature as a problem of non-commuting limits or singular perturbation in multiparameters homogenization problems: the limits of large wavelengths (compared with the pitch of the array) and large material contrast do not commute [64]- [67].

Numerical illustration: a phononic crystal for reduction of backward scattering
Modelling of the invisibility cloak has been carried out using the commercial finite elements package COMSOL. Figure 6(b) shows the 2D map of the displacement field when a plane wave is incident from above on the phononic crystal whose geometry is given in figure 6(a). Note that when η is small, there are many more ellipses in the cloak. A comparison between the ideal homogenized coating and the structured one shows that even if not perfect the structure actually works as expected: the forward scattering is clearly reduced; thanks to the reflectionless and absorbing nature of the cloak. The displacement field within the homogenized cloak is characterized by the shear matrix where θ = 0[π/4]. In figures 6(c) and (d), we plot the periodic potentials 2 used in the numerical computation of the second entry of the matrix [µ hom ] by solving equations (13)- (22) with i = j = 2. We note that the resulting anisotropy of the homogenized shear matrix in (25) depends upon the orientation of the ellipse: when the ellipse is rotated by an angle θ = π/4, the gradient of the periodic potential 2 clearly rotates by the same amount in figures 6(c) and (d), in accordance with (25). It might be worthwhile to optimize both orientation and eccentricity of the ellipse to reduce the forward scattering by the cloak.

Conclusion
In this paper, we have proposed an alternative route towards cloaking of square obstacles using a structured square acoustic cloak displaying a transversely anisotropic effective shear matrix. Whereas the ideal cloak obtained through the singular geometric transform in [38] convincingly cloaked a square region, it is fair to say that our structured cloak is not as efficient as we would have expected. But it is shown to reduce the forward scattering for freely vibrating and clamped obstacles to certain extent in some cases. The limitation in the cloaking effect could be attributed to the presence of corners in the outer boundary of the cloak, since the circular counterpart of our structured cloak was numerically shown to work over a finite range of frequencies in [68], whereby Neumann boundary conditions were also considered (in the context of transverse electric waves propagating through an infinite conducting micro-structured cloak). However, this will require further theoretical and numerical investigation before our conjecture is verified which lies beyond the scope of the present paper.
Finally, we have also shown through homogenization of an array of elongated ellipses that in the case of anisotropic absorbing layers, the cloak clearly improves the forward scattering, while it fails to reduce the backscattering. We thus conclude that in any event, it seems likely that a circular cloak might improve the cloaking mechanism compared with a square (or polygonal) counterpart, for the former cloak has smoother inner and outer boundaries. The circular design is therefore more promising for practical purposes: 'think outside the square' might be an adequate moral for this invisible story.