Dynamics of Phononic Materials and Structures: Historical Origins, Recent Progress, and Future Outlook

A periodic medium is a material or structural system that exhibits some form of spatial periodicity, which can be in the constituent material phases, or the internal geometry, or the boundary conditions. The study of periodic materials and structures has a long history in the field of vibrations and acoustics with origins that can be traced back to Newton's first attempt to describe the propagation of sound in air [1] and Rayleigh's early study of continuous periodic structures [2]. The topic has grown to become of fundamental importance in condensed matter physics due to the role that atomic vibrations (and electronic structure) play in determining the properties of crystals. The notion of a “phonon” has emerged in the context of vibrations in a crystal lattice. Formally defined as a quantum of vibrational energy in an elastic medium—which may be interpreted as a discrete particle-like quantity of sound in a solid—the term has also been associated with classical wave-like vibrations and acoustics, mainly in the context of periodic media. Subsequently, it has grown to become customary to refer to a periodic material, or structure, even if in the size scale of a large engineering system, as a “phononic material,” or a “phononic structure.” As such, there has been an exquisite connection, especially in the past few decades, between phonon physics and dynamics of periodic materials and structures. Perhaps what stands out the most is the concept of a band diagram—a diagram that represents the relationship between frequency (or energy) and wavenumber, along multiple directions. In physics, this diagram represents the backbone of electronic structure theory, credited for forming a basis for the classification of all crystals into metals, semiconductors, and insulators. In mechanics, a band diagram is precisely a representation of the dispersion relation describing the nature of free wave propagation in an elastic (or acoustic) medium.

Looking closely at the study of periodic systems in the past half-century, we find that researchers in vibrations and acoustics, and more broadly from the mechanics community at large, have conducted a considerable amount of work on key theoretical foundations, concepts, and analysis techniques that are relevant to periodic systems in other nonmechanics disciplines. Arguably the two most motivating applications in mechanics, going back to the 1950s, and extending through the 1990s and beyond, have been composite materials (which conveniently have been modeled as periodic materials) [3,4] and aircraft structures (which naturally exhibit some degree of periodicity emanating from the presence of ribs introduced primarily for strengthening and stiffening purposes) [5,6]. Other areas related to mechanical and civil engineering include multiblade turbines [7–9], impact resistant foam/cellular materials [10–12], periodic foundations for buildings [13,14], and multistory buildings and multispan bridges [15].

In recent years, starting in the early 1990s, the field has experienced a resurgence with the introduction of phononic crystals [16–21], and, roughly a decade later, with the study of acoustic metamaterials [22], which may be also broadly considered as examples of phononic materials. A phononic crystal is a composite or nonuniform material consisting of one, two, or more material phases (solid and/or fluid) arranged periodically in space. This is not much different from periodic composite materials studied earlier in the engineering literature, with the only difference that more rigor is applied in the treatment of symmetry by borrowing concepts from crystallography [23]. The already existing connection between periodic materials and phonon physics has received a substantial upgrade with the introduction of phononic crystals. An acoustic (or elastic) metamaterial, on the other hand, is also (usually) a periodic material except it has the added feature of exhibiting local resonance, or possibly other traits, that can give rise to “unusual” dynamical behavior, such as negative refraction and negative group velocity. Metamaterials also have a foundation in electromagnetics, and the disciplines of acoustics and elastodynamics have been steering a parallel path in the development of this highly promising class of materials.

While applied research in phononic materials has been abundant in recent years, inspection from a fundamental mechanics perspective, and particularly from the standpoint of dynamical systems, is needed to advance the field in new directions. For example, techniques already developed for the treatment of damping and nonlinearities in structural dynamics provide a rich resource for application to the problem of wave propagation in periodic materials and structures. Similarly, numerical and experimental approaches originally developed for the characterization of conventional materials and structures provide indispensable tools that may be employed toward better understanding and exploitation of phononic systems.

The goal of this review article is to bring the topic of phononic materials to the attention of the Journal of Vibration and Acoustics and the Applied Mechanics Reviews communities, and to encourage interest in the field among researchers and practitioners from across all related disciplines. To this end, the article begins with an overview of historical developments underlying the early years of research in periodic materials and structures, and follows with an in-depth literature review of progress that came about with the advent of phononic crystals and acoustic/elastic metamaterials with special consideration given to aspects pertaining to dynamics, vibrations, and acoustics. A detailed technical review is then presented on fundamental concepts in phononic materials theory, numerical techniques, and the treatment of damping and nonlinearity, as well as experimental research encompassing the field. Finally, a brief outlook is projected toward the future on the basis of the progress achieved to date.

According to Ref. [24], the first work on a periodic one-dimensional (1D) lattice is the result of Newton's attempt to derive the formula for the velocity of sound in air. One can speculate that because of lack of knowledge of differential calculus, a natural approach would have consisted in discretizing the continuous media supporting the propagation of sound into a series of lumped masses connected by a lumped spring. This process, which is illustrated in Fig. 1, leads to a periodic system, which is known to have a more complex wave behavior than the original continuous system.

The same system was later investigated by John Bernoulli and his son Daniel, who, in a series of studies starting in 1727, demonstrated that a system of N masses is characterized by N modes of vibration and associated frequencies, and essentially formulated the principle of superposition. Subsequent studies included the estimation of the velocity of wave propagation along one axis of a cubic lattice structure as a function of wavelength. These efforts have led to the first instance of observation of wave dispersion from the estimation of the phase velocity's dependency on frequency. These results were later expanded and explained in 1881 by Lord Kelvin [25] among others, and the reader is referred to the classic books by Brillouin [24] and Kittel and Kittel [23] for a detailed historical perspective on the topic.

Periodic layered composites are the motivation for one of the first studies, as in Ref. [26], where a 1D periodic arrangement of different materials is investigated. The study of this class of assemblies found much interest in the 1960s and 1970s with the advent of composite materials, and the necessity to assess their integrity through nondestructive evaluation (NDE) methods based on the propagation of elastic waves. This has motivated studies of layered composites and the development of equivalent, or continuum models, as in the case of Refs. [3,27]. A continuum theory was developed therein to describe the dynamic behavior of a laminated composite by expressing the displacements of the reinforcing and the matrix layers as two-term expansions of the layers' midplanes. The accuracy of the theory was estimated by evaluating the dispersion curves for waves with polarizations normal and parallel to the layering and comparing them with exact curves. Variational formulations for the analysis of harmonic wave propagation in periodic elastic composites were presented by Nemat-Nasser in Ref. [4], which is the first paper of a series presenting studies on layered configurations of increasing complexities [28–32]. An early review of variational methods for the analysis of periodic composites can be found in Ref. [33], while continuum theories have also been proposed by Nayfeh and coworkers in the same period [34]. Their study of wave propagation in layered composites, again mostly motivated by NDE applications, has continued through the years, and notable contributions can be found in Refs. [35,36] along with studies on guided waves in layered composites as in the case of Refs. [37–39]. In several studies, the transfer matrix approach and other, similar matrix-based calculations are applied to investigate dispersion properties and the propagation of harmonic waves through the thickness of composites of general in-plane anisotropy. Mixture theories for the estimation of equivalent stiffness and mass properties of layered structures have also been pursued by other authors; see, for example, Refs. [40,41]. The ability of equivalent, mixture, or continuous theories to predict the dispersion properties of heterogeneous media is still the object of intense investigation, and has been addressed in some of the multiscale studies performed in recent years. Examples are the studies by Hulbert and coworkers [42–44], and by Fish and others [45–47] who have considered the 1D periodic waveguide as a test case for the initial assessment of the developed approximation techniques based on multiple scale expansions. Periodic configurations of highly mismatched materials have formed the basis for most of the periodic systems considered by the community, and are still a common trait in the design of “phononic crystals.” One of the early investigations in this particular area is found in Ref. [18], where two-dimensional (2D) periodic materials and their complete bandgap properties are discussed as a parallel to what was previously observed in “photonic crystals” (see Sec. 1.3).

