Breaking the mass law for broadband sound insulation through strongly nonlinear interactions

Sound transmission through panels is governed by the well-known mass law in the mid-frequency range. This paper reveals a possibility of breaking this density-dominant law through strongly nonlinear interaction, while broadening the bandwidth for effective sound insulation. For this purpose, a basic model is established, and corresponding exact analytical methods for bifurcation and stability analyses are proposed. Influences of four typical types of nonlinear interactions on the wave insulation are analytically and numerically investigated. We find that, by introducing strongly nonlinear interactions at appropriate locations, the nonlinear model can not only break the barrier imposed by the mass law, but also entails broadband sound insulation by 2–3 times relative to the optimal linear model. Meanwhile, the sound insulation valley due to the coincident effects can also be eliminated. With bifurcation and effective mass, we clarify that the enhanced wave insulation of the strongly nonlinear models arises from the broader band of super mass induced by strongly nonlinear local resonances, which depends on the bifurcation of periodic solutions. The proposed models and the findings provide a solid basis and new possibilities for wave insulation in complex nonlinear structures and nonlinear acoustic metamaterials.


Introduction
Insulation of low-frequency and broadband sound and elastic wave is of great significance for a wide range of applications, such as protecting buildings from earthquake waves, reducing noise of airplane cabin (insulating fluid and aeroengine noise) [1,2], or conceiving various wave manipulation devices. It is well known that sound insulation of a monolayer plate/wall generally follow the mass-density law in the mid-frequency range beyond the stiffness-dominant region and before the coincidence frequency [3], i.e. the transmission loss T L = 20log(ωM/2ρc), where ω denotes the angle frequency; M the surface density of the wall; ρ and c the medium density and sound speed, respectively. This law stipulates that a heavy panel is needed to achieve good sound insulation. More specifically, improving the sound insulation by 6 dB requires doubling the surface density of the panel. For a finite plate with boundary, the insulation of very low-frequency sound depends on the structural stiffness. Between the stiffness-dominant band and the density-dominant band, a coincidence frequency appears, where the transmission loss becomes zero. It is an eternal pursuit to find ways to effectively insulate low-frequency and broadband sound with small weight cost. Many linear ways have been attempted, but the challenge remains [4].
From a different perspective and analogous to the mechanical nonlinear energy sink (NES) [5,6], Bellet et al [7] proposed a nonlinear acoustic absorber that can effectively dissipate the resonant or transient response of an acoustic cavity or tube. Such nonlinear acoustic absorber can be made of a viscoelastic membrane or loudspeaker attaching to an orifice of a cavity [7][8][9][10][11][12]. The problem can be simplified to a two-degrees-of-freedom (2DoF) nonlinear mechanical system [8]. However, the way to insulate acoustic waves (instead of absorbing waves) using NES has not been clearly understood. Figure 1(a) shows a plate (membrane or wall) with surface-mounted local resonators. The plate is immersed in acoustic media I and II. While sound insulation of LAMs can be well apprehended with frequency-domain finite element methods, the nonlinear system under investigation requires special treatment for analysis. As shown in figure 1(b), to find nonlinear ways that may overcome the 'narrow band' and 'concomitant insulation valley' disadvantages of acoustic metamaterial, we establish an analogical discrete model that can not only reflect the essential wave insulation properties of the nonlinear systems, but also can be efficiently solved in time domain as well as with approximate analytical solutions. The model consists of a linear phononic crystal I, a 3DoF nonlinear system, another linear phononic crystal II and a perfect match layer (PML). Phononic crystals I and II can simulate the low-frequency wave propagation (much lower than the bandgap frequency) in the homogenous media I and II, respectively (as further explained below). The nonlinear system, inserted between phononic crystals I and II, consists of a primary mass m B , two coupled local resonators m 1 and m 2 , and a boundary spring. The primary mass m B models the wall for wave/sound insulation, and the spring k B models the influence of the boundary stiffness. Therefore, this model can represent either a NES or a metamaterial design, because one generally considers a periodic unit cell of metamaterial.

