Wave attenuation and trapping in 3D printed cantilever-in-mass metamaterials with spatially correlated variability

Additive manufacturing has become a fundamental tool to fabricate and experimentally investigate mechanical metamaterials and phononic crystals. However, this manufacturing process produces spatially correlated variability that breaks the translational periodicity, which might compromise the wave propagation performance of metamaterials. We demonstrate that the vibration attenuation profile is strictly related to the spatial profile of the variability, and that there exists an optimal disorder degree below which the attenuation bandwidth widens; for high disorder levels, the band gap mistuning annihilates the overall attenuation. The variability also induces a spatially variant locally resonant band gap that progressively slow down the group velocity until an almost zero value, giving rise to wave trapping effect near the lower band gap boundary. Inspired by this wave trapping phenomenon, a rainbow metamaterial with linear spatial-frequency trapping is also proposed, which have potential applications in energy harvesting, spatial wave filtering and non-destructive evaluation at low frequency. This report provides a deeper understanding of the differences between numerical simulations using nominal designed properties and experimental analysis of metamaterials constructed in 3D printing. These analysis and results may extend to phononic crystals and other periodic systems to investigate their wave and dynamic performance as well as robustness under variability.


FE numerical validation of the wave trapping
The numerical FRF results as a function of space for the samples MM1 and MM4, which were used to validate the experimental wave trapping phenomenon induced by the variability, are shown in Figure S2. In these simulations, the measured material properties were included in the FE model to compute the dynamic response. These numerical results are in agreement with the experimental observations (see Figure 5 and Figure 6), which validates the method for spatially varying material property estimation using the cube specimens. Also, the numerical model confirms that for the MM4 sample the wave trapping is created for forward wave propagation (i.e., excitation at interface 0) and only wave attenuation is produced for backward wave propagation (i.e., excitation at interface 15). Figure S2. Numerical FRF as a function of space for the samples: MM1 (a, c) and MM4 (b,d). Excitation location at interface 0 (a-b) and at interface 15 (c-d). Legend: the colors represent the FRF magnitude in dB.

More about the rainbow metamaterial
In this section, more results about the proposed rainbow metamaterial (E(n) = E 0 [1 + α(n ∆ − L/2)] with α = 0.5) are presented. The FRF from interface 15 to interface 0 or vice-verse (i.e., G 15,0 = u 0→15 / f 15 , which is equal to G 0,15 = u 15→0 / f 0 ) is presented in Figure S3(a), where the resonance peaks inside the vibration attenuation zone create the "jagged" profile. In Figure S3(b), the FRFs as a function of space are presented for excitation at interface 15, as discussed in the paper, the wave trapping isn't created even with spatially varying band gap and only wave attenuation is observed because the group velocity doesn't slow down until an almost zero value at the upper band gap boundary. Finally, by exciting the system at interface 15 with the same angular frequencies from Figure 7(e-h), the backward wave propagation isn't trapped at the upper band gap boundary and only vibration attenuation is observed in the displacement fields, Figure S3(d-g).

Karhunem-Loeve expansion
A random field H(x, p) can be defined as a collection of random variables indexed by a continuous parameter x ∈ D, where D describes the system domain. In other words, for a given position x 0 , H(x 0 , p) is a random variable, and for a given outcome p, H(x, p) is a realization of the field. There are several methods available in the literature for generating random fields S1 and the Karhunem-Loeve (KL) expansion is a special case of series expansion using random variables and deterministic spatial functions, which are orthogonal and derived from the covariance function. A homogeneous random field with a finite, symmetric and positive definite covariance function, defined over a domain is given by S2 where H 0 is the random field mean value, X j are random variables, l j and f j (x) are eigenvalues and eigenfunctions, solutions of the Fredholm integral equation of the second kind from the covariance function S2 and the p dependence has been omitted. The eigenvalues and eigenfunctions can be ordered in descending order of magnitude of the eigenvalues to truncate the series in Equation (S1) to a finite number of terms N KL , chosen by the accuracy of the series in representing the covariance function, rather than the number of random variables S3 . As a rule of thumb, N KL can be chosen such that l j /l 1 < 0.1, and it depends on the correlation length of the random field. The longer the correlation length the more rapidly the eigenvalues decrease, meaning that fewer terms are needed to accurately represent the series.
If H(x) is a Gaussian random field, X j are always independent zero mean, unit standard deviation Gaussian random variables. However, if the random field is not Gaussian, it is not possible to use the KL expansion to directly generate H(x) because X j have unknown joint PDF. Some approaches have been proposed to overcome this issue, amongst them an iterative scheme S4, S5 , which uses directly the KL expansion and can simulate both stationary and non-stationary random fields. Moreover, if the target CDF is approximately Gaussian, only one iteration might be enough to achieve convergence. In general, the eigenproblem can only be solved numerically and normally involves some procedure for discretizing the random field S1, S6 . However, for some families of correlation functions and specific geometries, there exist analytical solutions of this integral equation S2 .
For the Monte Carlo analysis, the material properties of the unit cell (ξ = E, ρ, ν) were chosen to be randomly varying according to a random field, i.e., ξ (x) = H(x), where H(x) is a Gamma distributed homogeneous random field with correlation function given by where τ is lag or the distance between any two points in the random field, and b is the correlation length. Then, an analytical solution for Equation (S1) is available S2 . The properties within each cell are constant, which is usually accurate for large correlation lengths b > 0.6L S7 , and given by H(x c ), where x c is the center of each cell. The Gamma distribution is chosen according to a Maximum Entropy criterion S8 , which is given by where Γ(x) = +∞ 0 t x−1 e −t dt, x > 0 is the Gamma function and the parameters of the distribution are given by a 0 = 1/δ 2 ξ and b 0 = ξ 0 δ 2 ξ , where δ ξ is the dispersion parameter. Equation (S1) is used to simulate the non-Gaussian random field S4, S5 in which the X j are firstly generated using Gamma independent random variables and the Latin Hypercube Sampling (LHS) scheme as a stochastic solver S9 .
By knowing the experimental spatial distribution of material properties, Figure S1, the correlation function and length of Equation (S2) can be estimated, which can be used to estimate numerically the material property distributions and, hence, to provide a significant physical insight in the stochastic analysis of the ensemble of metastructures.

4/5
The experimental mass density and elastic modulus are highly correlated, as presented in the scatter plot of Figure S4(a-b). Due to the manufacturing process, they depend on a single parameter, which is the porosity. Moreover, the Poisson coefficient is only slightly correlated with these two parameters, Figure S4(a,c). By using Equation (S2) as a proposed correlation function, the correlation length was estimated as b = 0.22L for all the three material parameters, as seen in Figure S4(d). The experimental correlation was calculated by using an unbiased FFT-based estimator. Because the correlation length is larger than three times the unit cell length (a/b ≈ 3.3), the material properties inside the unit cell can be considered constant. The correlation function given by Equation (S2) with b = 0.22L, and the experimental standard deviations for mass density, elastic modulus and Poisson coefficient were applied on the KL expansion to generate numerical spatial profiles of material properties, and 10 random samples are shown Figure S5.

Comparison between experimental and numerical FRFs
Finally, the dynamic response of the metamaterial beam samples are shown in Figure S6 for all the 10 experimental measurements as well as for the FE simulation by using the KL estimation for the 10 material distributions presented in Figure S5. A similar dynamic behavior can be observed, which is also shown in the validation and discussions of Figure 8(a-b), where 250 material distributions were used for the sake of convergence in the Monte Carlo simulation.