The attenuation performance of locally resonant acoustic metamaterials based on generalised viscoelastic modelling

Acoustic metamaterials are known as a promising class of materials interacting with acoustic and/or elastic waves. Band gap formation is one of the most spectacular phenomena that they exhibit. Diﬀerent ways to broaden the attenuated frequency ranges are still being actively explored. It turns out that material damping through intrinsic viscoelastic material behaviour, if accurately tailored, may contribute to the enhancement of the performance of a properly designed acoustic metamaterial. In this study, a locally resonant acoustic metamaterial with periodic multicoated inclusions with viscoelastic layers is investigated. Multiple attenuation regimes obtained with the considered geometry are joined for a certain level of viscosity of the coating layer. The analysis is performed using a generalized Maxwell model, which allows for an accurate description of nonlinear frequency dependent elastic properties, and thus is widely used to model the behaviour of many polymeric materials in a realistic way. The study reveals that variation of the material parameters of the rubber coating with respect to frequency inﬂuences not only the position of the band gaps but also the eﬀectiveness of the wave attenuation in the frequency ranges around the band gaps.


Introduction
Acoustic metamaterials are a novel class of materials possessing some unusual properties, uncommon or non-existent in nature.The mechanism of low-frequency wave attenuation based on the local resonance [1] attracts scientific attention due to various potential applications, e.g., noise insulation [2] or energy harvesting [3].Unlike phononic crystals (PCs), which is another class of materials capable of forming band gaps (frequency ranges where no wave propagation occurs), the design of locally resonant acoustic metamaterials (LRAMs) does not require periodicity and allows generating subwavelength band gaps [4].The number of occurring band gaps can be controlled through the microstructural design, for instance, by exploiting multicoated inclusions proposed by Larabi et al. [5].
Such a microstructure has been studied based on 1D model of dual-resonators in [6], and has been parametrically optimised in [7] and [8], leading to a reduced distance between band gaps.
However, the main limitation in terms of application of LRAMs is still the fact that the attenuated frequency ranges (even if they are multiple) are rather narrow.In the literature, a few solutions dealing with this limitation have been proposed, among which optimizing the metamaterial topology [9], coupling of the effect of local resonance with Bragg scattering [10,11] and using resonators with distributed resonant frequencies [12,13].
More recently, also the potential influence of material losses on broadening attenuation regions, has started attracting researchers' attention [14,15].Such a solution seems to be particularly promising considering the fact that material damping is an intrinsic feature of polymeric materials typically used in LRAMs.
Studies on damped periodic structures started with the works of Mead [e.g. 16] in his analysis of a one dimensional periodic chain of masses with lossy springs, and by Mukherjee and Lee [17] who have investigated transient effects in damped laminates.Until recently, available studies mainly focused on PCs.In a number of works investigating one-and two-dimensional PCs, broadening of band gap regions due to the presence of viscoelastic components has been observed [18,19].It has been shown by Oh et al. [20] that lossy PCs are more effective than homogeneous viscoelastic media in terms of energy dissipation.
Moreover, in some cases where more advanced material models have been considered (e.g.generalised Maxwell), a shift of the attenuation regions, due to the frequency dependent storage modulus, has been reported [21].These predictions have been confirmed by the experimental analysis conducted by Merheb et al. [22].Currently, the studies of damped PCs are taken even further, for instance by including defect modes [23].
The influence of damping on the periodic material's performance is not restricted to the band gap width and position.The entire band structure changes significantly, as reported in various studies [e.g.24,25].Based on the complex band structure representation, Moiseyenko and Laude [24] have found that the losses in a phononic crystal have stronger impact on the real part than on the imaginary part of the wave number.However, as shown in [26], high damping ratios can actually support the wave propagation within attenuation regions by decreasing the imaginary part of the wave number.In the work of Hussein and Frazier [27], also the phenomena of branch overtaking and branch cut-off in the band structure have been observed and studied.
It turns out that the influence of material damping on the performance of LRAMs differs significantly from the case of PCs.First of all, not only the losses in the matrix material but also within the rubber-coated inclusions should be considered.Manimala and Sun [28] have shown, considering three types of viscoelastic models (Kelvin-Voigt, Maxwell and Zener), that tailoring the damping within the resonators instead of relying on the dissipation in the matrix material might be beneficial for broadening the attenuation spectrum.This has also been confirmed in the works of Wang et al. [14] and Krushynska et al. [15], using a locally resonant acoustic metamaterial with a single coated inclusion.
The viscoelastic behaviour of the rubber coating has been shown to have a critical impact on the material performance in comparison with the damped matrix.Moreover, in case of a LRAM, the imaginary parts of wave numbers are predominantly influenced by the material losses and, as a result, the attenuation peaks related to the local resonances are smoothed.The dissipative effect that leads to this response has also attracted attention of researchers.As a consequence, the notion of metadamping has been introduced by Hussein and Frazier [29] as an enhancement of material dissipation due to the presence of local resonance.Later studies [18] on both viscoelastic PCs and LRAMs have concluded that the effect of damping on the band gap size is actually more pronounced in case of PCs.However, the experimental and theoretical studies performed by Zhao et al. [30,31] have shown a wide absorption range at low frequencies for a composite polymer slab with embedded local resonators, due to the dissipative mechanisms in the coating material.
So far, the analyses of viscoelastic LRAMs typically assume simple linear viscoelastic models, like the Kelvin-Voigt model, in order to describe the viscoelastic behaviour of the constituents [30,14].With this model the material properties become complex, but the real part of the modulus (related to the elastic response) is still constant and only the imaginary part (associated with viscous behaviour) changes linearly with frequency.On the other hand, using the generalised Maxwell model, which has not been used extensively for locally resonant acoustic metamaterials, allows for a realistic variation of both terms with respect to frequency and as a consequence more realistically describes the material behaviour.This is important since most polymeric materials have properties that are known to be frequency dependent [32].
In this paper, the multiple band gaps obtained with a particular microstructure with coaxial multicoated inclusions have been joined using viscoelasticity of the coatings.Such a concept has recently been introduced in [8] and was studied based on a simple linear viscoelastic model in 1D (for longitudinal wave polarisation only).In the present study, an advanced viscoelastic model is used to describe the behaviour of the rubber coating, thus providing a more realistic insight into the influence of frequency-dependent material parameters on the locally resonant acoustic metamaterial performance.To this aim, a 2D analysis of a locally resonant acoustic metamaterial based on complex dispersion diagrams and power transmission spectra is conducted.The impact of viscoelastic material properties is studied in detail using the generalised Maxwell model for the coating layers in the inclusions.Furthermore, it turns out that a correct material model might be crucial if the focus is on exploiting viscoelasticity for joining band gaps.First, it is observed that due to the variation of elastic parameters with frequency, the band gap regions exhibit a shift.Depending on the elastic properties of the considered material, such a shift can be significant in some cases, which means that purely elastic predictions may not be sufficient to determine the band gap location.Therefore, while designing such a metamaterial, the dependence of the soft coating material behaviour on frequency should be verified.Secondly, if material damping is used for the purpose of joining band gaps, the target frequency range of wave attenuation (within the distance between the band gaps) needs to overlap with the region where the loss tangent level is sufficiently high.Otherwise the effect of bridging may not occur.
The paper is organised as follows, first, the modelling approach based on Bloch theory for obtaining complex dispersion diagrams and the finite element calculation of the power transmission spectra are described.Next, in section 2.2, details on the considered geometry and material properties are given.Finally, in section 3 the results of simulations are presented and discussed.Section 5 summarises the main conclusions of the paper.