Model
The phononic crystals I and II comprise periodic mass, m 0 , and linear springs, k 0 . The motion equation of the nth cell writes:  where u n denotes the displacement of the nth m 0 . Replacing k 0 with k 0I and k 0II distinguishes crystals I and II. By introducing the Bloch theorem, we obtain the dispersion equation (see figure 2) of this mono-atomic chain as: in which q symbolizes the wave vector and ω = 2πf denotes the angular frequency. At q = π, one gets the initial frequency of bandgap of the chain, f gap = k 0 /m 0 /π . The wave speed, c = dω/ dq, can be solved by equation (2). At q = 0, long-wave speed is c 0 = k 0 / m 0 . For the low-frequency wave of frequency f ≪ f gap , the wave speed decreases only slightly with respect to c 0 and can be considered constant. Therefore, at frequencies much lower than the bandgap frequency, phononic crystals I and II can simulate low-frequency wave propagation with constant wave speed in continuous media I and II, respectively. In this case, equation (1) can be regarded as a difference equation. Furthermore, u is a function of propagation distance, x, and time, t. Therefore, equation (1) can be rewritten as a partial differential equation: which is the linear sound equation in a homogeneous medium like air. Due to dispersive effects, note that the above equivalent process between equations (1) and (3) is valid only for f ≪ f gap . In the present case, if we define a limiting frequency f eq in order to comply with c(f eq ) = 0.8c 0 , the derivation dω/ dq yields f eq ≈ 0.61f gap .
In the 3DoF nonlinear system, the primary mass m B connects to the foundation through the nonlinear spring k B u B + k BN u 3 B , where k B and k BN denote the linear and nonlinear stiffness coefficients, respectively; and u B denotes the displacement of m B . m 1 is coupled to m B through the nonlinear force 3 , and m 2 interacts with m 1 by the nonlinear force k 2 (z − y) + k 2N (z − y) 3 , where k i and k iN denote the linear and nonlinear stiffness coefficient, respectively, i = 1, 2; and y and z symbolize the displacement of m 1 and m 2 , respectively. Therefore, m 2 is indirectly coupled to m B . Only the typical cubic nonlinearity is considered here. The equations of motion read: in which, F denotes the total force from the three boundaries. In the numerical model, where u LN is the displacement of the last m 0 in crystal I; u R1 is the displacement of the first m 0 of crystal II. An incident wave is emitted from the left terminal of the chain I. Let us define the incident, reflected and transmitted wave fields as u i (x, t), u r (x, t), u t (x, t), respectively. For a continuous system, taking equation (5) as a difference equation at the boundary x = 0 yields: This equation describes the interaction between the wave field and the wall structure.
The following parameters are used in the subsequent numerical analyses: . Moreover, f gap = 318 Hz and f eq ≈ 195 Hz. The incident wave amplitude is u ia = 40 µm. Nonlinearity is tuned by changing the nonlinear stiffness coefficients. We specify k 1 > k 2 to ensure effective energy transfer from m B to m 2 .

Numerical methods
Wave propagation in this chain can be solved by time-domain numerical integral method, but we must model the non-reflection boundary. For that, the length of phononic crystal I should be longer than c 0 × T, where T is the total time length of the incident wave packet. Here, phononic crystal I consists of 180 cells, and crystal II has 60 cells. Crystal II can be short because an optimized perfect match layer (PML) is connected to it. The PML consists of 120 periodic cells with damping. The motion equation of PML oscillators is Here, the damping coefficient d 0 = 0.0015. Moreover, to eliminate the influence of the boundary reflection and to distinguish the incident wave and reflection wave in crystal I, we use a ten-period sinusoidal wave packet as the input signal. The numerical integral time is 0.8 s. Typical wave forms are shown in figure 3, in which different wave components can be clearly identified. Furthermore, we extract all peaks of the transmitted wave packet, and use their average value u tav to calculate the transmission loss T L = 20log 10 (u ia /u tav ).

Analytical method for linear model
For linear model, when k BN = k 1N = k 2N = 0, we adopt the effective mass of the metacell to calculate the wave insulation. In this case, taking force F in equation (4) as an unknown variable, the effective mass of the linear model reads: where The total mass m t = m B + m 1 + m 2 . Local resonators will not only present the negative mass m eff,L < 0, but also the unusual bands for |m eff,L /m B | > 1 or |m eff,L /m t | > 1, which are defined as the 'super mass' relative to the primary mass m B and the total mass m t , respectively. Super mass is paramount to wave insulation. Then, replacing F with m eff,LüB in equation (7) gives the interaction equation: To solve this equation, let us specify By substituting these variables into equation (8), one has, Moreover, based on the continuity of displacement at the incident and transmission boundaries, one obtains boundary conditions: Here we consider the case when q 1 = q 2 = q. Combing the two equations gives the transmission loss: Moreover, by replacing m eff,L with m B or m t , we can study the influences of mass and boundary stiffness k B on the wave insulation, which is termed as mass law. Equation (12) indicates that T L = 0 dB for k B = ω 2 m eff,L . This frequency is termed as the coincident frequency, f co .
Furthermore, considering wave propagation in a 'medium' consisting of the periodic 3DoF cells in figure 1 The solution of equation (13) denotes u = u a e i(ωt±q eff x) , where q eff denotes the equivalent wave vector of such medium, If the imaginary part |Im(q eff )| > 0, wave attenuates quickly, thus giving the corresponding frequency range referred to as bandgap.

