Gradient index phononic crystals and metamaterials

1 Introduction

Phononic crystals and acoustic metamaterials enable to achieve innovative properties for the propagation of mechanical waves (air-borne sound waves, water-borne acoustic waves, water waves, elastic waves, surface acoustic waves and Lamb waves) inapproachable in natural materials. Phononic crystals consist of periodic arrangements of scatterers in a given matrix, and they were firstly proposed in 1993 [1], [2]. The attention that phononic crystals received originally was due to the existence of a phononic Bragg band gap where the propagation of acoustic waves was forbidden. In this regime the wavelength λ of the acoustic field is comparable to the periodicity a of the lattice, λ≈a, and to be observable it is also required that the thickness D of the bulk phononic crystals be at least four or five periodicities, D>>a. In the subwavelength range, λ>>a, it is possible to find local resonances in the scatterers, and the designed structures exhibit hybridization band gaps which give rise to novel effects such as negative mass density or negative elastic modulus. In this exotic regime the structure behaves as a special type of materials called “acoustic metamaterials”, which were firstly proposed in the seminal work [3] in 2000. Over the past two decades, dramatically increasing efforts have been devoted to the study of acoustic artificial structured materials driven by both fundamental scientific curiosities with properties not found previously and diverse potential applications with novel functionalities [4], [5], [6], [7], [8], [9], [10], [11], [12], [13], [14], [15], [16].

While interesting, most of the extraordinary properties of metamaterials are in general single-frequency or narrow-band, since outside the resonant regime metamaterials behave as common composites. However, non-resonant phononic crystals in the low-frequency regime behave as homogeneous non-dispersive materials whose effective parameters can be easily tailored, and in this regime gradient index (GRIN) acoustic materials, or GRIN devices, are easily doable. These devices were firstly proposed for acoustic waves by Torrent et al. in 2007 [17] and for elastic waves by Lin et al. in 2009 [18], and they allow to manipulate acoustic waves to enforce them to follow curved trajectories. They are characterized by a spatial variation of acoustic refractive index, which is designed by locally changing the geometry of units. For instance, the effective index obviously depends on the filling ratio of scatterers, namely, the size of scatterers, which is initially proposed for GRIN control of acoustic waves [17], [18], [19]; instead, the variation in the lattice spacing while keeping the size of scatterers is also achievable [20], [21]; for triangular shape of scatterers, rotating the angles of the triangular shape can also affect the effective acoustic velocity [22]; for lead-rubber pillared metamaterial plate, by changing the height of lead layer in pillars, the effective mass density is tuned, resulting in a change in effective phase velocity following a given law [23]; for phononic crystal plates, the effective velocity of antisymmetric Lamb wave is directly related with the thickness of plates, which can also further affect the effective velocity of symmetric mode [24], [25], [26]; and gradient effective index can also be achieved from coiling-space structures [27], [28], [29], [30]. In addition to the geometric parameters, the effective refractive index can also be tuned in relation to the elastic properties of the scatterers, for instance, the choice of materials [18], [31], or to an external stimulus such as electric field [32] or temperature [33].

The attained effective refractive index of the GRIN device can be smaller or larger than that of the background medium, which corresponds to phase advance [20] or delay approaches in wave propagation. Nevertheless, the phase delay approach is mostly employed as higher refractive index enables to design more advanced functionalities and reduce the entire width of the devices. Comparing to the wavelength, the thinner the device, the larger will be the highest required index. If the width of the whole device is downscaled to the subwavelength regime, the GRIN device becomes a metasurface [34], [35]. However, to find a high enough effective refractive index to design GRIN metasurfaces is a big challenge that needs new technologies in material science. Recently, it has been reported that soft porous materials [36], [37] with soft-matter techniques can achieve a relative refractive index higher than 20 [38], resulting in soft GRIN metasurface [39]. GRIN devices can also be designed with negative index of refraction at frequencies lying in the first negative slope of the acoustic band structure [21]. However, these devices would have a narrow band, and it may limit potential applications.