In the area of structural engineering, Cremer and Leilich [48] were among the first to investigate harmonic flexural wave motion along a 1D periodic beam. The periodicity was introduced by either simple supports or by point masses located at regular intervals. The periodic beam is an example of a monocoupled periodic system whereby each unit cell communicates with its neighbors through a single degree of freedom. Its wave characteristics at any frequency are therefore described by a single wave mode and by a single pair of equal and opposite propagation constants. Other examples of earlier studies are reported in Ref. [49] and in the works of Heckl [50,51] and Ungar [52]. Therein, investigations were conducted on beam-plate systems, 1D beam structures with various periodic features and 2D beam grillages. Most of these systems and their propagation characteristics can be conveniently investigated through the formulation of a transfer matrix. Its repeated application through the number of cells composing a finite periodic system leads to the estimation of natural frequencies and mode shapes. In addition, the analysis of a single unit allows the estimation of the dispersion relations according to a procedure briefly summarized later in this paper (see Sec. 2.2). Among the pioneering efforts on the formulation of the transfer matrix, the work of Lin and McDaniel [53] is worthy of note, along with the work of Faulkner and Hong [54], which provides examples that can be useful to researchers who want to begin investigations in this area. Of interest is also the work of McDaniel [55,56], Thomas [57], and Williams [58] on rotationally periodic structures, which is an area that appears to have developed on an independent path, mostly as a result of its strong practical relevance and associated specific technological problems. Prime examples of rotationally periodic structures are bladed disks in turbomachinery, where several characteristic dynamic properties of periodic structures have been observed and linked to large responses and localization of vibrations, and have been identified as causes of failure [59,60].

Many noted contributions in this area are attributed to Mead and his coworkers at the University of Southampton, which has been a hotbed of activities starting from the 1960s. An excellent review of their work has been produced by Mead himself in Ref. [61], and will not be repeated here for the sake of brevity. Among the Southampton studies, we would like to mention the work of Gupta [62,63], which relates the modal content of beams, stiffened panels, and ribbed plates to their dispersion properties. The investigation of modal properties of monocoupled and multicoupled periodic systems can also be found in Refs. [64,65]. Reference [64], in particular, relates the bounding frequencies of the propagation zones of monocoupled systems to the natural frequencies of the individual unit cells and of a finite assembly. The work is extended to multicoupled systems in Ref. [65], where the effect of damping is also considered. The forced response of periodic systems is studied in Refs. [66,67], where a response function is found for an infinite, uniform, 1D structure subjected to an array of harmonic forces or moments. The case of random convecting loading is addressed in Refs. [68,69]. The application of numerical techniques for the prediction of dispersion and overall dynamic response is presented in Ref. [70], where the finite element (FE) implementation of wave periodicity conditions as defined by Bloch's theorem is illustrated and applied to simply supported beams and to complex built-up structures. The formulation of a real-valued eigenvalue problem for dispersion estimation, which is suitable for its solution through built-in modal solvers in commercial FE packages, was later presented by Thomas [57] and subsequently by Åberg and Gudmundson [71]. The application of the FE method for the analysis of periodic structures is an area of continued interest, as indicated, for example, by the work of Mace and coworkers [72,73]. In particular, the periodicity of the structure, or the consideration of a uniform structure as the limiting case of a periodic assembly, can lead to significant advantages in terms of computational cost, or to the ability to exploit commercial FE codes for the derivation of the dispersion properties of complex structural components, such as stiffened plates, sandwich panels, or cylindrical shells [73]. These approaches have been successfully extended to the study of 2D structures, as in Ref. [74]. The analysis of 2D periodic structural components, with a focus on their directional properties in addition to the stop and pass behavior, traces back to the early studies of Langley (see, for example, Refs. [75–77]) and continues to draw interest in regards to the influence of topology and geometry on the wave propagation characteristics of structural lattices [78–83]. The presence of forbidden zones of response, the observation of “caustics” [84], and the relation between wave directionality and the equivalent anisotropic elastic properties of the lattice are of relevance to the design of impact absorbing structures, acoustic waveguides, vibration isolation devices, and solids with tailored mechanical properties.

As mentioned earlier, a remarkable trait of phononic materials/structures research is the broad interest across disciplines. While research in composite materials and periodic structures has evolved among mechanicians the notion of artificial periodicity has also been appealing to the electromagnetics and photonics communities (see Ref. [85] for a review of historical developments in these areas). A key milestone occurred in 1987 when Yablonovitch [86] and John [87], almost simultaneously, proposed the concept of what was later referred to as a photonic crystal—a multiphase dielectric material ordered with symmetries resembling those of atomic-scale crystalline materials (as dictated by the rules of crystallography). Only a few years later, the same concept emerged for acoustic or elastic waves, which was subsequently called a “sonic crystal” or phononic crystal [16–21]. The public had a luminous introduction to the concept when a large periodic array of cylinders was put on display in Madrid to represent a sound sculpture [88].

In 1992, Sigalas and Economou [16] considered in-plane (mixed longitudinal-transverse) elastic wave propagation in a 2D phononic crystal, albeit composed of fictitious materials. In the following year, they repeated their study using gold cylinders embedded in a beryllium host and incorporated out-of-plane shear waves [17]. Kushwaha et al. [18,19], on the other hand, studied out-of-plane waves in a periodic array of nickel alloy cylinders in an aluminum alloy matrix. An analysis of in-plane wave propagation in a 2D phononic crystal consisting of carbon cylinders as inclusions and an epoxy resin as a matrix material was conducted around the same time by Vasseur et al. [20]. The problem of a three-dimensional (3D) phononic crystal consisting of solid sphere inclusions in a solid host was also considered by Economou and Sigalas in 1994 [89], and later by Psarobas and coworkers [90] and Sigalas and Garcia [91]. Meanwhile, there has also been significant interest in sonic crystals consisting of multiple phases of liquids and/or gases. The acoustic band structure of these systems were first reported by Sigalas and Economou [17], followed by other studies, such as the theoretical investigations by Kushwaha and Halevi [92,93] and the experimental investigations by Sánchez-Pérez et al. [94] (audible range) and Montero de Espinosa [95] (ultrasonic range). Regardless of the phase state and choice of materials, two distinct features were identified as critical for the design of phononic crystals: the unit-cell topology and the lattice symmetry. It was shown that for elastic waves a “cermet topology,” in which a stiff/dense material forms the matrix and a compliant/light material is used as inclusions, generally produces the largest bandgaps [96]. More recently, this concept was taken further by virtue of full-fledged unit-cell topology optimization [97,98]. As for lattice symmetry, it was revealed that phononic crystals patterned over a hexagonal lattice typically display larger bandgaps compared with those based on square lattices [99]. The specific choice of the constitutive material phases and the possibility of interweaving more than one lattice also provide means for improving bandgap characteristics [100].

In the late 1990s, the concept of phononic crystals was extended to free surfaces involving 2D semi-infinite domains consisting of a periodic array of inclusions extending indefinitely into the depth in the direction normal to the surface [101–104]. In this configuration, Rayleigh-type Bloch waves propagate along the surface. Among numerous investigations, the case of anisotropic material properties was considered in this context [101,104]. Similar to free surfaces, phononic crystals were also studied and realized in the form of plates of finite thickness [105–117]. The interest in this class of problems has been mostly in the study of Rayleigh–Lamb waves in plates with cylindrical inclusions composed of either a second solid phase [106,107,109] or a fluid, such as air or even simply a vacuum [112–114]. The possibility of using spheres as inclusions in plate-based phononic crystals has been studied as well [108]. Contribution of phononic crystals to silicon MEMS [118] and standard MEMS [119] has emerged from this line of research. Further analysis and physical realization of phononic crystal plates with unit-cell sizes at the micrometer scale created a platform for radio-frequency sensing applications [112–114]. For additional detail, the reader is referred to the review article by Olsson and El-Kady [115], and to other recent reviews on phononic crystal plates, or 2D phononic crystals in general, which can be found in Refs. [116,117].