Modeling approach
The dynamic characteristics of a material can be obtained through the study of harmonic wave propagation, which is typically analysed based on its dispersion relation: the relationship between frequency ω and wave number k.Some extraordinary properties of locally resonant acoustic metamaterials are depicted in such diagrams, already in the range of real wave numbers where the presence of band gaps may be captured [1,4].
However, due to the continuity principle for the dispersion curves, branches within the band gap ranges are present in the domain of complex wave numbers, wherein the imaginary part of the wave number is often used as a measure of wave attenuation.Therefore, by considering complex band diagrams (real frequencies and complex wave numbers) information on spatial wave propagation as well as spatial attenuation of elastic waves is obtained.The analysis of the complex bands of dispersion diagrams may also contribute to a better understanding of band gap formation since tracing the evanescent Bloch waves becomes possible [14].

Complex dispersion diagram: formulation based on Bloch theory
The classical way for obtaining a dispersion relation is based on Bloch theory.The fundamental theorem for wave propagation in a periodic, infinite material states that the wave field in such a medium is also periodic [33,34], and as a consequence, the analysis of such a structure can be restricted to a single unit cell.The application of Bloch theorem [35] allows to characterize the wave solution in an irreducible Brillouin zone in the following way [34, p. 204 where u denotes the displacement field, ũ the Bloch displacement function with spacial periodicity (the same as the periodicity of the structure), x the position vector, k the wavevector, ω frequency and t time.The dispersion relation can then be easily obtained [34].The Bloch solution (1) is substituted into the equation of motion with the additional assumption that the materials are isotropic and linear elastic.The strong form of the resulting problem is further converted into a weak form and discretised by means of the finite element method.The resulting eigenvalue problem is given by the equation: where K and M denote the stiffness and mass matrices, respectively.Since the eigenvalues λ = ω 2 are calculated with respect to the wave number k, this approach is typically called the ω(k) formulation.
In order to obtain the complex band structure with complex wave numbers and real frequencies, a reformulation of the eigenvalue problem (2) is needed, for which the formulation proposed in [36] is adopted.The analysis is restricted to in-plane wave polarizations, where the wave vector components k x and k y are linearly related as k y = αk x , with α being a real constant.The eigenvalue problem can then be written in a first-order form: where K, M denote the new stiffness and mass matrices, respectively, and the periodic displacement is given as u = [k x u, u] T .Equation ( 3) is solved with respect to real frequencies, for a single unit cell, assuming periodic boundary conditions.It is known as the k(ω) approach.
The k(ω) approach has been used here to incorporate viscoelasticity.Since the governing equations may be formulated with respect to frequency using the Laplace-Carson transform, Laplace transform [37] or other, accounting for frequency dependent mechanical properties becomes straightforward.Following such a transformation, the frequencydependent linear elastic constitutive relation reads: where κ and G denote the bulk and shear moduli, respectively, and ε d stands for the deviatoric part of the strain tensor.According to the elastic-viscoelastic correspondence principle [38], the dispersion relation in the case of a viscoelastic material can be obtained by replacing the elastic material parameters with their frequency-dependent counterparts in the constitutive equation ( 4), parameters κ and G have been replaced by complexvalued κ = κ(ω) and G = G(ω).A detailed description of the viscoelastic material model used in this work is given in section 2.2.

Model geometry and material parameters
A locally resonant acoustic metamaterial with coaxial multicoated inclusions is considered.The geometry of the periodic unit cell of the material is shown in Figure 1.
Infinite tungsten cylinders are coated with two layers of rubber alternating with one layer of tungsten.These multicoated cores are embedded in an epoxy matrix.The geometry and dimensions of the investigated unit cell (Figure 1) have been chosen in order to trigger multiple low frequency band gaps located close to each other, where the band gaps are formed due to the vibration of the heavy masses (in phase motion for the first band gap, out of phase for the second band gap).The filling fraction (the ratio between the area of the coated inclusions and the area of the entire cell) for the given geometry equals 60%, which is at a level where the interactions between adjacent inclusions are not yet affecting the material performance [4].
In order to achieve the elastic wave attenuation within the low frequency range, heavy tungsten masses are used.Tungsten has a high mass density, as well as further advantages, such as a low environmental threat in comparison with lead, typically used in the literature.The properties of silicone rubber used in this work are based on [39].The compliant nature of the medium has been considered as well as the fact that rubber is an essentially (nearly) incompressible material, characterised by a high ratio between bulk and shear moduli.The material properties of all components of the metamaterial used in the elastic modelling are listed in Table 1.
For the viscoelastic analysis, frequency dependent material parameters are used for the coating constituents.The tungsten cores and epoxy matrix are assumed linear elastic.
The former is motivated by the fact that metals behave in a nearly elastic manner for small stresses and low temperatures [41], whereas the later is based on the results reported in [15].Rubber, which is the material chosen for the coating layers, exhibits relaxation in both bulk and shear properties, as most polymers.However, as pointed out by Lakes [41], the change of the storage shear modulus is orders of magnitude higher than the bulk modulus variation.Also, the loss peak in the case of shear is much broader.Therefore, it is reasonable to assume that the bulk modulus remains real and constant κ = κ, thus frequency independent [42].The viscoelastic behaviour of the rubber is then described by the dependence of the shear modulus G = G(ω) on frequency.In this study, the generalised Maxwell model has been chosen, which adequately captures experimentally observed stress relaxation curves [43].The mechanical analogue of this model consists of a series of dashpots and springs arranged in a parallel manner (Figure 2).The relaxation where n is the number of modes included in the model, G i is a shear modulus contribution, τ i denotes the relaxation time defined as τ i = η i /G i , with η i the dashpot viscosity.For a viscoelastic solid, the viscosity of one of the dashpots is infinite, which results in the constant elastic term G ∞ in equation ( 5).An instantaneous shear modulus is then defined Finally, the complex shear modulus is formulated in the frequency domain as a combination of the storage modulus G and the loss modulus [44], where A useful parameter describing the damping property of a viscoelastic material (here rubber) is the loss tangent given by: which provides a measure of the ratio of energy lost to energy stored in a cyclic deformation [43].
In the present work, the generalised Maxwell model is first applied in the form of a single-mode model, i.e. considering one relaxation time only, in order to study the impact of the frequency of the loss tangent peak on the material performance.The equilibrium modulus G ∞ = 50 kPa has been kept constant, as well as the instantaneous shear modulus G 0 .The assumed maximum loss tangent is approaching tan δ = 0.1, which has been reported as a representative value for silicone rubber [41].Five relaxation times have been considered separately, namely τ i = {10 −1 s, 10 −2 s, 10 −3 s, 10 −4 s, 10 −5 s}, each with the same shear modulus contribution G i = 10 kPa.The frequency dependent material parameters are visualised in Figure 3.Note that with increasing relaxation times, the loss tangent peak is shifted to higher frequencies, whereas the storage modulus growth with increasing frequency is less rapid.

Power transmission spectrum analysis
Since the dispersion diagrams are determined for an infinite periodic material, the wave propagation in a finite size model is further evaluated.To this aim, transmission spectra are typically calculated in order to obtain the amount of wave attenuation.
The finite element method simulations have been performed using COMSOL Multiphysics, following the standard approach given in e.g.[45,46] in the frequency domain.
The model consists of 8 unit cells arranged in a bar with homogeneous matrix material attached at both sides (Figure 5).Boundary conditions assigned to the structure are as follows: in the homogeneous matrix material on the left a unit amplitude sine signal of given frequency is applied by means of imposed displacements (in x and y direction for longitudinal and shear waves, respectively); on the left and right sides so-called lowreflecting boundaries [47] have been used in order to mitigate the reflection from the free boundaries.In the frequency domain, the low-reflecting boundary condition is described by: where n is the normal vector to the boundary and d is a diagonal tensor containing mechanical impedances.The following damping coefficients have been adopted: The point of departure for the power transmission analysis is the time-harmonic averaged energy balance equation, which following Cerveny and Psencik [48] can be written in a general complex form (neglecting body forces) as: where Re(σ : ε * ) the potential energy, Pd = −ω 1 2 Im(σ : ε * ) the dissipated power; v stands for velocity, the bar symbol is a time-harmonic average and an asterisk denotes a complex conjugate.The use of the complex form of the energy balance equation in this study is dictated by: (1) the complex constitutive law related to the viscoelastic model and (2) the frequency domain analysis which results in a complex wave field (the exponential time factor e −iωt has been assumed).By separating the real and imaginary parts of equation (10) one obtains: It can be noted that in the literature, the quantities Re( P) and Im( P) are called active and reactive powers, respectively [48,49,50].Since the active power carries the information about the energy flow [50], the subsequent study is based on equation (11).
In order to assess the performance of a finite size metamaterial, the transmitted active power 1  2 Re(σ • v * • n) spatially integrated over the output section P out (see Figure 5), is compared with the active power integrated over the output section calculated for the reference case of homogeneous epoxy P 0 .Power transmission is therefore given by: In addition, for the purpose of evaluating the contribution of the dissipative mechanisms, the dissipated power for the generalised Maxwell model can be expressed as: spatially integrated over the viscoelastic surface.

Result
In this section, the results of the analysis are presented.Complex band structures have been obtained for the ΓX border of the Brillouin zone (Figure 1) according to the procedure described in section 2.1.In order to assess the performance of the locally resonant acoustic metamaterial, attenuation diagrams have been introduced, where the attenuation factor is given by the angle: which qualitatively indicates the regions of wave attenuation: the full wave propagation occurs for an angle φ = 0.The study is completed by transmission spectra considering longitudinal (x-polarised) and shear (y-polarised) wave propagation separately.

Elastic case
The lossless case is considered as a reference for the viscoelastic metamaterial.In the studies reported in the literature, where due to the simplicity of the viscoelastic material models, the loss tangent was linearly increasing with frequency and the real part of the modulus was constant [14,8], it was sufficient to present only a single elastic reference case.However, adopting the generalised Maxwell model results in a frequency dependence of the storage modulus.This indicates at least two reference cases: linear elasticity with a shear modulus of the rubber coating equal to either the equilibrium modulus G ∞ or the instantaneous modulus G 0 .These cases are marked as elastic G ∞ and elastic G 0 , respectively.
Figure 6 shows the complex band structure and the projections of the dispersion curves on the real and imaginary planes calculated for the case elastic G ∞ .Purely real, purely imaginary and complex bands can be distinguished in the diagram, forming two band gaps in the considered frequency range: the first between 677-815 Hz and the second between 906-1092 Hz.For comparison, the classical real-valued dispersion curves obtained via the ω(k) approach [51] are also depicted in Figure 6 with solid grey lines and show perfect agreement with the real branches of the complex band structure (aside from the horizontal branches associated with torsional modes which are typically badly resolved with the k(ω) approach [15]).It is worth noting that the mechanism opening both band gaps is local resonance of the coaxial heavy masses, which has been thoroughly described by Larabi et al. [5] and Chen et al. [8].This aspect distinguishes these two band gaps from the ones obtained for a single resonating mass [4], where the second band gap is related to the vibration inside the coating layer.Note that the continuity of the branches in the dispersion diagram is preserved: the real branches surrounding the band gaps are connected by bands passing through the complex domain.Moreover, the bands coinciding in the symmetry points are sharply bent.
In Figure 6, a distinction between wave types has been introduced, following Wang et al. [14].The colour varying from red to blue denotes the change of wave polarisation from longitudinal to shear, respectively.Note, that in the dispersion diagrams, a very low level of mode conversion in an infinite material is observed.The longitudinal and shear waves show different behaviour, not only on the real plane projections, where the slope of the curves corresponds to the velocities of the propagating waves, but also inside the The results obtained for the dispersion diagram of the elastic G 0 case are similar (and therefore not shown here) except for a shift of the band gap ranges to higher frequencies, namely 742-893 Hz and 992-1197 Hz.
In Figure 7, attenuation diagrams are presented, where the colours ranging from blue to red are associated with the lowest value of the attenuation factor φ ranging from 0 to π/2, respectively, calculated with equation (15).Ranges where waves have a propagative character are coloured in blue (Im(k x )=0), a deep red colour denotes band gaps with purely imaginary branches (Re(k x )=0 and φ → π/2 ) and light green to red regions are associated with complex shear wave branches passing through the band gaps (Figure 6).
In other words, in the elastic case, the band gap regions can be distinguished by the absence of blue colour, and thus in Figure 7, the mentioned shift of the band gaps to higher frequencies between elastic G ∞ and elastic G 0 cases is clearly depicted.The power transmission spectra used for the evaluation of the results obtained in the previous steps on a finite size structure are presented in Figure 8. Within the band gap ranges predicted by the Bloch approach (shaded regions), transmission dips can be observed for both longitudinal and shear waves.Also the shift of the band gaps between the two elastic cases is again visible.An apparent correspondence exists between the imaginary projection of the dispersion diagram and the transmission spectra: transmission of longitudinal waves exhibits a rapid dip whereas the transmission abatement for shear waves is more gradual within the band gaps.Note that the regions outside the band gaps where the normalised power transmission exceeds 1 and characteristic ripples occur can be associated with the eigenfrequencies of the analysed finite structure.Therefore, their location depends on the chosen size of the sample as well as on the applied boundary conditions, not exactly matching the impedance of the metamaterial.In turn, the location of band gaps is not dependent on the considered set-up, however, the number of unit cells influences the depths of the transmission dips.

Single mode Maxwell model
In Figure 9, the complex band structure as well as its projections to the real and imaginary planes are presented for a selected viscoelastic case with a single (intermediate) relaxation time τ = 10 −3 s.As already indicated in the elastic study, the colours of the curves are associated with wave polarisation.As a reference, the dispersion curves obtained for the elastic case with G 0 are also presented (in black).This choice is motivated by the attenuation regions being located in both cases within the same frequency ranges.
The viscoelastic bands for frequencies approaching the elastic band gap openings are bent out and finally deviate from the reference ones.Moreover, the previously sharp band corners at the symmetry points are significantly smoothed, which is typical for a band structure of viscoelastic metamaterials [14]; this occurs at both band gap regions.It should be noted that material damping mostly affects the two cusps of the imaginary part of the wave number of the longitudinal wave, whereas only a slight decrease is visible in the case of the shear wave.Moreover, the wave attenuation in between the former elastic band gap regions is confirmed by the non-zero values of Im(k x ).It is worth mentioning that in the viscoelastic studies, a band gap in the strict sense does not exist.Therefore, the notion of an attenuation range is used further.
The attenuation diagrams for single mode Maxwell models with different relaxation times are presented in Figure 10.With an increasing relaxation time in the single mode Maxwell model, a shift in the position of the attenuation ranges can be observed.This behaviour, which has been observed also in studies of other PCs and LRAMs [52], can be explained by the variation of the value of the rubber shear modulus G with frequency and has been demonstrated for a one-dimensional case with a single inclusion in [28] based on a discrete Zener-type oscillator.Note that a shift of approximately 100 Hz can be observed for a relatively small variation of storage modulus in the range of 50-60 kPa.In fact, even larger variations of storage modulus, in the frequency range considered, can be encountered in case of silicone rubber as shown in experimental studies [53,22].Since the level of viscosity changes for every frequency, the magnitude of the attenuation angle in between and around band gaps varies.From this point of view, the favourable effect of viscosity is most pronounced in Figure 10b (τ = 10 −4 s), where the attenuation angle between the two attenuation regions exceeds π/6.
Based on the attenuation ranges observed in the analysis of the attenuation spectra, two cases have been selected for a transmission study.The transmission spectrum analysis has been performed for the relaxation times τ = 10 −4 s and τ = 10 −3 s.The results are shown in Figure 11 with a distinction between longitudinal and shear waves.As a reference, the results for the elastic cases are also shown.In line with the predictions from the dispersion analysis, the effect of wave mitigation in the region between the attenuation ranges is more pronounced for the case with a relaxation time τ = 10 −4 s for both wave polarisations.However, this effect is in general more pronounced in the case of a shear wave, which can be explained by the lower group velocity of this wave polarization observed in the region in between band gaps (Figure 6).It has been reported [54] that flat bands are more affected by losses.The material damping contributes also to smoothing of the transmission dips in the longitudinal wave spectrum.In this case, attenuation in the ranges between the elastic band gap regions is present, but rather low in magnitude.On the other hand, the transmission spectrum for the shear wave shows that a low level of material damping does not significantly affect the transmission dips within the band gap regions; in fact, it contributes to the wave attenuation in a much wider frequency range.With increasing loss tangent, wave mitigation around the band gap regions also increases, although, minor smoothing of the transmission dips is observed as well.In both spectra, the shift of the attenuation regions in terms of frequency due to a frequency dependent value of the storage modulus G of rubber is visible.The presented results are in agreement with the predictions from the dispersion analysis.
The dissipative mechanism occurring in the viscoelastic rubber coatings is further investigated using dissipation power spectra.Based on the levels of the dissipated power presented in Figure 12, the viscous dissipation in space can be assessed.Note, that the dissipated power is divided by the output power P 0 obtained for the reference epoxy case, in order to account for the increase of the input power with frequency.First of all, the presence of two maxima in the dissipated power graphs at the band gap openings for the models with both relaxation times and both wave polarisations can be observed.Secondly, most of the dissipation takes place within the first two unit cells, which can be inferred from the small spacing between the plots of the dissipated power for increasing cell numbers.Moreover, it turns out that the total dissipation at the openings of the bang gaps is higher for the viscoelastic case τ = 10 −3 s, even though the values of the loss tangent for this model are lower at the discussed frequencies.This can be explained by higher strain values in the presence of local resonance.A significant level of dissipation can also be noticed in between the former elastic band gap regions for both wave polarisations, with higher values obtained for the viscoelastic case τ = 10 −4 s, which is associated with the effect of joining the attenuation regions.2).A significant shift of the attenuation ranges can be observed for material A compared to B and C, which largely results from the significant difference in the value of G 0 between these materials.In the case of material B (Figure 13b) the favourable effect of joining attenuation ranges occurs (with an attenuation factor exceeding π/6).In Figure 13a, due to the small loss tangent value of material A, the attenuation regions are clearly separated, which can be observed also in Figure 13c, despite the fact that the model of material C differs, with respect to material B, by only one relaxation time (Table 2).