GRIN phononic crystals and metamaterials can be applied to various types of waves in a long frequency limit, such as surface water waves [40], [41] for a few hertz, air-borne sound waves [19], [42], [43], [44], [45], [46], [47], [48], [49], [50], [51] for 10³–10⁵ Hz, water-borne acoustic waves [20], [31], [52], [53], [54] for 10⁴–10⁶ Hz, Rayleigh waves [55], [56], [57], [58], [59] for 10–10⁸ Hz and Lamb waves [23], [24], [25], [32], [60], [61], [62], [63], [64], [65], [66], [67], [68], [69], [70], [71], [72] for 10³–10⁸ Hz, among others, with functionalities like focusing [18], [19], waveguiding [70], [73], mirage [74], beam splitting or deflection [41], [62], [63], cloaking [55], [64] or energy harvesting [71], [72], [75]. It is in the domain of elastic waves in plates where the most interesting applications at the nanoscale are found, since from the technological point of view low-dimensional materials are more interesting than bulk materials. The excellent compatibility of phononic crystals and metamaterials with the nanoelectromechanical systems has been proven for applications in wireless telecommunications, sensing or thermal control, among others [76], [77], [78], [79], [80], [81]. However, the propagation of elastic waves in plates is in general composed of three polarizations, which travel at different speeds and for which refractive devices designed for one of these polarizations will not work for the other two. Recently, the designs of GRIN devices based on phononic crystal plates have demonstrated that the simultaneous control of all fundamental modes propagating in thin elastic plates in a broadband frequency region is possible [25], [26], [63].

The objective of this review is to provide a comprehensive picture of the evolution of the domain of GRIN devices for mechanical waves, to demonstrate the wide variety of applications at the nanoscale that these devices offer and to present the challenges to be accomplished to further develop this field. The paper is organized as follows: after this introduction, we will review the fundamental ideas of the homogenization of sonic and phononic crystals in Section 2 and of phononic crystal plates in Section 3, which offers efficient tools to obtain the effective elastic properties and, consequently, their effective refractive index. Then we will review GRIN devices for bulk acoustic and elastic waves in Section 4 and for flexural waves in Section 5; the advanced full control of polarizations in elastic plates with multimodal GRIN devices will be reviewed in Section 6. The last section will present some conclusions and future challenges.

2 Homogenization of sonic and phononic crystals

Graded materials with specific variations of the refractive index are obviously not found in nature, and they have to be artificially engineered. The inclusion of scatterers in a given matrix is an excellent way for the realization of artificially graded materials, since the average behavior of the scatterers is to modify the effective velocity of acoustic waves in the matrix, and this effective velocity can be tuned by means of the size of the inclusions. When these inclusions are arranged in a regular array we call the composite a sonic (fluid matrix) or phononic (solid matrix) crystal, and the computation of the effective sound velocity is made by a set of mathematical tools called homogenization theories. Roughly speaking, the homogenization region is the range of frequencies at which the propagating field cannot distinguish the individual scatterers of the composite and perceives the structure as a uniform material with some effective parameters (stiffness constant, mass density, viscosity, etc. are the parameters to be computed by means of the most adequate homogenization method).

Homogenization is an old problem, and a great bibliography is available in this realm (see for instance [82]). Our aim here is not to review all these methods but to present some examples of these methods applied to the specific case of sonic and phononic crystals.

2.1 Homogenization of sonic crystals by multiple scattering theory

As mentioned before, a sonic crystal consists of a periodic distribution of solid or fluid inclusions embedded in a fluid matrix. When the operating wavelength is larger than the typical distance between inclusions (lattice constant) the field cannot distinguish individual scatterers and perceives the structure as a homogeneous material with some effective parameters (mass density, compressibility, speed of sound etc). There exists a vast literature about homogenization of periodic and random materials [82], but it has been shown that multiple scattering theory is a versatile method in the case of sonic crystals. For instance, for the simple case of a two-dimensional arrangement of fluid inclusions (parameters labeled with “a”) in a fluid background (parameters labeled “b”), the effective bulk modulus B_eff and mass density ρ_eff are given by

(1)1Beff= fBa+1−fBb, ρeffρb=ρa(Δ+f)+ρb(Δ−f)ρa(Δ−f)+ρb(Δ+f)

where f is the filling fraction of the crystal (area of the inclusions divided by the area of the unit cell) and Δ represents a modification of the filling fraction due to the multiple scattering processes [83], [84].

For airborne propagation any solid material will have a much higher density than that of air, and the extreme impedance mismatch blocks the penetration of the acoustic wave into the solid inclusions. Therefore, those solids can be assumed to be “rigid” inclusions, and their mass density and compressibility can be set as infinite. In this case, if the filling fraction is small and the multiple scattering interactions can be neglected (Δ can be regarded as 1) the effective mass density and velocity can be simplified as [83]

(2)ρeff=ρair(1+f)1−f,  ceff=cair/1+f

and we recover the effective mass density and effective sound velocity derived by Barryman [85] and the heuristic model [86], respectively.