Analytical method for nonlinear model
For the nonlinear model, we specify g = y − u B , s = z − y to simplify the equation (4) as It is challenging to solve these nonlinear equations with boundary conditions. We propose a complete nonlinear method based on the first-order harmonic balance method.

Complete nonlinear method
First, we can directly solve the interaction equations by specifying the variables as below: High-order harmonics are not taken into consideration. By substituting equation (16) into equation (15), one obtains six equations with the first-order harmonic balance, .
Moreover, by substituting equation (16) into equation (6) and specifying q 1 = q 2 = q, one obtains: Based on the boundary conditions described by equation (11), one has: Substituting equation (19) into equation (18) yields At last, replacing F 1 and F 2 in equation (20) by equation (17), we can solve the six unknown variables A, B, G 1 , G 2 , S 1 , S 2 for any given incident wave amplitude u ia . Specifying we can re-write the algebraic equations into matrices as: in which, .
The equations can be analytically solved. We note that there is only one real number but several complex number solutions. Bifurcations occur when different branches of periodic solutions cross. Given the fact that these periodic solutions rely on stability analysis [46], we adopt the theory proposed by Fang et al in [47] to analyze the stability of a solution based on harmonic balance method. The Jacobian matrix J for a periodic motion is In which the notation [[x]] stands for the square matrix built upon the vector x; the expression A • B denotes the element product of two matrices A and B. If the real part of an eigenvalue of J is positive, the corresponding solution is unstable. Otherwise, the solution is stable. Note that a real number solution of equation (21) may be unstable, and a complex solution may be stable.
Furthermore, the transmission loss is: For linear cases, results from equations (23) and (12) are the same. The nonlinear effective mass is As confirmed, for linear systems, solutions from equations (23) and (12) are equal.

Equivalent linear method
In addition, because the stiffness coefficients between the oscillators m B , m 1 , m 2 are tuned by nonlinear interactions, we can build an equivalent linear model of the nonlinear model, i.e. the equivalent linearized approach. For this, supposing the equivalent linear stiffness are k 1eq = γ 1 k 1 , k 2eq = γ 2 k 2 , k Beq = γ 3 k B , where γ i is an approximate factor. By substituting k i with k ieq in equation (7), one gets the effective mass of the equivalent linearized system as: Then we can obtain the wave insulation of the linearized system with equation (12).

Wave insulation of linear model
Before investigating nonlinear models, we recall the properties of wave insulation of linear structural or metamaterials, as shown in figure 4. Two other cases, m eff,L = m B and m eff,L = m t , are also depicted in figure 4. Numerical results agree with the analytical curves, which confirms the validity of analytical method. There are three coincident frequencies, f co1 < f co2 < f co3 . Owing to the boundary stiffness k B ̸ = 0, waves in 0 ∼ f co1 are greatly insulated because |Im(q eff )| > 0 in 0 ∼ f co1 , as indicated in figure 4(b). This is the stiffness-dominant band. The range f > f co1 is termed as mass-dominant band, in which T L increases by 6 dB when the mass is doubled. Two insulation peaks appear near the locally resonant frequencies, f r1 and f r2 , and the peaks exceed the mass law limit corresponding to m eff,L = m t . Figure 4(c) and equation (12) indicate that the high transmission loss appears in the band for k B − ω 2 m eff,L >> 2k 0 ω/ c 0 , which means m eff,L → ±∞. Therefore, the insulation peaks near f r1 and f r2 are generated by the mass amplification arising from the local resonances. However, as stated in introduction, in linear model, every narrow insulation peak is accompanied by a broad insulation valley near the coincident frequency. T L in valleys is much smaller than curves for m eff,L = m B and m eff,L = m t .

Cases Parameters
Case A: remote nonlinear local resonator