As research in phononic crystals continued to grow, other configurations were proposed, for example, by Torrent and Sánchez-Dehesa, who considered a radial phononic crystal for which a frequency band diagram was calculated for wave motion in radial coordinates [120]. More subtle symmetries have also been considered, such as quasicrystal [121] and fractal [122,123] phononic crystals—in analogy to similar studies concerned with photonic crystals. Furthermore, the field also witnessed advances in the type of theoretical investigations carried out. For example, Laude et al. [124] and Romero-García et al. [125] devised techniques to compute evanescent Bloch wave modes associated with phononic crystals alongside the usual propagating modes. The effects of periodicity truncation rendering phononic crystals finite in their spatial extent have also been examined, e.g., see [126–129], which builds on earlier studies on the topic by Gupta et al. [62,63].

The realization of phononic crystals in different forms, e.g., bulk, semi-infinite surfaces, plates, beams, has opened up a range of opportunities for realistic applications in engineering and applied physics. The concept of a bandgap naturally lends itself to applications involving vibration and acoustic filtering and control. For example, vibration isolation [130,131] or minimization [128,132], as well as noise isolation and control [94,133], may be realized with profound effectiveness. Similarly, nondestructive testing of materials for defects may be conducted by utilizing phononic crystal phase-surface characteristics for the emission of probing “wave beams” [80].

One of the remarkable functions that was realized using phononic crystal is waveguiding and localization by means of defect-engineering [134–139]. This has been achieved by introducing line or point “defects” into an otherwise perfect phononic crystal. The defect is incorporated by simply removing one or more inclusions (or other features that the characteristic unit cell possesses) in a manner as to introduce an anomaly in the form of a given pattern. In the defected region, the system will enable mode localization within the phononic crystal's bandgap frequencies, and this in turn allows for waveguiding or localization in a cavity. The latter phenomenon provides the basis for signal sensing in the microscale phononic crystals mentioned above [112–114]. These systems hold promise for radio-frequency communications due to the practicality of direct integration over CMOS circuitry with low insertion losses and their ability to provide frequency selection with high quality factors. In addition, the principle of defect-based elastic wave confinement has been utilized toward the realization of wave demultiplexers—passive devices that enable waves at certain frequencies to travel in different directions [140–142]. Further development along this track has started to produce early realizations of digital logic circuits based on phononic crystals. For example, phase differences associated with lack of collinearity between the wavevector and the group velocity (as revealed by equifrequency contour diagrams) have been used toward the creation of basic logic gates [143,144].

Utilization as sensors for the determination of liquid properties has been yet another application for phononic crystals [145,146]. The underlying concept is based on introducing a line defect that acts as a slit cavity, which in itself is also part of the fluidic system. The cavity mode causes a distinct transmission peak within the bandgap of the phononic crystal where the frequency of maximum transmission depends on the speed of sound of the confined liquid. In other work involving microfluidics, Cooper et al. [147] have demonstrated a simple interface between a piezoelectric surface acoustic wave device and a disposable microfluidic chip, patterned with a phononic crystal structure. They have illustrated that microfluidic manipulation, such as centrifugation of blood, may be performed on a disposable phononic chip. Employment of phononic crystals for mass sensing has also been demonstrated, e.g., Nardi et al. used a hypersonic surface phononic crystal to realize ultrahigh sensitivity mass sensing at operational acoustic frequencies that exceed 100 GHz [148].

Another active line of applied research is concerned with wave collimation and refraction. The interest here has been primarily in acoustics, although the ideas are also applicable to elastic waves. The function of sound collimation has been studied by Chen et al. [149], Espinosa et al. [150], Christensen et al. [151], and Shi et al. [152]. Refraction and focusing of sound waves have been a major topic due to the potential significant impact on the fields of subsurface imaging and NDE. Phononic crystals can be used to induce two types of refractions: positive and negative. Positive refraction is carried out by shaping the phononic crystal boundaries in order to refract an incoming wave [153–155]. The realization of negative refraction on the other hand is more subtle as it is based on tuning the phononic crystal to exhibit opposite phase and group velocity signs at the frequency of interest [156–163]. Use of phononic crystals for focusing has also been realized using nonlinear solitary waves generated in an array of precompressed metallic beads [164].

Another intriguing application of phononic crystals is wave rectification, whereby a device is designed to allow acoustic or elastic waves to travel along one direction but not in the opposite direction. This principle promises to be of great benefit for vibration mitigation but also could be used as a building block in acoustic/elastic logic gates in analogy to electronics. While not yet perfectly implemented, some promising approaches have been proposed. For example, Liang and coworkers utilized 1D layered phononic crystals composed of two segments, a linear segment to filter an incoming wave and a nonlinear segment to enable energy to transfer to a different mode at a slightly different frequency [165,166]. Nonlinearity was used once again, now in the context of bifurcation phenomenon within a linear chain of precompressed beads, to generate a rectification effect [167]. Two-dimensional linear phononic crystals were employed as well to realize an acoustic diode, but with a mechanism based on diffraction and modal conversion at bandgap edges [168]. A major challenge that is yet to be tackled is to achieve rectification in a perfectly asymmetric fashion, without alteration of frequencies and along a straight line and a broad frequency band.

The coupling of acoustic or elastic waves with other types of waves, such as electromagnetic signals provides yet another avenue for promising phononic crystal applications. For example, a crystal may be designed to exhibit simultaneous phononic and photonic bandgaps coinciding at the same range of wavevectors [169–174]. This can be used to induce an interaction between the two types of waves in a defect-type cavity [175]. This coupling could in turn be exploited for signal manipulation or as a mechanism to slow down light. The simultaneous phononic/photonic bandgap concept is now often linked with the term “phoXonic crystal,”² and one of the first investigations of the concept can be found in Ref. [169]. The concept of wave coupling has been utilized in the context of acousto-optics as well [176]. A recent review article discusses developments in this area as well as other applications of phononic crystals and related concepts [177].

All the work reviewed thus far on phononic crystals (and other types of periodic media) has been based on the concept of passive control of wave propagation. While very practical in its own right in many applications, a natural extension is the development of techniques for active control or “tunability.” Such a prospect certainly holds promise as it implies an “on the fly” capability to change the location and size of bandgaps, as well as the shapes of the dispersion curves themselves. A simple approach for realizing active control of phononic crystals involves anchoring the unit-cell inclusions on axes and thus allowing them to be rotated on command [178]. As such, the orientation of the inclusions with respect to the incoming waves may be tuned thus altering the band structure during operation. A variety of innovative approaches have been proposed to achieve the same objective but by other means. Piezoelectricity, for example, has been utilized with noticeable effectiveness [179–182]. Other approaches involve the use of light-induced substrate potential in colloidal phononic crystals [183] or the use of electrorheological materials in which tuning is done via an external electric field [184]. Magnetism, in its part, has also been considered by several researchers as a platform for realizing phononic crystal tunability through alteration of material elastic properties. For example, Wang et al. [185,186] employed piezomagnetism and examined magnetoelectroelastic phononic crystals, and Robillard et al. [187] considered magnetoelasticity. Temperature control has provided yet another approach for tunability. This was explored as early as 2000 by Ruzzene and Baz [188] in the context of shape memory alloys whose properties are temperature dependent. Another phononic crystal study involving temperature investigated its effect on the bandgaps of surface and bulk acoustic waves [189]. Similarly, Jim et al. [190] considered the possibility of thermal tuning of a phononic crystal composed of a ferroelectric ceramic and epoxy. Finally, we also report on efforts that utilize finite elastic deformation [191] and lattice-level mechanical instabilities [192], both being approaches that take advantage of large elastic deformation.