As an example of application, a circular cluster of 151 wooden cylinders arranged in a hexagonal lattice embedded in air was homogenized by multiple scattering, and the resulting effective parameters as a function of the filling fraction are shown in the left panel of Figure 1. We can clearly observe, especially from the insets, that the simplified approach by Eq. (2) for effective mass density is valid when the filling fraction is no larger than 0.7; meanwhile, the maximum limit of filling fraction is 0.6 for effective sound velocity. Experiments verified the homogenization method, and it was shown that for wavelengths larger than 4 times the lattice constant, the cluster could be considered as a homogeneous cylinder with fluid-like properties [83].

Figure 1:

Effective elastic parameters of sonic crystals.

(A) Effective mass density and effective sound velocity for a circular cluster consists of 151 rigid cylinders embedded in air. The black dots stand for the results from the multiple scattering simulations in Eq. (1), and the dashed lines are values obtained from Eq. (2) [83]; (B) Impedance-velocity diagram of sonic crystals consist of square lattice of solids 1 and 2 cylinders in air. Insert shows the unit cell [87].

The theory can be generalized to more cylinders per unit cell. If we consider two types of cylinders, with materials labeled as 1 and 2, embedded in a square lattice, their fluid-like parameters being (ρ_a1, B_a1) and (ρ_a2, B_a2) and the filling fractions being f₁ and f₂, respectively, the effective parameters can be obtained as

(3)1Beff=1−fBb+f1Ba1+f2Ba2, ρeff=1+f1η1+f2η21−f1η1−f2η2ρb,  ceff=Beffρeff

where f=f₁+f₂, η_i=(ρ_i−ρ_b)/(ρ_i+ρ_b).

The right panel of Figure 1 shows a “phase diagram” of the effective impedance and effective sound velocity for several selected pairs of materials under the low filling fraction condition (neglecting multiple scattering interactions), for different pair of materials and backgrounds (Pb-Fe, Pb-Al and Fe-Al in water and aerogel-rigid cylinders in air). Any point on each curve is calculated sweeping over the total filling fraction. The overlapped corner point (1, 1) in this diagram corresponds to the case of f₁=0 and f₂=0. Notice that the effective impedance matches that of the background when the sonic crystal is made of aerogel and rigid cylinders embedded in air, making the sonic crystal transparent to air possessing a different refractive index. This interesting property was used in ref [17] to design a completely transparent GRIN lens.

Multiple scattering allows the homogenization of finite clusters; therefore, it is a versatile method to study the behavior of the effective parameters when disorder is introduced in the lattice [84] or even the effect of the cluster’s size, as it was done in Ref. [88], where it was found that small clusters with specific number of scatterers presented the same effective parameters as that of the infinite medium.

3 Homogenization of phononic crystal plates

In a homogeneous thin elastic plate, the propagation of flexural waves (antisymmetric Lamb mode wave) can be approximately described by the bi-Helmholtz equation (assuming time harmonic dependence of the field)

where ρ_b, h_b and Db=Ebhb3/12(1−vb2) are the mass density, the thickness and the rigidity of the plate, respectively, E_b being the Young’s modulus and v_b the Poisson’s ratio. W is the out of plane displacement field. Multiple scattering can also be applied for the homogenization of distribution of scatterers in elastic plates, and frequency-dependent effective parameters were obtained in [61]. The expressions derived there can also be applied in the low-frequency limit to obtain the following effective parameters

where the subscript “a” (“b”) means the parameters for inclusions (background). The effective phase velocity can be obtained by the ratio of the wavenumber between the effective medium and the background, giving

As an example, let us consider circular inclusions in a thin aluminum plate (h_b=0.1a) and let the inclusion’s material be chosen as a hole (empty inclusion), lead and rubber, whose parameters can be found in Ref. [61]. The effective parameters as given by Eqs. (10) and (11) as a function of the filling fraction are plotted in the left panel of Figure 2. The effective mass density is the volume average approach as for bulk waves. We see that the effective rigidity and the Poisson’s ratio for the case of rubber behaves like that of a hole, since the ratio of the Young’s modulus between the rubber and the aluminum is very low. For the effective phase velocity, the difference between the rubber and the hole is more evident. In the right panel of Figure 2 is shown the dispersion curves from finite element methods (red dots) and the effective parameters in Eq. (10) (blue lines) for a triangular lattice of hole/lead inclusions in a thin aluminum plate with the radius of inclusion being 0.3a. We can clearly observe that the agreement between the numerical simulation and the effective theory is quite good for the antisymmetric (A₀) and symmetric (S₀) Lamb modes. This is remarkable since we have developed the theory for the antisymmetric mode only; however, there is a connection between the speed of these two modes [26] if a wave propagates from plate h₁ to plate h₂, nA2=nSh1h2, which allow us to relate both, even if the propagation of the S mode has not been included in the theory.