Multi-mode Maxwell model
Therefore, to ensure the effect of bridging the attenuation ranges, the loss tangent of the rubber coatings should at least reach the value of 0.1.Note that the attenuation regions are broadened also from the top and bottom in all instances, thus, the total width of the wave mitigation range increases in comparison with the linear elastic cases.This effect is most pronounced for the metamaterial with a viscoelastic rubber coating of material B.
The power transmission spectra for these materials (Figure 14) are in good correspondence with the predictions given by the dispersion analysis.The location and character of the transmission reduction show both the shifts due to growth of the storage modulus G and the desired attenuation in the region between the elastic transmission dips.Note that in the case of materials B and C, the shift of the band gaps in comparison with material A is more pronounced, due to the stronger difference between the values of the storage moduli G .Furthermore, as demonstrated in the previous section, although the mechanism of joining attenuation ranges in the case of a longitudinal wave can be observed, in terms of the transmission level, it is rather modest.For shear waves, a moderate effect of bridging the transmission dips is present even for materials A and C with rather low loss tangents, and the effect is particularly pronounced for material B that is characterised by a higher material damping.Moreover, material damping for this wave polarization does not significantly affect the band gap regions as for a longitudinal wave, where the transmission dips at the opening of the band gaps range are smoothed.
Based on the dissipated power spectra shown in Figure 15, the losses occurring for leads to a shift of approximately 100 Hz.The loss tangent level, which is nonlinear in the Maxwell model, strongly influences the attenuation performance of the metamaterial around the band gap frequencies.Therefore, taking into account the typical nonlinear frequency dependence of both the storage modulus and the loss tangent of the coating is essential in order to localise the band gap frequencies and correctly assess the influence of material damping.
In the presented analysis, dissipated power spectra demonstrate the efficiency of the proposed viscoelastic metamaterial in terms of energy dissipation.A significant amount of energy is dissipated at frequencies between resonance attenuation ranges which is the underlying mechanism for the effect of band gap bridging.For instance, the dissipation level within the bridged frequency intervals for the metamaterial modelled with the rubber with almost constant loss tangent approaching 0.1 (material B) was only 6 times lower than the highest dissipation peaks.
It has been also demonstrated, that for the purpose of bridging the wave attenuation ranges, the viscosity of the rubber coating needs to be sufficiently high in the frequency range of band gaps formation.Otherwise, the effect of bridging separate attenuation regions is not pronounced.For the investigated case, the loss tangent level should reach at least 0.1, which is a moderate value for silicone rubber.Moreover, the shift of the attenuation ranges should be taken into account as well, if the mentioned bridging is aimed at.
Finally, a material with an optimal value of viscosity should be chosen since material damping may decrease the level of wave attenuation inside the band gap regions (mostly due to the influence on local resonance).The optimal level of viscosity can be found with the aid of the analysis presented in this work by comparing the dissipated power curves for different loss tangents.For overdamped cases, the dissipation level will decrease due to the annihilation of resonance effects.