While the above applications appear at a variety of scales, they all may be practically viewed as macroscopic in nature (e.g., may be modeled using continuum theory). The concept of a phononic crystal, however, can in principle be also realized at the atomic scale whenever the conditions for coherent phonon transport exist. Under such conditions, the dispersion characteristics of both the constitutive material lattice (defined by a primitive cell) and the phononic crystal lattice (now defined by a nanoscale supercell) contribute to the physics of the problem, which mostly pertain to transport properties, such as the thermal conductivity. This viewpoint in itself has opened up an entire track of research that overlaps with the emerging discipline of nanoscale heat transfer [193,194]. The reader is referred to Refs. [195–198] and [199–201] for, respectively, theoretical and experimental investigations of “nanophononic crystals.” More recently, a new analogy between the acoustic/elastic regime and the thermal transport regime has been introduced through the proposal of a “nanophononic metamaterial” whereby thermal conductivity reduction is achieved by local resonance [202]. These two concepts of nanoscale phononic materials provide promising opportunities for advancing a broad range of applications in heat management, such as those discussed in Ref. [203].

As mentioned earlier, the scope of interest in periodic materials received yet another surge in the early 2000s, and this has not been any less significant than the one brought about by phononic crystals in the early 1990s. The new phononic material that emerged has been denoted acoustic, or elastic, metamaterial, in analogy to an electromagnetic metamaterial. The prefix “meta” in this context is obtained from the Greek language to indicate the notion of “after” or “beyond.” This is linked to the understanding that metamaterials exhibit properties that are generally viewed to go beyond what we expect to find in naturally occurring, or conventional, materials. Examples of such unusual properties include the generation of subwavelength bandgaps and negative refraction. A key feature of acoustic/elastic metamaterials is the presence of local resonance, or possibly other traits that may bring about similar effects. While periodicity is not necessary for an acoustic/elastic metamaterial, its presence enables intrinsic, unit-cell-based, description of dynamical properties and wave propagation characteristics as well as order in the introduction of the locally resonant elements.

A seminal paper that has introduced a conceptual realization of an acoustic metamaterial was published in 2000 by Liu et al. [22]. In this paper, a 3D array of lead spheres coated with a 2.5-mm layer of silicone rubber was stacked in a simple cubic arrangement within an epoxy matrix. The lattice constant of the periodic material system was 1.55 cm. Upon excitation with acoustic waves, it was revealed that the medium exhibited a bandgap at a wavelength well below the regime corresponding to bandgap generation based on spatial periodicity of the impedance mismatch areas, which leads to “Bragg scattering.” In addition, the localized resonant structure was shown to cause the material to behave as if its effective elastic constants were negative at certain frequency ranges (precisely at frequencies around the bandgaps associated with the local resonances). Interestingly, the form of a locally resonant bandgap has appeared earlier in the Sigalas and Economou [16] paper although not for the purpose of exploring a new bandgap mechanism and the results were based on fictitious materials.

The use of local resonance as a mechanism for the generation of bandgaps and the achievement of other intriguing subwavelength properties has effectively opened up a new major thrust in phononic materials research. After the initial proposal of a three-phase composite material [22,204], numerous other configurations for locally resonant periodic materials have been proposed. These include the use of binary materials with the “inclusion” material being very soft [205,206], hollow cylinders or spheres [207], split rings or spheres [208,209], beams or plates with suspended masses [210–212], plates or surfaces with pillars [213–216], and structures with inertial amplifiers [217].

The notion of a locally resonant bandgap has been a source of studies concerned with analysis of the nature of the bandgap and how it differs from the more familiar Bragg scattering bandgap. Goffaux et al. [218], for example, have shown that the transmission spectrum corresponding to a locally resonant bandgap displays fano-like features. Xiao and coworkers [219,220] provided closed-form analytical expressions for the prediction of bandgap edges, and Liu and Hussein [221] provided a criterion for the transition between Bragg scattering and local resonance as the resonator properties are varied for beam structures with a periodic array of suspended masses.

Like their electromagnetic counterpart, acoustic/elastic metamaterials enjoy the allure of exhibiting unique dynamical behavior at long wavelengths, which is attributed to the presence of subwavelength bandgaps. This trait has made this type of materials an excellent candidate for the realization of desired homogenized properties at the long-wave end of the spectrum. Consequently, there has been a growing interest in deriving effective properties of locally resonant phononic materials.