Figure 2:

Effective elastic parameters of phononic crystal plates.

(A) The effective mass density (a), rigidity (b), Poisson’s ratio (c) and phase velocity (d) with respect to the filling fraction for different cylinder inclusions in an aluminum plate [61]; (B) Dispersion relations for a phononic crystal plate made of a triangular lattice of circular holes and lead inclusions with a radius of 0.3a. The red dots show the band structures obtained by the finite element method, while the blue lines show the dispersion curves calculated with homogenized parameters from Eq. (10) [26].

Actually, there are three fundamental plate modes as it is seen in the right panel of Figure 2, with the third one named shear-horizontal mode (SH₀), which needs an additional effective stiffness component c66eff to compute its dispersion relation. The theory developed for flexural waves allows the calculation of the velocity of the S₀ mode only due to the special relationship between the velocities of A₀ and S₀ modes; however, a similar relationship is not found, or at least obvious, for the SH mode. However, the homogenization theory for bulk phononic crystals presented in Section 2.2 offers the full components of the effective stiffness matrix, and it can be used to compute the effective velocity of shear waves in a phononic crystal plate. This connection is not obvious in principle; however, it can be explained by a two-step homogenization procedure.

Figure 3, upper left panel, shows a two-dimensional phononic crystal made of cylinders embedded in a bulk matrix. The phononic crystal plate at the upper right panel can be regarded as a finite “slice” taken from the bulk phononic crystal. In the low-frequency limit, the phononic crystal behaves as a homogeneous material with effective parameters and tetragonal anisotropy (square lattice) or transversal isotropy (triangular lattice), and this symmetry is maintained for the phononic crystal plate, as illustrated from the upper panel to the lower panel in Figure 3A. Thus, the effective phononic crystal plate can also be considered as a finite “slice” of the effective homogeneous material with the same symmetry. Such homogenization procedure has a non-trivial implication in accordance with symmetry.

Figure 3:

Homogenization of phononic crystal plates for all three fundamental Lamb modes.

(A) Illustration of the homogenization relationship between phononic crystal and phononic crystal plate [26]. (B) Dispersion comparison of a phononic crystal plate between the finite element method (black dots) and the homogenization theory with effective parameters from Eq. (8) (red lines) for the plate thickness as 0.1a (left), 0.5a (middle) and a (right). The phononic crystal plate consists of a triangular lattice of hole-gold shell unit cell whose inner radius is 0.2a and outer radius is 0.4a [25].

As an example, a phononic crystal plate with the unit cell shown in Figure 3B is considered for three different thicknesses, consisting of hole-gold shell unit cell in triangular lattice. The inner radius is 0.2a, while the outer radius is 0.4a. The dispersions obtained from homogenization method with effective parameters in Eq. (8) are compared to those from finite element methods, as displayed in Figure 3B. In the low-frequency limit, these effective phononic crystal plates behave like transversely isotropic plate. It is clear that the agreements between the finite element method and homogenization method are excellent for S₀ and SH₀ modes for all cases of thickness before the mode branches deviate, as their dispersions do not depend on the plate’s thickness. However, the dispersion of A₀ mode relates with the plate’s thickness. For thinner plates, the dispersion of fundamental A₀ mode keeps the parabolic shape throughout the Γ-J direction of the first irreducible Brillouin zone, so that there is an excellent agreement in the whole range of Γ-J for plate’s thickness as 0.1a. For thicker thickness such as 0.5a and a, the agreement for A₀ mode goes well before the A₀ mode branch deviates.

As a conclusion, the homogenization procedure for phononic crystal plates illustrated in Figure 3A works well for low-frequency limit and thin plates.

4 GRIN devices for bulk acoustic and elastic waves

This section is devoted to the review of GRIN devices for waves in bulk materials; therefore, no interfaces are assumed and surface or guided waves are excluded in this section (they will be considered in the next one). In principle we could divide the devices into two types, depending on the nature of the matrix which can be fluid or solid. This division makes sense since in the former shear waves does not exist in the background and the field is described entirely by a scalar pressure field, while in the latter shear waves with different propagation velocities have to be as well considered. However, for the fluid matrix an additional division can be made, since the physics involved in the scattering process is different if the background is air or water. While in air almost all solid materials are acoustically rigid due to the impedance mismatch between air and any solid, in underwater acoustics this impedance mismatch is of the same order, so that we have to consider the elastic nature of the inclusions. For this reason, we divide this section into three subsections, corresponding to the study and design of GRIN devices for air, water and solid backgrounds.