Summary and conclusion
In this paper, the dynamic behaviour of locally resonant acoustic metamaterials consisting of multicoated coaxial inclusions of rubber and tungsten embedded in epoxy has been investigated.The concept of merging the band gaps by introducing material viscoelasticity to the rubber coating has been analysed in detail for longitudinal and shear wave polarization.For this purpose, the generalized Maxwell model, which accounts for a realistic (nonlinear) frequency dependence of both real and imaginary components of the complex elastic properties, has been adopted.It should be noted that this model, although it is the most general viscoelastic model, has not been yet exploited in the literature for the analysis of locally resonant acoustic metamaterials, and in particular, with the purpose of bridging band gaps.
The analysis presented in this paper shows that the frequency dependent behaviour of the polymer coating material influences the response of the metamaterial in a twofold manner.First of all, the variation of the storage modulus is reported as being responsible for the shift of the band gaps in the frequency spectrum.Secondly, the loss tangent value determines the attenuation performance of the metamaterial.These two observations are crucial if the effect of band gap bridging is aimed at.The shift of attenuation ranges needs to be taken into account in order to adequately estimate the frequency regions where the attenuation is being enhanced with the aid of viscosity.Next, the level of loss tangent in between band gaps needs to be sufficiently high for bridging to occur.Finally, the optimal variation of loss tangent values should be found in order to favourably combine the effect of viscous attenuation in between band gaps with the resonance-based attenuation inside the band gaps.
The findings presented in this paper may help to understand and optimise the performance of locally resonant acoustic metamaterials with viscoelastic constituents.It is clear that more realistic modelling is a necessary step towards tuning metamaterials' behaviour.Moreover, taking advantage of molecular methods for designing polymers with desired properties [37] and controlling their viscoelastic behaviour (e.g.relaxation time) may result in manufacturing of metamaterials with a fully optimised response.Based on the present study and accounting for currently available routes towards active control of material properties, for instance by varying the location of loss tangent peak using electric current or magnetic field [55], some new ways of enhancing metamaterial performance can be envisioned.