Using the coherent potential approximation method [222], Li and Chan [223] derived the effective elastic modulus and effective density of a periodic composite material consisting of soft rubber spheres suspended in water (a choice motivated by the rubber's much slower speed of sound compared with that of water). They showed that at certain frequency ranges near the local resonance frequency, the effective bulk modulus and density become negative. Furthermore, their results have shown that there is a narrow frequency window at which both properties are simultaneously negative. The authors discussed how double negativity relates to rather exotic properties, such as negative refraction. They also correlated the observed phenomenon with the type of resonance modes in the local unit-cell structure. Specifically, it was explained that a negative elastic modulus is associated with a monopolar resonance, and a negative density is associated with a dipolar resonance. Upon reflection on the intriguing phenomenon of negative effective properties, one may interpret these as dynamical conditions whereby the restoring and inertial forces act in the direction of motion (rather than opposite to it) for the cases of effective elastic modulus and effective density, respectively.

Other studies followed with the aim of calculation, analysis, and/or utilization of negative effective properties, demonstrating negative effective elastic modulus [224–227], negative density (or mass) [228–232], or both [233–236]. A natural extension to these is the possibility of anisotropy in the effective mass or density [155,237–241], which is yet another intriguing concept in its own right. While most effective property calculation techniques are based on unit-cell analysis, it is also possible to compute these properties from the reflection and transmission coefficients that may be obtained by simulations or experiments [242]. Furthermore, while most studies have focused on actual, realistic models of locally resonant acoustic or elastic metamaterials, much can be learned by exploring the problem in the context of a simple “mass-in-mass” lumped parameter model as done by Huang et al. [243].

The topic of effective properties of acoustic/elastic metamaterials cannot be viewed in isolation from the classical theory of homogenization (e.g., see Ref. [244] and subsequent high-order homogenization formulations for solids [42–47], and homogenization of rigid bodies in a fluid [245]). However, in order to adequately obtain frequency-dependent homogenized properties, it is necessary to fully generalize the theory of homogenization to the dynamic regime and to do so in a manner that allows for an accurate regeneration of the dispersion relation at any wavenumber and frequency. This problem, commonly referred to as “dynamic homogenization,” has attracted significant attention in recent years [246–249].

As we discuss above, acoustic and elastic metamaterials draw their power from their low frequency, subwavelength properties. A nonexhaustive list of applications includes vibration and acoustic shielding, acoustic absorption, subwavelength focusing for fine-scale imaging applications, and acoustic/elastic cloaking. Below we briefly review some of the applications that have been studied. Reference [250] also provides an overview of state-of-the-art acoustic metamaterial applications.

An obvious application of a periodic material with a subwavelength bandgap is acoustic and vibration isolation, especially if this can be done over a relatively broad band of frequencies [251–253]. One may anticipate this concept finding application to sound proofing in buildings, vibration minimization in vehicle structures and home appliances, shock mitigation in armored gear, etc. In addition to isolation, where reflection is desired, researchers also considered sound absorption as an objective that also requires the minimization of reflection, e.g., see Refs. [212,254,255]. Such functionalities, would be even more effective if they can be realized in a tunable manner, which is in line with the goal of active control of phononic crystals. This problem has been pursued recently by means of fluid cavities separated by piezoelectric boundaries [256,257] and plates with a periodic array of piezoelectric transducers [258].

Another major category of applications of acoustic/elastic metamaterials is subwavelength focusing. The key principle here is utilization of operating points in the subwavelength portion of the medium's frequency spectrum whereby the isofrequency contour lines have opposing convexities across a straight path in the Brillouin zone (BZ) [208,259–265]. This concept provides a powerful alternative to conventional imaging and nondestructive inspection technologies because it essentially breaks the conventional diffraction limit.

Perhaps, the most intriguing of applications of acoustic/elastic metamaterials is cloaking, whereby an acoustic (or an elastic) wave is steered around an object rather than be scattered by it. A common approach for this problem is based on transformation acoustics (or elasticity) theory, whereby a target material distribution around the object of interest is mathematically synthesized from the point of view of heterogeneity and anisotropy [238,266–274]. Phononic crystals and acoustic/elastic metamaterials are viewed as promising candidates for practically realizing the unusually extreme material properties that emerge from this transformation process [239,271,275–284]. This is an open area of research that promises to continue to attract the attention and imagination not only of scientists, but also of the general public as a whole. For further reference, the reader may consult a review article [285] and a journal special issue [286] dedicated to the subject.

The resulting relation Ω = Ω(μ), obtained by imposing a range of wavenumbers and solving for frequency, is plotted in Fig. 2, where some of its important properties can be highlighted as discussed in the following.

Periodicity. The dispersion relation is symmetric about μ = 0 and periodic with period 2π. The behavior of the relation in the fundamental period $μ∈[-π,+π]$ therefore fully describes the characteristics of harmonic plane waves in terms of the relation between their wavelength and the frequency of wave motion. This space/wavenumber duality is a feature common to all periodic media. The range of wavenumbers defining the fundamental period of the dispersion relation, i.e., $μ∈[-π,+π]$⁠, is known as “first Brillouin zone” (FBZ) [24], which is generally defined by the reciprocal lattice, the dual of the direct lattice describing the geometry of the periodic domain. The symmetry of the dispersion relation provides the opportunity to fully characterize the dispersion properties of the lattice by only considering half of the FBZ, i.e., $μ∈[0,+π]$⁠, which defines the “irreducible Brillouin zone” (IBZ). The periodicity also indicates that for a given value of Ω, the dispersion relation is verified for μ + 2mπ, where m is any integer, i.e., for wavenumbers $k+2mπ/d$⁠. The FBZ, corresponding to m = 0, is bounded by a maximum wavenumber k = π/d. This corresponds to a minimum wavelength equal to λ = 2d, which defines a wave where two neighboring masses move out-of-phase. Conceptually, all waves of wavelength λ = 2d/(2m + 1) satisfy this condition, which supports the wavelength indetermination discussed above. However, for the particular case of the discrete 1D lattice, the cases for $|m|>0$ describe wave motion with several intermass oscillations, which are not justified by the lumped parameter model and the associated degrees of freedom (DOF). This issue is, however, of great importance and introduces complexity in the analysis of periodic structures defined by unit cells with multiple degrees of freedom or of continuous periodic structures.

which leads to the following relation between wave amplitudes $u˜(μ)=u˜(μ')$⁠, i.e., $u˜(μ)=u˜(μ+2mπ)$⁠. Equations (2) and (3) express the plane-wave solution for harmonic propagation. In this regard, the wave amplitude is the amplitude of motion for frequency ω and wavenumber κ = μ/d, which defines the spatial spectral harmonic content of the wave at the specified frequency and wavenumber. Therefore, Eq. (8) expresses the periodicity of the spatial spectrum of the propagating wave, which is fully consistent with the spatial sampling of the continuous rod as postulated by Shannon–Nyquist theorem.

Dispersion. The nonlinearity of the dispersion relation indicates that the periodic lattice is dispersive in contrast to the continuous rod, which is nondispersive. This is another effect of the discretization process. In the long wavelength limit, the dispersion relation behaves linearly and matches that of the continuous rod, where Ω ≈ μ (Fig. 3). The phase velocity in the lattice is far from being constant and is instead given by

Comparison of the dispersion relations for the continuous rod and the lattice, showing the matching long wavelength behavior and of the variation of phase velocity in terms of frequency, is shown in Fig. 3. The obtained expression for the phase velocity demonstrates how the attempt to use a system of the kind of Fig. 1 to estimate the velocity of sound, leads to a complex result, while the application of partial differential calculus to the continuous system provides a much simpler value for the wave speed.

Propagation and attenuation bands. The dispersion relation solved in terms of frequency for assigned wavenumber provides solutions in a limited range, i.e., for Ω ∈ [0, 2]. For Ω > 2, there are no corresponding real-valued wavenumbers, and harmonic plane waves do not propagate. The range of Ω > 2 therefore corresponds to an attenuation range, also denoted in the literature as a “stop band” or, when this phenomenon occurs over a limited frequency range, a “bandgap.” In contrast, for Ω ∈ [0, 2], plane-wave propagation occurs without attenuation. The amount of attenuation encountered within the stop band can be quantified by solving the dispersion relation at a given frequency in terms of the wavenumber. This allows the estimation of both the real and imaginary parts of the wavenumber as functions of frequency. Figure 4 shows the imaginary part of the wavenumber, which is equal to zero in the propagation range and different than zero in the stop band. In this case, the imaginary part quantifying wave attenuation increases with frequency, which indicates that the rate of spatial decay for a wave propagating through the lattice also increases with frequency.

The plots of Fig. 5 show the frequency-domain variation of the steady-state response of the last of N = 10 masses of a finite lattice excited by a harmonic force, f₀ applied to the first mass. A series of resonant peaks (10) within the pass band, and strong response attenuation occurring within the stop band, are noteworthy features of the frequency response plot.

The spatial attenuation is further demonstrated by plotting the dynamic deformed shapes at frequencies just below and right above the stop band boundary at Ω = 2, for a lattice of N = 100 masses. Within the pass band all masses undergo motion, which is an indication that waves propagate through the chain, while within the stop band, the amplitude of motion decays through the lattice. The magnitude of such decay, as discussed above, increases with frequency as clearly demonstrated by the deformed dynamic shapes of Fig. 6.

A simple extension to the spring-mass lattice problem is the diatomic lattice of Fig. 8, which features two different masses m₁ and m₂. In addition to being an interesting academic exercise, this case has physical relevance as it describes wave propagation in the [111] direction for NaCl and other crystal lattices [23].

where $un=[u2n,u2n+1]T$ contains the degrees of freedom of the nth unit cell of the lattice. Under the assumption of lattice periodicity, the notation can be simplified by letting K_n = K₀, where K_m, with m = −1, 0, +1, are matrices defining the spring interactions within a representative unit cell and with its neighbors. Because of periodicity, these matrices are identical for all cells and are thus independent of the cell location within the lattice. Of note is the fact that the description of the system in terms of its mass and stiffness matrices suggests the extension of the techniques to more complex systems, with multiple degrees of freedom resulting, for example, from a FE discretization.

Two different approaches can be followed in order to estimate the dispersion relation of the media. The first technique generalizes the approach presented for the spring-mass system presented above, whereby the propagation constant μ is imposed to solve the resulting expression in terms of frequency. This, under the fundamental assumption of absence of dissipation, in general, leads to a simple, linear eigenvalue problem which has a convenient solution. The result of this procedure, however, does not provide any information regarding attenuation, and only predicts the real dispersion branches. This technique is often denoted as the “inverse method.” Extensions of the techniques to the analysis of systems with dissipation are an active area of research, as described in upcoming sections (see Sec. 4.1). A second technique, known as the “direct method,” leads to the dispersion relations by imposing frequency and solving for the wavenumber. The most popular implementation of the direct method formulates the “transfer matrix” of a unit cell, and evaluates its eigenvalues. The direct method is easy to implement for general 1D periodic system and can also be applied to FE discretizations of complex systems. Its main drawback consists of its difficult extension to the analysis of 2D systems, which is still an open research problem. The two techniques will be briefly described in the following.

The direct method and transfer matrix. The implementation of the direct method generally relies on the formulation of the transfer matrix of the lattice, which relates the generalized displacements and forces at the left and right cell boundaries (see Fig. 10). The description of the behavior of the unit cell in Fig. 10 in terms of the dynamic stiffness matrix $D(ω)=K-ω2M$ is a convenient starting point to relate boundary forces to boundary displacements. In expanded form, the dynamic stiffness relation can be written as

is the transfer matrix.

A variant to the case of the 1D monatomic lattice of Fig. 1 consists of an addition of a second-order oscillator of spring constant k_R and mass m_R, connected to each base mass (Fig. 12). This leads to a configuration characterized by internal resonances centered at the resonant frequency of the oscillators, i.e., $ωR=kR/mR$⁠. This produces a bandgap in the vicinity of ω_R. In contrast to the bandgaps discussed for monatomic and standard diatomic lattices, the internal resonance gap is not related to wavelength, and therefore in principle can be produced at arbitrarily low-frequencies/long wavelengths. The bandgaps observed in monatomic and standard diatomic lattices instead only occur at wavelengths which are of the order of the unit-cell size, and therefore are strictly constrained by the unit-cell dimensions. These are the result of interference generated by the interaction of incident and scattered waves at the unit-cell boundaries, which is the phenomenon often referred to as Bragg scattering. The wavelength-independent nature of internal resonance bandgaps has important implications, which have motivated several studies for practical implementation of the concept [22], along with numerous studies devoted to the interpretation of the physical meaning of the locally negative slope of the dispersion relations caused by the internal resonance, which is an indication of negative group velocity [224]. Interpretations in terms of negative mass can be found, for example, in Refs. [240,243,296]. The reader is referred to Sec. 1.4.3 for a detailed overview of the literature on locally resonant metamaterials, and to other parts of this review where the notion of local resonance has been considered.

The frequency of resonance Ω_R can be selected, or “tuned,” to obtain attenuation over a desired band, without any constraint imposed by wavelength. Of course, in practice, low-frequency attenuation requires large mass values and/or low stiffness, which poses a challenge for the implementation of this concept. However, the application of the internal resonator concept for systems more complex than the one illustrated in this section has been broadly pursued by the research community, and several experimental demonstrations are also available. Examples of configurations and experimental results are presented in Sec. 5 of this paper.

The spring-mass lattice can be extended to obtain a 2D periodic structure that develops over a plane (Fig. 14). The unit cell consists of a mass undergoing motion perpendicular to the plane of the structure and connected to its neighbors through two springs that apply a force proportional to the relative out-of-plane motion of two neighboring masses. This example, which was originally introduced in Ref. [76], can be described through a single degree-of-freedom unit cell, which is sufficient to uncover additional interesting properties of periodic systems. According to the configuration of Fig. 14, wave motion is investigated by expressing the equation of a generic unit cell at location m, n

The strong directivity is further illustrated by computing the harmonic response of a finite lattice at the chosen frequencies. A finite 81 × 81 lattice in free–free configuration is excited by a concentrated load applied at the center of the lattice. The steady-state harmonic response is computed by inverting the dynamic stiffness matrix of the assembly, and is plotted in the form of the color map of the response magnitude. Proportional damping is added in order to reduce the effects of boundary reflections and better visualize the response directivity. The response plots for the two considered frequency values are shown in Fig. 18, which demonstrate the difference in dynamic behavior, and, in particular, the strongly directional response of the lattice at $Ω=2$⁠.

The occurrence of wave “beaming” due to limited directionality associated with specific frequencies of wave propagation is illustrated in this simple case, and has been broadly investigated in the literature. Directional wave propagation in crystal structures of various kinds of anisotropy has been studied extensively and documented both numerically and experimentally [76]. Strong focusing of waves in specific directions is associated with the occurrence of “caustics” [84,297], which correspond to closed circles in the group velocity plots as shown in Refs. [80,298,299]. Directional wave propagation has been investigated in periodic structural components, such as beam grillages, hexagonal, re-entrant and chiral lattices, and periodic plates among others. In addition, the analysis of dispersion surfaces and group velocity directions has provided the foundation for the design and implementation of phononic crystals with negative acoustic refraction and focusing capabilities [156,160,163,260,261,300, 301,302].

allows us to obtain the frequency band structure along any desired path in, or at the boundaries of, the BZ. In Eq. (36), $u˜$ is the unit-cell Bloch displacement function that is periodic over the unit cell, x = {x, y, z} is the position vector, and $κ={κx,κy,κz}$ is the wavevector.

This section reviews a variety of computational techniques for analyzing response and designing continuous phononic systems. Included are descriptions of the plane-wave expansion methods (PWE), finite difference (FD) methods, the FE methods, and the multiple scattering (MS) method. The literature reviewed is not intended to list exhaustively those papers using each method, but instead to document early literature employing the method and extending it, and to provide the proper technical review necessary for a reader to begin using each method.

The general approach used in the plane-wave expansion method for continuous systems is to expand the sought-after solution field and the material properties in a Fourier series, and to then invoke orthogonality of the basis functions to solve individually each of the introduced solution coefficients. A very early description of this method can be found in the beginning of Sec. 2 of Elachi's review article [85], although the method is not given a name therein. Phononic crystal research in which this approach has been taken date back to the early 1990s [17,18] with work on acoustic and elastic composites.

where $r=xx∧$ denotes position, $κ=kx∧$ denotes the wavevector, and $∑Gl,∑Gm$ denote sums over all reciprocal lattice vectors (e.g., over all l, m). The transformed density and stiffness are denoted by $ρ(Gm)$ and $τ(Gm)$⁠, respectively. The transformed displacement is given by $uκ(Gl)$ for any wavevector $κ$⁠. Note that Eq. (38) embodies Bloch's theorem.

Used in arriving at Eq. (41) is the fact that for each $Gn$⁠, the only nonzero terms in the inner product are those satisfying $Gl+Gm-Gn=0$⁠, leading to the condensing of the summation over $Gm$⁠. By truncating the expansions for $Gl$⁠, and evaluating Eq. (41) for a similarly truncated set of $Gn$⁠, a well-posed eigenvalue problem results for eigenfrequencies $ω(κ)$ and eigenvectors $[uκ(G-N)uκ(G-N+1)…uκ(GN)]T$⁠, where 2N + 1 denotes the number of retained terms. The band structure is then obtained by sampling $κ$ in the first Brillouin zone. As in any approximation technique, a convergence criterion should be used to determine an appropriate value for N. Convergence is discussed in Ref. [303], where it is pointed out that it may be appropriate to expand the inverse of $ρ(x)c2(x)$ to achieve faster convergence.

The finite difference method is a classical method that has seen wide application in almost every discipline where partial differential equations (PDEs) govern behavior. It has one critical advantage over PWE methods in that it can easily analyze periodic media composed of multiple states (e.g., solid and fluid) [304]. This includes a large class of materials with voids. The earliest use in a phononic system known to the authors is that by Yang and Lee [305]. They found the band structure of a two-layered one-dimensional composite (what is now referred to as a bimaterial rod) assuming harmonic time dependence and using a spatial finite difference approximation of the resulting governing equation. The authors also enforce a Floquet (or Bloch) condition on the meshed domain, effectively limiting the computation to a single unit cell. Solution of a discrete eigenvalue problem at each frequency yields the dispersion curves.

A more frequently encountered finite difference approach is one in which a PDE is discretized in both space and time to compute the time-domain response. Frequency content is then acquired from Fourier transforms (FTs) of the time-history data. This subapproach is appropriately termed the finite difference time-domain (FDTD) approach. Note that if the Floquet (or Bloch) condition is not imposed in arriving at the PDE, the computational domain must include many unit cells, an excitation source, and a choice of boundary conditions. For two-dimensional phononic crystals, such an approach was first presented by Garcia-Pablos et al. [304]. By varying the source's excitation frequency, attenuation versus frequency can be characterized and used to identify bandgaps. In contrast, Tanaka et al. [306] presented an FDTD approach in which the PDE is formulated with the Floquet (or Bloch) condition, thus requiring simulation only over a unit cell. By varying the wavevectors through the Brillouin zone, Fourier transforms of the time-history data can be used to generate band structures.

Substitution of Eq. (43) into Eq. (42) yields the desired n by n complex eigenvalue problem. Sweeping through μ = ka and solving the eigenvalue problem at each wavenumber generates the band structure.

where M denotes a mass matrix, K denotes the stiffness matrix, and q stores the nodal degrees of freedom. Note that these equations are arrived at by considering a unit cell in isolation from other unit cells, and thus neighbor forces do not appear on the right-hand side; other treatments can be found in Refs. [82,307,308] where neighbor forces are included, but lead to the same final form. As can be expected, the finite difference equations Eq. (42) take this same form if rearranged with the discretization length h incorporated into the stiffness matrix; one notable difference is that in the FE treatment, the mass matrix is not necessarily diagonal. Since the one-dimensional case was discussed for finite difference modeling, the two-dimensional case will be discussed here for the finite element case.

In arriving at Eq. (44), all elements lie wholly within the unit cell and all boundary nodes lie on the perimeter of the unit cell (i.e., not inside of the unit cell). Thus, all unit-cell energy is accounted for, and no energy associated with M or K belongs to a neighboring unit cell. However, an excess number of nodes exist in q. Specifically, nodal degrees of freedom at the top and left sides of the unit cell (see Fig. 19) need to be condensed out since they belong to neighboring unit cells.

where $M∧(μx,μy)=T¯TMT$ and $K∧(μx,μy)=T¯TKT$⁠. Sweeping through all values of μ_x and μ_y along the irreducible Brillouin zone and solving the eigenvalue problem generate the band structure.

where $U˜$ is the discrete Bloch vector. This equation may be solved for different values of $κ$ to generate the band structure. More details on the FE discretization of Eq. (49) are available in Ref. [309].

The MS method appears to be first proposed for evaluating band structure in periodic elastic composites by Economou and coworkers starting in 1999 [310,311] and by Liu and coworkers in 2000 [312]. As the name of the method implies, it is particularly useful in phononic systems composed of scatterers in a host medium (in contrast to layered phononic systems). The method can handle multiple media states (e.g., elastic scatters in a fluid, or air spaces in an elastic medium) and can accurately analyze high-contrast problems where the other methods (PWE, FD, FE) may exhibit slow convergence.

An excellent review of the multiple scattering method and its use in determining both band structure and frequency-dependent transmission coefficients can be found in Ref. [313]. Reviewing the method, even in simple cases, requires a lengthy development beyond the scope of the present review. For the purposes herein, we only present the general ideas that guide the method.

Note that a scatterer located at R_n generates a pressure field $pnsc(r)$ that in itself arises from incident pressure waves from neighboring scatterers. If the single scatterer problem has been solved (i.e., the transmission coefficient is known for an incident spherical wave on the scatterer), a linear algebraic system can be established that determines amplitude coefficients for each scattered wave leaving a scatterer at R_n. For an acoustic medium, the propagation away from the scatterer is then the product of the amplitude coefficients with spherical Hankel and spherical harmonic functions whose arguments contain r – R_n, as per a standard Green's function approach. Note that the problem as described is infinite in size, so computed solutions require truncating the number of neighbors used to describe the waves scattered at R_n. Finally, if band structure is desired, Bloch's theorem can be invoked to relate amplitude coefficients for a scatterer in one unit cell to those of an analogous scatterer in a reference unit cell.

As thoroughly discussed in Secs. 1 and 2, the study of wave propagation in periodic media in general utilizes Bloch's theorem, which allows for the calculation of dispersion curves (frequency band structure) and density of states. Due to crystallographic symmetry, the Bloch wave solution needs to be applied only to a single unit cell in the reciprocal lattice space covering the first BZ [24]. Further utilization of symmetry reduces the solution domain, even more, to the IBZ. As mentioned in Sec. 3.1, there are several techniques for band structure calculations for phononic crystals and acoustic/elastic metamaterials (which are also applicable to photonic crystals and electromagnetic metamaterials). Some of the methods involve expanding the periodic domain and the wave field using a truncated basis. This provides a means of classification in terms of the type of basis, e.g., the PWE method (Sec. 3.1.1) involves a Fourier basis expansion and the FE method (Sec. 3.1.3) involves a real space basis expansion. The pros and cons of the various methods are discussed in depth in the literature [314].

Regardless of the type of system and type of method used for band structure calculations, the computational effort is usually high because it involves solving a complex eigenvalue problem and doing so numerous times as the value of the wavevector $κ$ is varied. Therefore, means for speeding up band structure calculations are very beneficial. Examples of techniques that achieve this speed-up include utilization of the multigrid concept [315], development of fast iterative solvers for the Bloch eigenvalue problem [316,317], extension of homogenization methods to capture dispersion across the spectrum [43,318,44], and formulation of multiscale finite elements as illustrated recently in Ref. [319].

In this section, we describe a model reduction method that is based on modal transformation [309]. This method, which is referred to as the reduced Bloch mode expansion (RBME) method, involves carrying out an expansion employing a natural basis composed of a selected reduced set of Bloch eigenfunctions. This reduced basis is selected within the IBZ at high-symmetry points determined by the crystal structure and group theory (and possibly at additional related points). At each of these high-symmetry points, a number of Bloch eigenfunctions are selected up to the frequency range of interest. As mentioned above, it is common to initially discretize the problem at hand using some basis of choice. In this manner, RBME constitutes a secondary expansion using a set of Bloch eigenvectors and hence keeps and builds on any favorable attributes the primary expansion approach might exhibit. The proposed method is in line with the well-known concept of modal analysis, which is widely used in various fields in the physical sciences and engineering. To follow, a brief description of the RBME process is given and its application in a discrete setting (e.g., using FEs) is demonstrated for a phononic crystal problem.

The starting point for the RBME method is a discrete generalized eigenvalue problem emerging from the application of Bloch's theorem. For demonstration, Eq. (50) will be used for this purpose. The Bloch eigenvalue problem is first solved at a reduced set of selected wavevector points (i.e., reduced set of $κ$-points). This provides the eigenvectors from which a reduced Bloch modal matrix, denoted Ψ, is formed. Several schemes are available for the selection, the simplest of which is the set of eigenvectors corresponding to the first few branches at the high-symmetry points, e.g., the Γ, X, and M points for a 2D square lattice as illustrated in Fig. 20 (more details on selection schemes are given in Ref. [309]). The matrix Ψ is then used to expand the eigenvectors $U˜$⁠, i.e.,

where M_R and $KR(κ)$ are reduced generalized mass and stiffness matrices. The eigenvalue problem given in Eq. (53) can then be solved for the entire region of interest within the IBZ at a significantly lower cost compared with using the full model given in Eq. (50).

To demonstrate the RBME approach, we analyze a linear elastic, isotropic, continuum model of a 2D phononic crystal under plain-strain conditions (the method is also applicable to periodic metamaterials). As an example, a square lattice is considered with a bimaterial unit cell (with ratios of Young's moduli and densities of E₂/E₁ = 16 and ρ₂/ρ₁ = 8, respectively). The topology of the material phase distribution in the unit cell is shown in the inset of Fig. 20 (the filling ratio of the inclusion is 0.1111). The unit cell is discretized into 45 × 45 uniformly sized four-node bilinear quadrilateral FEs, i.e., 2025 elements. With the application of periodic boundary conditions, the number of degrees of freedom is n = 4050. Figure 20 shows the calculated band structure (a) and density of states (b) using two-point expansion, that is, the selection is carried out at the Γ, X, M points in $κ$-space. In the calculations, eight modes were utilized at each of these selection points. As such, a total of 24 eigenvectors (m = 24) were used to form the Bloch modal matrix. The results for the full model are overlaid for comparison indicating excellent agreement, despite a reduction of model size from 4050 to 24 degrees of freedom. For models with a larger number of degrees of freedom, and a calculation with high $κ$-point sampling, two orders of magnitude or greater reduction in computational expense will be achieved (as shown in Fig. 20(c)).

The presence of damping is unavoidable in materials and structures. Its consideration therefore improves the predictive ability of any model of a phononic system [320]. From a characterization perspective, damping can be identified from experiments [321]. From a design perspective, damping may be introduced intentionally to provide stronger attenuation and/or more control of the wave propagation characteristics, as will be demonstrated in this section. In both the prediction and design scenarios, the treatment of damping within the framework of dispersion analysis provides a fundamental description of the effects of dissipation on the intrinsic properties of a phononic material or structure. This is more powerful than consideration of the effects of damping on finite systems where only extrinsic (i.e., structural level) response behavior may be observed.

In this section, we consider the problem of damped phononic materials and structures. We first review relevant studies from the literature in Sec. 4.1.1, and follow in Secs. 4.1.2 through 4.1.4 with a detailed technical review on the treatment of viscous damping in the calculation of the dispersion curves. The technical review will consider (1) a 1D diatomic chain like the one depicted in Fig. 8 except with viscous damping elements included between the masses; (2) a 1D diatomic chain incorporating local resonators similar to the one depicted in Fig. 12 but with damping included; and (3) a 2D continuum model of a damped phononic crystal discretized by the finite element method.

There are numerous studies in the literature that consider the treatment of damping in the context of periodic materials or structures. Some of these studies are concerned with finite structures [322–324] or the transmission [325] or Green function [326] response in a portion of a periodic material or structure. In this review, the interest is in the examination of materials and structures mathematically modeled as infinite systems since this provides intrinsic wave propagation properties as discussed above. Research in this area may be divided into two categories. The first seeks a wave dispersion solution for prescribed harmonic motion and the second is concerned with obtaining a dispersion relation for free wave propagation. In practice, the first problem is useful for comparison with experiments on long phononic structures, or portions of a phononic material, subjected to forced harmonic excitation. The second problem is suited for investigation of the natural dynamic characteristics of phononic materials and structures or when the loading conditions are of an impulsive nature.

A key distinction between harmonic and free wave propagation lies in whether the frequency is restricted to real values (harmonic) or allowed to be complex (free). The wavenumber (or wavevector) can be either, depending on whether the interest is in both propagating and evanescent modes (for which the wavenumber has to be complex) or in only propagating modes (for which it is sufficient to only consider a real wavenumber).

The study of harmonic wave propagation in damped periodic materials and structures in the context of a linear eigenvalue problem goes back to Mead [6] in which a monocoupled periodic chain with a complex stiffness connecting each unit was considered. Langley [327] investigated a similar model only with the damping introduced in the inertial term. Merheb et al. [328] and Zhao and Wei [329] also incorporated damping in a linear eigenvalue problem setting but considered more complex frequency-dependent damping models that carried information on the time history of the dissipation mechanism. The approach by Merheb et al. [328] was based on the FD time-domain method, whereas Zhao and Wei [329] used the PWE method. Recently, Manconi and Mace [330] used finite elements to treat frequency-dependent damping in the context of harmonic wave propagation in laminated panels. On a different track, treatment of harmonic wave propagation through the setting up of a quadratic eigenvalue problem has been considered for both 1D models [331–333] and multidimensional models [334–336]. The advantages of the quadratic eigenvalue problem route are twofold: (1) it easily enables the incorporation of frequency dependency in the material properties and (2) it readily provides the wave solution for both spatially propagating and spatially decaying modes.

Calculation of dispersion relations for free wave propagation in damped phononic materials and structures may be traced as far back as the 1975 paper by Mukharjee and Lee [337] that studied laminated media with structural damping and where the finite difference method was used to obtain dispersion relations for various levels of damping. Sprik and Wegdam [338] and Zhang et al. [339] later considered similar problems—in which damping was introduced by adding an imaginary component to the elastic modulus—using the plane-wave expansion method. Recently, Hussein [340] and Hussein and Frazier [341] derived dispersion relations for viscously damped phononic materials. The treatment of temporally nonlocal viscoelasticity was also considered [342]. In all studies concerned with free wave propagation, a complex frequency is assumed in the application of Bloch's theorem. This provides the usual frequency dispersion curves—now affected by the presence of damping—and also curves associated with the temporal attenuation for each Bloch mode. A wavenumber-dependent damping ratio diagram may be extracted from the latter to explicitly determine the loss factor for each Bloch mode [340–342].

To demonstrate the treatment of damping in free wave propagation analysis, we start with the simple 1D diatomic lattice model studied in Sec. 2.2, and add linear viscous damping elements between each of the masses (Fig. 21). We provide details of a procedure for obtaining the dispersion curves as described in Refs. [340–342]. Considering periodicity, the set of homogeneous equations describing the motion of each mass in the unit cell is obtained as follows:

where $uαj$ is the displacement of mass α in an arbitrary jth unit cell. In this notation, a unit cell and its neighbors may be identified by j + r, where r = 0, −1, +1 denote the present, previous, and subsequent unit cells, respectively.

Note that $ωs(κ)=Abs[γs(κ)]$⁠. For the special case of proportional (Rayleigh) viscous damping, the value of ω_s is equal to the undamped frequency [343].

To get some insight into the effects of damping on the wave propagation characteristics of phononic materials and structures, we generate dispersion curves for the 1D chain as shown in Fig. 21. In addition to customary frequency dispersion curves—which for our damped model we obtain from Eq. (63)—we may also plot curves for the Bloch mode damping ratios using Eq. (64). For demonstration, a specific set of material parameters are used: r_m = m₂/m₁ = 9, r_c = c₂/c₁ = 0.5, r_k = k₂/k₁ = 1, and $ω¯=k2/m2=149.07 rad/s$⁠. A damping parameter, β = c₂/m₂, is varied to show the dependence on the damping intensity. In Fig. 22, a frequency band diagram appears on the top and a damping ratio band diagram appears on the bottom. While these results are specific to the values chosen for the parameters r_m, r_c, and r_k, they provide a basic insight on the effects of damping on the elastodynamic behavior of a periodic material or structure.

We observe in Fig. 22 (top) that as the damping intensity $β/ω¯$ is increased, the optical branch drops while the acoustic branch experiences little change—this in turn leads to a reduction in the size of the bandgap. It also leads to a branch-overtaking phenomenon [340–342] whereby a higher branch overtakes a lower branch as damping is raised. Furthermore, we observe that the descent of the optical branch takes place at slightly faster rates at low wavenumbers compared with high wavenumbers. The damping ratio diagram in Fig. 22 (bottom) follows a consistent trend with an indication that the effect of damping is slightly more significant at low wavenumbers. In looking further at both diagrams in tandem, we observe that the optical branch cuts-off for long wavelengths that are overdamped, and the cut-off wavenumber is precisely at ζ₂ = 1 (which corresponds to the critical damping ratio) [341,342]. This cut-off phenomenon implies that under high damping conditions, there is a possibility for wavenumber bandgaps to arise as illustrated in Fig. 23, which displays a bandgap map for both frequencies and wavenumbers.

The most striking contrast is that in the latter case, the drop rate of the optical branch is significant at high wavenumbers compared with low wavenumbers. This is reflected in the rapid rise of the dissipation factors at high wavenumbers as shown at the bottom of Fig. 25. We also observe that at low wavenumbers, the damping ratio values corresponding to the optical modes for the conventional lattice exceed those of the locally resonant lattice for a given damping intensity $β/ω¯$ and vice versa at high wavenumbers. With regard to the damping ratios of the acoustic branch modes, these are similar at low wavenumbers but are higher for the regular lattice at high wavenumbers.

The qualitative and quantitative differences revealed in the dissipation diagrams as shown at the bottom of Figs. 24 and 25, respectively, have recently been examined in more detail [344]. In order to yield some insights into the difference in the overall damping capacity, a comparison was made between the same two models except with a slight alteration in the parameter values in order to allow both systems to exhibit an identical value of long-wave speed. It was shown that despite both systems having the same values of prescribed damping in the dashpots, and being statically equivalent, the locally resonant lattice (i.e., the metamaterial) exhibits more dissipation across the spectrum. This phenomenon, which has been termed metadamping, is indicative of an emergence of damping due to the presence of local resonance and may be utilized toward the design of materials that need to be simultaneously stiff and dissipative.

where C is the conventional elasticity tensor, $αd$ is the velocity-dependent dissipative component of the stress tensor.