4.1 GRIN devices for bulk acoustic waves in air

In this section, we will discuss the flat GRIN lens for focusing whose surface is flat instead of being curved as in traditional lenses.

The index profile of GRIN flat lenses keeps constant along the propagating x axis while it changes along the vertical y axis as [18]

(12)n(y)=n0sech(αy),  α=1hcosh−1(n0nh)

where h is the half-height in y axis of the lens, n₀ is the refractive index at y=0, and n_h is the index at the upper or lower edges y=±h. Such index profile is designed with the aim of presenting low aberrations [18]. The position of focal point x_f at the right side of the lens can be derived from the differential equation satisfied by the energy rays [19].

As shown in Figure 4A, a GRIN lens was fabricated with nine columns of metal rods whose radii are varied locally by an inverse design as follows: first make a sweep of the filling fraction for sonic crystals made of rigid cylinders to plot the diagram of effective index and the filling fraction; then calculate the required index at a given position along the y axis of the GRIN lens; next, choose the corresponding filling fraction of the rigid cylinder from the diagram in the first step in accordance with the required index in the second step; finally, calculate the required radius from the selected filling fraction in the third step. Figure 4B shows the pressure distribution fields after the GRIN sonic crystal at a frequency corresponding to a wavelength of about 3.8a numerically (upper panel) and experimentally (lower panel). Similar focal spots and side lobes are evidently seen. To quantitatively characterize the focusing effect, the pressure along the x and y axes are also plotted in Figure 4C, showing very good agreements. In this approach with rigid cylinders in GRIN sonic lenses, the effective impedance of the lens is not matched with the background, so that the maximum pressure at the focal spot is limited. As indicated in Figure 1B in Section 2.1, the mixed lattice of aerogel and rigid cylinders can also tailor the effective index while keeping the effective impedance matched with the background, so that zero reflectance can be achieved to make the lens transparent. The first GRIN sonic Wood lens proposed in Ref. [87] has a higher pressure focalization at the focal spot.

Figure 4:

Sonic GRIN flat lens.

(A) Schematic view of a GRIN lens with the ray trajectory. For the impinging plane sound wave, the rays will bend and focus at x_f position at the right side of the lens. (B) Normalized sound amplitude at 4.5 kHz mapped at the right side of the GRIN lens from the multiple scattering method (upper panel) and the experiment (lower panel). (C) Normalized sound amplitude along the x and y axes crossing the focusing point. Blue lines stand for simulation, and black dots stand for experiment [19].

4.2 GRIN devices for bulk acoustic waves in water

Different from the impedance mismatch between the air and solid scatters, the mass density and acoustic velocity are in the same order for most scatterers and water, so that the scatters cannot be regarded as rigid anymore. From Eq. (1) in Section 2.1, one can also calculate the effective parameters of phononic crystals in water.

For the same type of GRIN flat lens as in the last section, similar phononic crystals made of radius-varied steel cylinders are proposed in water [93]. In the upper left panel in Figure 5, the profile of radius is calculated to meet the required index profile like Eq. (12), the maximum index being at the center and minimum at the edges in vertical axis. A plane wave with center frequency at 20 kHz is excited, corresponding to wavelength as 4 times lattice constant around the cutoff frequency in the low-frequency limit. A clear focusing spot is observed from the experimental pressure intensity field. From the multiple scattering theory(MST) simulation, a similar position focal spot is also obtained as shown the dot “+” in subpart (c) of upper right panel of Figure 5. It is worthy to note that this GRIN lens also supports the source off the central axis, which was demonstrated numerically and experimentally [93].

Figure 5:

GRIN flat lens with phase delay and advance approaches.

Upper left panel: Schematic view of GRIN flat lens with radius-varied steel cylinders embedded in water. Upper right panel: Normalized pressure amplitude from the source (a), measured normalized pressure intensity after the GRIN lens (b) and the same pressure intensity map using MST (c) [93]. Lower left panel: Schematic view of GRIN flat lens made of the same hollow air-filled aluminum tubes with anisotropic lattice spacing. Lower right panel: Sound pressure intensity at 13.05 kHz from 2D MST (a), experiment (b) and 3D Rayleigh-Sommerfeld approximation (c) [20].