Figure 1 :
Figure 1: Geometry of the unit cell (dimensions are in mm) and the reciprocal cell in k-vector space with the shaded area indicating the irreducible Brillouin zone

Figure 2 :
Figure 2: Schematic of generalised Maxwell model

FrequencyFigure 3 :
Figure 3: Frequency dependent parameters of the single-mode Maxwell models with different relaxation times used to simulate a viscoelastic rubber coating: (a) storage modulus, (b) loss tangent.

FrequencyFigure 4 :
Figure 4: Frequency dependent parameters of three viscoelastic materials (A, B and C) described by the five-mode Maxwell model: (a) storage modulus, (b) loss tangent.

Figure 5 :
Figure 5: The numerical set-up for transmission computations

1 ReFigure 6 :
Figure 6: Dispersion diagrams for the linear elastic case with G ∞ (a) 3D band structure (b) 2D projections.The band gap regions are shaded.Colours represent wave polarisations from shear (blue) to longitudinal (red).Grey solid lines depict the real-valued band structure obtained via the ω(k) approach.

Figure 7 :
Figure 7: Attenuation diagrams for the purely elastic metamaterial with different shear moduli of the rubber coating: (a) G ∞ and (b) G 0 .The values of the attenuation angle φ are depicted by colours ranging from blue (0) to red (π/2).