Wave insulation of nonlinear models
Four types of nonlinear interaction are considered, where nonlinearity is imposed at different positions. As listed in table 1, they are distinguished by the nonlinear boundary, remote nonlinear local resonator (nonlinearity at m 2 ), close nonlinear local resonator (nonlinearity at m 1 ), and completely nonlinear model (nonlinearity acts at all positions).

Case A: remote nonlinear local resonator
The bifurcation analyses of model A with typical variables, k 2N = 3 × 10 10 and k 2N = 1 × 10 13 N m −3 , are shown in figure 5, alongside the wave insulation properties shown in figure 6. We combine figures 5 and 6 to clarify the phenomena and the underlying physical mechanisms. Similar figures are presented for other cases. As shown in figure 5, there is a saddle node bifurcation point, marked by f B , in the wave insulation curve T L and the generalized effective mass |m eff /m t | curve. For f < f B , there are two analytical branches: the branch with real number solution is stable (red circles); the other complex solution is unstable (black circles). Therefore, the numerical result will follow the first stable branch for f < f B . For f > f B , there are three branches, 1, 2, 3. The real number solution on branch 1 becomes unstable, while the complex branch 2 is stable. Therefore, there may be a jump from branch 1 to branch 2 at f = f B .
For clarity, besides the numerical results, only the real number analytical solutions are present in figure 6. Analytical curves agree well with numerical results. The numerical responses mainly follow with the real number solution. The influences of the nonlinear coefficient k 2N on the transmission loss are shown in figure 6(b). Changing k 2N first alters the frequency of the resonator m 2 , whose frequency is lower than m 1 . When increasing k 2N from 1 × 10 10 to 1 × 10 13 N m −3 , the first insulation peak always exists but its frequency increases gradually from 44 to 66 Hz. The second insulation peak first shifts a little to higher frequencies, but a moderate nonlinearity (k 2N = 3 × 10 10 ) changes the trends. For k 2N = 3 × 10 10 , the insulation peak near 180 Hz, as predicted by analytical solutions, is not generated in numerical T L because this peak corresponding to the unstable solution in figure 5(a). Instead, as predicted by the bifurcation jump in figure 5(a), the numerical T L features a jump near 140 Hz. As shown in figure 6(a), further increasing k 2N > 3 × 10 10 will diminish the second insulation peak until it completely disappears, so that only a low frequency insulation peak is left. The equivalent linearized theory shows that this peak is formed when k 2 → +∞, indicating that it happens because enhancing the hardening nonlinearity merges the two  Moreover, the frequency ranges, relating to the shifting and disappearing of the wave insulation peaks in figures 5(a), (c) and 6, are in accordance with the varying trend of the super mass |m eff /m t | > 1 in figures 5(b) and (d). Except for the insulation in the stiffness dominant band for f < f co1 = 20 Hz here, we can always obtain T L ∝ |m eff /m t | along the stable branches. Therefore, the amplified super mass due to local resonances remains the essential mechanism for the high sound insulation of nonlinear model A.
However, in Case A, the nonlinear interaction between m 1 and m 2 cannot offer broader-band wave insulation than the equivalent linearized model (see the black curve in figure 6(a)). This means the nonlinear model A cannot provide better insulation capability than a proper linear system with same mass. Moreover, the spectra of the transmitted waves at the insulation peaks in figure 7 show that many high-order harmonics are generated by strong nonlinearity.