Equation (12) gives a picture of special index profile for focusing with low aberration, while it does not define that the maximum index at the center n₀ should be larger than that of the background n_b; in another word, Eq. (12) can also be achieved by phase advance approach. In the lower panel of Figure 5, phononic crystals made of air-filled aluminum tubes are arranged in water with anisotropic lattice spacing. The thickness of the aluminum shell is 1/10 of the radius that is designed to have impedance matching with respect to water [94], a higher velocity but lower mass density acting as metafluids. The required filling fraction profile to meet the index distribution is implemented by using anisotropic lattice spacing [95]. From the lower left panel of Figure 5, the lattice space in vertical axis is large at the center while it reaches the minimum at the upper and lower edges. The minimum velocity at the center is still larger than that of the background water, resulting in a phase advance GRIN lens. A good comparison is found in the lower right panel of Figure 5 among the experimental measurement and two predicted numerical simulations. It is worth mentioning that such GRIN flat lens is broadband, which can be implemented in water from experimental point of view as the short pulse from the transducer in water is normally broadband.

4.3 GRIN devices for bulk elastic waves

The effective parameters for bulk elastic wave in 2D phononic crystals can be found in Eq. (8) in Section 2.2. Here, another approach based on the analysis of the lowest band among the dispersion curves is introduced. Let us consider a phononic crystal made of square lattice solid inclusions embedded in an elastic matrix. For small anisotropic ratios, the effective index for a bulk elastic wave can be obtained as [18]

(13)neff=nΓX+nΓM2,  nΓX=cbcΓX=cbdω/dkΓX,  nΓM=cbcΓM=cbdω/dkΓM

where c_b is the acoustic velocity in the background matrix, c is the group velocity of the bulk elastic wave mode, and ΓX and ΓM are two orientations in the first Brillouin zone. For an epoxy matrix phononic crystal, the diagram of the first band of shear-vertical (SV) mode is plotted in two ways: (i) steel inclusions with different filling fractions (different radii) as shown in the left panel of Figure 6A; and (ii) fixed filling fraction but with different material inclusions as shown in the right panel of Figure 6A. The elastic parameters for different solid materials can be found in Ref. [18]. Both approaches can change the slope of the first SV mode leading to a variation of effective index by Eq. (13).

Figure 6:

GRIN phononic crystals for bulk elastic waves.

(A) Left: The first shear-vertical (SV) mode of bulk elastic wave along the ΓX for varying the filling fractions of steel cylinders embedded in epoxy. The arrow indicates the increasing of the filling fraction. Right: Same band diagram but for phononic crystals made of different solid material inclusions in epoxy with a fixed filling fraction as 0.126. The arrow indicates the decreasing of stiffness. (B) The finite difference time domain (FDTD) simulated SV wave propagations in a phononic crystal whose inclusions are made of different materials embedded in epoxy at two different frequencies. The calculated and simulated focal lengths as a function of frequency are also plotted as solid line and dotted line, respectively [18].

From the dispersion curve variations, for a given band the effective index will change locally along the reduced frequency, and for different bands the dispersive properties of the effective index are slightly different. Therefore, if a GRIN flat lens is designed based on a given frequency, it may have a small deviation in focal length as the effective or fitted gradient coefficient α in Eq. (12) changes, as seen in the right panel of Figure 6.

5 GRIN devices for flexural waves

The control of elastic waves in plates provides the most interesting applications such as on-chip acoustic propagation, thermal conductivity control and opto-mechanic interactions in micro/nanoscale for GRIN devices, which can be attractive for the communities of physics, materials and nano science. For 2D thin plates, or membranes, the elastic waves propagate as surface waves due to the constraints of the plate’s two boundaries, resulting in antisymmetric, symmetric and shear-horizontal polarizations. Among those modes, the fundamental antisymmetric mode, in another word flexural mode, is widely studied since it is mainly characterized by out-of-plane displacement component in plates that can be easily detected by means of optical methods.

Figure 7A shows the foci of the GRIN flat lens connected to a linear phononic crystal waveguide for flexural waves in a piezoelectric plate based on deep reactive ion etching. The well-designed GRIN lens focuses a plane incident wave as a spot at the interface between the GRIN lens and the phononic waveguide so that waveguiding is demonstrated. This would have potential for developing active microplate lenses [96]. In Figure 7B, GRIN flat lens is fabricated with silicon pillars erected on top to behave as a metalens to demonstrate the focusing behind the lens in the near-field beyond the diffraction limit. The dipolar resonant pillars and the flexural waves can exhibit polarization coherency to enhance the evanescent waves in order to help the focal spot include better information [68]. The conception of mixed local resonance and GRIN lens can also be applied to other types of waves. At nanoscale, a gradient optomechanical/phoxonic crystal in one-dimensional semiconductor slab is fabricated with the minimum radius of the hole as only 100 nm [78], as displayed in Figure 7C. This phoxonic crystal demonstrates a full phononic band gap at the gigahertz domain. The optomechanical coupling originates from moving interfaces and the photoelastic effect. The acoustic mode inside the band gap exhibits high mechanical quality factor making the optomechanical phonon coherent manipulation possible. Figure 7D illustrates that the confinement of elastic waves in the GRIN devices can be used for energy harvesting [71], [72]. Piezoelectric energy harvester disks are deposited on a GRIN flat lens and a Luneburg lens at a position where the incident plane wave is focused as a spot. Comparing to the energy harvested in the background plate, the GRIN devices can generate output electrical power an order of magnitude higher.