FrequencyFigure 8 :
Figure 8: Normalised power transmission spectra for (a) longitudinal and (b) shear wave polarisation for the elastic cases with G ∞ and G 0 .The corresponding band gap regions are shaded in blue and grey.

1 ReFigure 9 :
Figure 9: Dispersion diagrams for the viscoelastic single mode Maxwell model with τ = 10 −3 s (a) 3D band structure (b) 2D projections.Colours represent wave polarisations from shear (blue) to longitudinal (red).The band gap regions for the elastic G 0 case are shaded and the corresponding dispersion curves are shown in black.

Figure 10 :
Figure 10: Attenuation diagrams for the metamaterial with a viscoelastic rubber coating modelled by single mode Maxwell models with different relaxation times.The values of the attenuation angle (equation (15)) are depicted by colours ranging from blue (0) to red (π/2).

FrequencyFigure 11 :
Figure 11: Normalised power transmission spectra for the metamaterials with a viscoelastic rubber coating modelled by a single mode Maxwell model with two selected relaxation times.The longitudinal (a) and shear (b) wave polarisations are presented separately.In the background, the transmission curves and band gap regions for the two elastic reference cases are shown.

Figure 13 presents
Figure 13 presents diagrams of attenuation angles as a function of frequency for three locally resonant acoustic metamaterials, where the viscoelastic behaviour of the rubber coating is modelled by three different five-mode Maxwell models: materials A, B and C (Table2).A significant shift of the attenuation ranges can be observed for material A

FrequencyFigure 12 :
Figure 12: Dissipated power calculated for subsequent pairs of unit cells (starting from the excitation side) for longitudinal (a) and shear (b) wave polarisations; green and blue colours denote the viscoelastic case τ = 10 −3 s and τ = 10 −4 s, respectively.

Figure 13 :FrequencyFigure 14 :
Figure 13: Attenuation diagrams for the metamaterial with a viscoelastic rubber coating of material A (a), material B (b) and material C (c) modelled by multi-mode Maxwell models.The values of the attenuation angle (equation (15)) are depicted by colours ranging from blue (0) to red (π/2).

Table 1 :
Elastic material properties of the metamaterial components (ρ denotes density, G shear modulus, K bulk modulus, c l and c s longitudinal and shear wave velocities, respectively)

Table 2 :
Relaxation times τ i and associated shear moduli G i of three viscoelastic materials (A, B and C) described by the five-mode Maxwell model