Case B: close nonlinear local resonator
In Case B, nonlinearity only acts between m 1 and m B , k 1N > 0, k 2N = 0. We combine the bifurcation analysis in figure 8 and numerical results in figure 9 to show the physical phenomena. As shown in figure 9(b), the first insulation peak nearly persists when increasing k 1N . The second insulation peak shifts upwards, but will not disappear as what happens in Case A. As shown in figure 8, the model gives two or three analytical  solutions for each given frequency, and the bifurcation depends on k 1N . For k 1N = 1 × 10 11 , all three branches after the bifurcation frequency f B are stable, but two of them lead to negative insulation T L that will not occur in practice. Therefore, as demonstrated by the numerical result in figure 9(a), the practical curve follows the highest T L with the real number solution. For k 1N = 1 × 10 12 , the real number solution in f < f B is stable, and the middle branch in f > f B is stable, indicating a jump at f = f B . For k 1N = 1 × 10 13 , the second peak is shifted to high frequency; the curve with higher T L is stable in most part of the frequency range, with T L increasing with the frequency in 60-220 Hz. This is also demonstrated in in figure 9(a). Still on this curve, stable and unstable solutions appear alternately near 170 Hz, so the numerical results fluctuate. Therefore, numerical results still agree well with analytical curves with stable and real number solutions.
More importantly, the solutions for both the first and the second insulation peaks remain stable in model B (the second peak for model A is unstable). Therefore, model B can provide two insulation peaks. Interestingly, both the bandwidth and transmission loss of the second peak are greatly increased along with the increasing k 1N . Meanwhile, the coincidence effect behind the second peak also disappears even with weak nonlinearity, k 1N = 1 × 10 11 . This happens because of the bifurcation in figures 8(a) and (b). As a result, T L in the insulation valley increases by more than 10 dB under strong nonlinearity. Furthermore, we check whether an equivalent linearized model can realize this effect. As shown in figures 8(c), (d) and 9(a), we tune k 1eq = 2k 1 to make the equivalent linear model approximately gets the same initial frequency for T L > 10 dB near the second peak of the nonlinear model k 1N = 1 × 10 12 . The band for T L > 10 dB of the linear model is 125-168 Hz, while it is 120-200 Hz for the nonlinear model with k 1N = 1 × 10 12 , 1.86 times broader. Therefore, the nonlinear model B can eliminate a coincident effect and present broader wave insulation with the same mass.  Furthermore, we analyze the responses at 180 Hz and the effective mass for k 1N = 1 × 10 12 , as shown in figure 10. The waveforms and spectra suggest that the incident and transmitted waves have the same main frequency, and high-order harmonics are not remarkable. Therefore, the frequency is nearly preserved. The high-order harmonics in the transmitted is suppressed because their frequencies locate in the bandgap of phononic I and II, i.e. 3f > f gap . As illustrated in figure 8, the peaks of |m eff /m t | are consistent with the insulation peaks, and their varying trends are also same when enhancing nonlinearity. Therefore, the broader and greater wave insulation of nonlinear model B stems from the broader band of super mass induced by strongly nonlinear local resonances. The mechanism is consistent with model A.

Case C: both close and remote nonlinear local resonators
In this section, we consider the nonlinear interaction appears at both the remote and close local resonators, and k 1N = k 2N > 0 is specified here. As shown in figure 11, this model has seven analytical branches, which result in many bifurcations. Fortunately, only the branch featuring the first insulation peak is stable in f < f B ; and only one branch stemming from the saddle node bifurcation point f B is stable, i.e. jumping happens here. As shown in figure 12, analytical curves agree with numerical solutions very well, including the jumping at f B = 162 Hz for k 1N = k 2N = 5 × 10 10 .
In whole ( figure 12(b)), while increasing the nonlinear coefficient from 1 × 10 10 to 1 × 10 13 , the first insulation peak gradually shifts from 47 Hz to 120 Hz; the second insulation peak shifts upwards first, but diminishes when k 1N = k 2N > 1 × 10 11 . The second peak diminishes because the periodic solution becomes unstable at this peak, as illustrated in figure 11(a), like in Case A. Properties of cases A and C indicate that the second peak will be suppressed under moderate nonlinearity if there is nonlinearity in resonator m 2 .
Under strong nonlinearity (descripted by k 1N = k 2N = 1 × 10 13 or 1 × 10 14 ), though there is only one insulation peak, the nonlinear model offers very broad and extremely high wave insulation that breaks through the mass law for m eff = m t . Taking k 1N = k 2N = 1 × 10 13 as example, the equivalent linearized model has T L > 10 dB in 80-126 Hz, but it is 72-200 Hz for the nonlinear model. The generalized bandwidth (f end /f start − 1) is broadened by three times by strong nonlinearity and the insulation valley is also eliminated. Solutions for k 1N = k 2N = 1 × 10 14 in figure 12(c) validate the enhanced wave insulation.    In figure 12(a), numerical curve for k 1N = k 2N = 1 × 10 13 jumps near 116 Hz. As shown in figure 12(d), for f = 115 Hz before the jump, the nonlinear model generates remarkable third-order responses. When the fundamental wave transmits into phononic crystal II, the third harmonic is localized inside the nonlinear resonators because 3f locates in the bandgap of the crystals. Behind the jump (117 Hz), the third harmonic in the oscillators is not as large as the one at 115 Hz. The main difference between the two transmitted waves lies in the responses in 200-300 Hz.
Furthermore, we increase the wave speed c 0 to 1500 by increasing k 0 to study the influences of wave speed. As shown in figure 13, a broadband wave insulation of strongly nonlinear model is achieved again. Theses analyses confirm the generality of the enhanced insulation capability of the strongly nonlinear model which embed proper nonlinear interactions. Moreover, a jumping of T L appears at the insulation peak again.