Figure 7:

GRIN devices for flexural waves.

(A) Acoustic waveguide by combining GRIN flat lens and phononic crystals in an 80 μm thick plate. SEM image shows samples in the rectangular dotted box with the radii of gradient holes varying from 15 μm to 25 μm [96]. (B) Acoustic metalens by combining GRIN flat lens and pillars (30 μm radius) to exhibit the near-field subdiffraction focusing [68]. (C) 1D phoxonic crystal with varied unit cells (the minimum radius of holes is 100 nm) to exhibit the optomechanical interactions [78]. (D) GRIN flat rectangular lens and Luneburg circular lens for energy harvesting [71], [72].

In Section 3, Eq. (10) gives a useful tool to calculate the effective parameters of phononic crystal plates whose thickness is fixed as a constant. Therefore, GRIN devices can be designed by varying the filling ratio, in other words changing the radii of holes [62]. Now let us revisit the equation of motion (9) which describes flexural waves using Kirchhoff-Love approximation. When a plane wave propagates in a thin plate with a wavenumber k, Eq. (9) gives the solutions for the dispersion relation and phase velocity as

from where one can see that the phase velocity is not only related to the elastic parameters but also depends on the plate’s thickness. This special property of flexural waves allows designing of GRIN devices just by locally varying the plate’s thickness without any phononic crystal with holes. The refractive index for flexural wave can be easily derived when all elastic parameters remain fixed as in [60]

In the upper left panel of Figure 8, a set of omnidirectional GRIN devices are displayed with their refractive index profiles, namely, Luneburg, Maxwell fish-eye, 90° rotating, Eaton and concentrator lenses [97], [98], as shown by the different colors. Implementing the refractive formula Eq. (15), their corresponding thickness variations are also shown in the panel below. To achieve higher refractive index n_eff, it needs smaller thickness h_eff. For 90° rotating, Eaton and concentrator lenses, the maximum index at the center of the circular lens is larger than 5, so that the relative thickness there is smaller than 1/25. The designed GRIN lenses based on thickness variation are performed with full wave simulations, and their corresponding wave functionalities are displayed at the bottom of Figure 8 [24].

Figure 8:

GRIN lenses for flexural waves based on plate’s thickness variations.

(A) Refractive index n profiles and corresponding plate’s thickness h/h_b variations for Luneburg, Maxwell fish-eye, 90° rotating, Eaton, concentrator lenses. (B) Wave performances of those lenses [24]. (C) Omnidirectional invisible lens for flexural waves based on thickness variation for two different frequencies [64].

Invisibility is always a challenging problem for wave control in physics. Acoustic cloaks in principle require inhomogeneous and anisotropic materials for their realization, which are quite complicated to fabricate in practice [99], [100], [101]. A special gradient refractive index profile for omnidirectional invisible lens is derived that requires an infinite index in the vicinity of the center [102], [103]. However, such singularity can be achievable for flexural wave by Eq. (15). Nevertheless, from a practical point of view, the central index can be chosen with a proper finite value in accordance with the step in design. A ray enters the lens and exits in the same direction by flowing a loop around the center inside the lens, thus producing an invisible effect. The multi-frequency performance of the invisible lens is simulated with cylindrical waves from a point source as seen in Figure 8C. [64]. It is worthy to note that acoustic cloak for flexural wave based on nonlinear transformation acoustics method is proposed and realized also by varying the plate’s thickness locally to fulfill the required rigidity profile [104].

6 Multimodal GRIN devices