Case D: completely nonlinear model
At last, we study the model D with nonlinearity in two local resonator and the nonlinear boundary. Based on the results above, we specify k 1N = k 2N = 1 × 10 13 and change k BN in case D. Bifurcation analysis is not repeated here. As illustrated in figure 14, the influences of boundary nonlinearity mainly lie in the stiffness-dominant region and it has little effect on the insulation peak near 120 Hz. Increasing k BN can greatly expand the stiffness dominant range. In cases B and C, only the coincident frequencies behind the resonant insulation peaks are eliminated by nonlinearity. Here, the coincident insulation valley between the stiffness and mass dominant bands are eliminated when k BN is large. As predicted by the analytical method, T L curve undergoes a drastic change at the junction between the two bands. The numerical results demonstrate these properties. However, compared with the equivalent linearized system by changing k B , T L induced by nonlinear model is smaller than the one in the linear model for f → 0.  Moreover, numerical curve T L in figure 14 jumps from 6 dB to 27 dB in 80-82 Hz. We analyze the waveforms and spectra for before (75 Hz) and after (80 Hz) the jump, as shown in figure 15. At 75 Hz, there are high-order or even chaotic responses but the three oscillators, m B , m 1 and m 2 , vibrate almost synchronously and their amplitudes are approximately equal, which means the local resonance in m 1 and m 2 is not formed. At 80 Hz, the motions are still synchronous but energy is localized in m 1 and m 2 . Therefore, the great insulation is generated by nonlinear local resonance.

Comparison of different nonlinear cases
Herein we briefly summarize the influences of different nonlinear interactions on the acoustic wave insulation. The optimal results for cases A, B and C are illustrated in figure 16: (i) When the nonlinear interaction only appears between m 1 and m 2 (case A), the wave insulation is not drastically improved with respect to the linear model. (ii) When strongly nonlinear interaction only appears between m B and m 1 (case B), the low-frequency insulation peak remains; the coincident insulation valley after the second insulation peak is eliminated; the bandwidth for high insulation is expanded by two times. (iii) When there are strongly nonlinear interactions in both pairs m B -m 1 and m 1 -m 2 (case C), the bandwidth for high insulation can be broadened by about three times. (iv) When strong boundary nonlinearity is introduced (case D), the coincident insulation valley between the stiffness and mass dominant bands are eliminated, but the transmission loss for f → 0 becomes smaller than the one offered by equivalent linear model.
(v) When nonlinearity exists at the boundary and two resonators, low-frequency and broadband wave insulation is achievable.

Conclusions
This paper reveals a possibility of breaking the mass law of wave insulation through exploring system nonlinearity, while broadening the frequency bandwidth for effective sound insulation. As illustration, a basic and typical 1D model, consisting of two local resonators and elastic boundaries, is established to simulate the acoustic/elastic wave insulation. Three nonlinear sources, placed in the boundary and two resonators, are respectively considered. Analytical methods for both linear and nonlinear model are proposed, and their accuracy is validated by numerical integration method. Bifurcation and stability of periodic solutions are performed. It is shown that the wave insulation efficacy closely depends on the location of nonlinearity. Influences of four typical types of nonlinear interactions are investigated with the underlying mechanisms explained. Nonlinear interaction changes the wave insulation properties in different ways. By introducing proper strongly nonlinear interactions, especially strong nonlinearity between the primary mass and the first resonator attached on it, the bandwidth for superior wave insulation can be expanded by 2-3 times relative to the optimal linear model with the same mass (the equivalent linearized model). Meanwhile, the insulation valley arising from the coincident effects can also be eliminated. Moreover, the mechanism underpinning the observed broadened bandwidth for large wave insulation is elucidated by combining the bifurcations, stability, effective mass and equivalent linearized model. We establish that the observed broader and superior wave insulation offered by the strongly nonlinear model stems from the broader band of super mass induced by strongly nonlinear local resonances.
This study provides theoretical foundation for acoustic/elastic wave insulation in complex nonlinear structures, which can be extended to the metamaterial. It sheds lights on the underlying physics while offering a new way to design high performance light-weight structures for acoustic wave insulation.

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.