From Figure 3B, the fundamental Lamb modes have another two polarizations in addition to flexural mode. In fact, excitation in the thin plate from a given source will generate all kinds of modes, which means a design for only one mode may hinder the full functionalities of a GRIN device for applications such as absorption and energy harvesting. As discussed in the last section, the flexural wave can be independently controlled by thickness variation for given elastic parameters. The S₀ and SH₀ modes mainly depend on the elastic parameters (c₁₁, c₁₃, c₃₃) and (c₆₆), respectively, for transversely isotropic materials. Therefore, it needs a homogenization theory to obtain the full effective stiffness matrix. As illustrated in Figure 3 in Section 3, the full homogenization method described in Section 2.2 for bulk phononic crystal can be also applied to the thin phononic crystal plate by considering the latter as a finite thin “slice” taken from the former. All the three polarizations of the Lamb waves are well described in the low-frequency limit. Taking the same gold-shell-hole unit cell, we make a sweep for the outer radius while also varying the inner radius for a given outer radius. Figure 9A shows a full phase diagram of all possible pairs of (n_S0, n_SH0) domain with this sweeping approach, which exhibits several interesting properties: (i) the dotted red line cutting the blue domain indicates that the effective indexes for the S₀ and SH₀ modes are identical, supporting the simultaneous control of these two modes; (ii) the maximum effective index for S₀ mode can reach up to almost 5, while it can be larger than 10 for SH₀ mode; (iii) the wide blue domain shows that it is possible to design a special (n_S0, n_SH0) trajectory that makes a GRIN device work as one type of lens for S₀ mode and as another lens for SH₀ mode. After defining the index profiles for SH₀ and S₀ modes, the index of A₀ mode can be designed by means of an independent thickness variation and such that it works similarly to the one or both of the SH₀ and S₀ modes or as another totally different lens. Thus, a full control of Lamb modes and advanced multimodal GRIN device becomes achievable.

Figure 9:

Multimodal GRIN devices to fully control all three fundamental Lamb waves.

(A) Diagram of all possible values for the effective indexes for (n_S0, n_SH0) by varying the inner and outer radii of the gold-shell-hole structure. (B) GRIN Maxwell fish-eye lens designed for all three Lamb modes simultaneously. (C) An advanced multimodal GRIN device works as Luneburg lens for the S₀ and SH₀ modes while as Maxwell fish-eye lens for the flexural mode [25].

In Figure 9B, the red dotted circle shows the area of the Maxwell fish-eye lens that is designed to work similarly for all the three fundamental modes. A point source is excited at the left border of the lens, and a focused point image is formed at the right diametric border for a given frequency. It should be mentioned that the wavelengths of the three modes are different since their speeds are not the same due to different dispersion properties. The design of the full control means the positions of point focusing are the same, but not the wave field patterns. Finally, an advanced multimodal GRIN device is displayed in Figure 9C: it works as Luneburg lens for the S₀ and SH₀ modes, while it behaves as Maxwell fish-eye lens for the A₀ mode. The full control method is also extended to GRIN flat lens and beam splitter [63].

7 Summary

Over the past two decades, we witnessed the rapid growth and development of phononic crystals and metamaterials as well as the recent emergence of metasurfaces to deepen the understanding of wave physics. Although the effective index is a homogenized result from infinite periodic structures, it is still proper to define a local index in a gradient structure, being a bridge between infinity and locality. The homogenization theories of phononic crystals/metamaterials and the design principles of GRIN devices are implemented in terms of wavelength, mainly in long wavelength limit λ>3~4a, a being the lattice constant. Therefore, the GRIN conception can be applied to the structure’s size from nanometer to meter and different wave natures propagating in solids, liquids and gases with free or constrained boundary conditions, among which the applications for elastic waves in thin plates at nanoscales show a promising potential in future for wireless telecommunications, heat conductivity, phonon sensor and so on. The complex properties of elastic waves in thin plates are also fully demonstrated to design more advanced and multimodal GRIN devices. The applicationsof GRIN phononic crystals and metamaterials are mainly related with energy redirection and focusing; however, more advanced devices could also be considered some subtype of GRIN devices, though theirfunctionality is not designed in the framework of “ray trajectory design”. This is the case of acoustic cloaks [100], [105], [106], [107], [108], parity-time-symmetric acoustic structures [109], [110], [111], one-way acoustic diode [112], [113] and acoustic rainbow trapping [114], [115], [116]. Nevertheless, several challenges remain to be answered in future: (1) When wave propagates from one medium to another one, from a uniform surface point of view its behavior is governed by the Snell’s law that involves the ratio of the refractive index between the two media. However, the wave energy transfer rate is governed by the impedance ratio between the two media. For some GRIN devices such as flat lens and transmitted-type metasurface, the impedance matching condition is not always easily fulfilled, which will limit the focusing energy level and the corresponding energy harvesting efficiency. (2) New refractive index profiles are needed to discover mathematical methods that can help to design updated wave functionalities. (3) New physics such as thermal conductivity behaviors and nonlinear optical/acoustic responses in gradient phoxonic crystals at nanoscale. It can be foreseen that new breakthroughs will be reported by the united efforts in the communities of physics, material science, engineering, nano science and mathematics among others.