Frequency stop-band optimization in micro-slit resonant metamaterials

An optimized unit cell design of a micro-slit resonant metamaterial is proposed to increase the size of the frequency stop-bands and to enhance the sound absorption at normal incidence. Micro-slit resonant metamaterials offer a compact and lightweight solution for low-frequency noise reduction, in contrast to traditional methods such as absorptive foams. A combination of numerical and semi-analytical solutions based on dispersion and absorption curves is presented. A novel algorithm allows for the decoupling of wave types from raw numerical data in the dispersion curves, without using the stiffness and mass matrices. A thorough optimization process of unit cell designs with genetic algorithms is performed. Focus is given to the ﬁrst frequency stop-band, located in the frequency range of application of micro-slit resonant metamaterials. The process shows that relatively large resonators (with respect to the unit cell total area) produce a larger ﬁrst frequency stop-band, whereas slit size has a negligible effect. The optimized design increases the ﬁrst frequency stop-band by 20% and the second stop-band by 25% compared to the literature standard. The absorption curves at normal incidence of acoustic waves are derived numerically for a rigid and elastic frame of the metamaterial backed by a cavity. These curves are validated by the JCAPL semi-analytical model. The optimized unit cell design shows a 9% increase in the ﬁrst peak of absorption coefﬁcient compared to the literature standard at a cavity depth of 30 mm, and an increase of 10% at a cavity of size 53 mm. Stop-band behavior does not inﬂuence sound absorption at normal incidence of acoustic waves in the frequency range of interest. (cid:1) 2021 The Author(s). Published by Elsevier Ltd. This is an open access article under the CC BY license (http:


Introduction
Acoustic noise reduction has been a topic of interest in the scientific community for years.As technology evolves, novel solutions are emerging to deal with problems as soundproofing of recording studios, aero-engine noise reduction, or preventing the propagation of structure-borne sound to protect high-tech machinery.For the past decades, heavy materials such as absorptive foams or porous materials have been the most common choice.These materials are proven effective for wavelengths up to a quarter of the thickness of the material [1].For this reason, foams require a large thickness to achieve low-frequency noise reduction.In the architectural sector, the impact of weight and thickness of the absorptive materials is not a primary concern.However, for highspeed trains or aircraft, designs are required to be compact and lightweight.Acoustic metamaterials offer a lightweight and compact solution for noise reduction in harsh environment with high heating and ventilation such as launcher fairings [2] or mufflers [3].Metamaterials are plates with periodic structures consisting of small repeated unit cells, with dimensions in the order of centimeters [4].
In recent years, acoustic metamaterials have gained a lot of interest [5][6][7][8][9].One particular property of interest is the use of resonators to create frequency stop-bands to achieve low-frequency sound reduction.A frequency stop-band is a range of frequencies where the free propagation of incoming acoustic waves is prohibited, because of a fano-type-like interference occurs between the incoming and re-radiated waves [10].In Claeys et al. (2013) [7], the potential of applying stop-bands to decrease the vibrational response of panels is discussed.Another property of interest is the inclusion of micro-slits to improve acoustic absorption.When the width of the slits is in the sub-millimeter range (micro-slits), viscous and thermal losses occur in the perforations and improve acoustic absorption.Because of the dimensions of the slits, the accuracy of the manufacturing process is of key importance.Firstly, the sharpness of the edges in contact with the slits can not be guaranteed.Rounded edges heavily affect the impedance of the micro-slits and the absorption.Secondly, the periodicity of the perforation pattern can be not perfect.Both these effects are discussed in Aulitto et al. (2021) [11].Several manufacturing techniques can be employed to realise accurate slits such as laser cutting or milling.Plates with micro-perforations (MPPs) and micro-slit plates (MSPs) backed by a shallow cavity have been introduced by Maa (2000) [12].As shown in recent works, MPPs and MSPs are efficient sound absorbers in the low-frequency range [13][14][15][16].Metamaterials can achieve similar effects by embedding microslits in the design of the unit cell.In each unit cell, the microslits and resonator are created by cutting out a resonant shape instead of positioning a mass-spring resonator on top [4,6].Another advantage of micro-slit metamaterials is that they are easier to manufacture than the original structure with resonators added on top.In the works of Ruiz et al. (2016) [6] and Zielin ´ski et al. (2019) [4], a numerical and experimental study is performed on the normal absorption of the unit cell shown in Fig. 1.The design can be optimized to improve the size of the stop-bands.The stop-bands can be improved by changing the shape of the resonator and by increasing the ratio between the resonant area and the total area of the unit cell.In these works, the presence of the stop-bands is assumed, based on the presence in previous findings on the original metamaterials.
In this paper, an optimized unit cell design based on the work of Ruiz et al. ( 2016) [6] and Zielin ´ski et al. ( 2019) [4] of a micro-slit resonant metamaterial is proposed with larger frequency stopbands and enhanced sound absorption at normal incidence.Furthermore, an elastic numerical model is described to derive the absorption curves for micro-slit resonant metamaterials.The software used for the numerical simulations is COMSOL Multiphysics V5.5 [17].In Section 2, the methodology and proposed data processing technique are discussed to derive the dispersion curves of the unit cells.In this novel algorithm, the stiffness and mass matrices are not used, and the procedure can be applied even when the bending waves are non-smooth, unlike the branch-tracking algorithm discussed in Magliacano et al. (2020) [18].In Section 3, the unit cell design is optimized with the use of genetic algorithms to maximize the size of the first frequency stop-band.In Section 4, the absorption curves of the proposed unit cell design are compared to the design currently used in literature using a combination of rigid and elastic numerical models and the semi-phenomenological JCAPL model [1,[19][20][21][22][23][24].

Methodology
The structure of resonant metamaterials is given by the repetition of the same unit cell.In Fig. 1, an example of such a unit cell is shown: a double-legged resonator (DLR) design with corresponding design variables and slit size as used in the work of Ruiz et al. (2016) [6] and Zielin ´ski et al. (2019) [4].
Floquet-Bloch theory is applied to reduce the computational cost for analysis of these materials [25].The theorem states that the response of a two-dimensional periodic system can be expressed in terms of the response of a reference unit cell and an exponential term describing the amplitude and phase change as the wave travels from one cell to the adjacent cell [26].As stated in Fok et al. (2008) [5], the unit cells are very small and have minimal crosstalk, leaving the individual resonator eigenfrequencies insensitive to lattice parameters and direction.As a result, to describe the behavior of the entire structure, only a single unit cell has to be analyzed.In this work, the dispersive behaviors of various unit cell designs are analyzed along the irreducible Brillioun contour (IBC) 0; 1; 2; 3 # ð0; 0Þ; ð0; LÞ; ðL; LÞ; ð0; 0Þ, where L is the length of the square unit cell [27][28][29].The IBC is the smallest contour in the wave space that captures all information, that is, the minimum and maximum eigenvalues, required to compute the frequency stop-bands for the unit cell.For a 2D periodic square unit cell, the IBC is shown in Fig. 2. The Floquet wavenumber k F = ½k x k y > is spanned along this contour by imposing wavenumbers in xdirection k x , and y-direction k y .Floquet boundary conditions are applied at the edges of the unit cell using the Floquet wavenumber.The Floquet wavenumber is represented by 60 discretizations and the eigenvalue problem is solved with the use of COMSOL Multiphysics [17].The output of the model is a matrix containing the eigenvalues along the IBC used to produce the dispersion curves.The Finite Element Method (FEM) model used to derive the dispersion curves is discussed in Appendix A and a validation of the model is given in Appendix B.

Decoupling waves algorithm
The raw data retrieved from the FEM model contains in-plane waves (longitudinal and transverse waves) and bending waves.The curves are intertwined since the sorting order of the eigenvalues is mixed in the matrix representation.In Fig. 3, this issue is visualized.In Magliacano et al. (2020) [18], a branch-tracking algorithm is discussed for periodic porous materials.The algorithm only considers the gradient between points and does not consider the euclidean distance between them, hence failing when the curves become non-smooth.For the unit cell discussed in the present work, the resonant element introduces additional dispersion curves above its resonance frequencies, that are non-smooth.This phenomenon is similar to the result shown for a two-dimensional infinite structure with a mass-spring system, as discussed in Claeys et al. (2013) [7].Consequently, the algorithm cannot be applied.As a solution, a new algorithm is designed.The proposed algorithm does not utilize the stiffness and mass matrices, allowing for fast computations.The proposed algorithm considers the gradient between points on the dispersion curves and the euclidean distance between points.Furthermore, the new algorithm removes in-plane waves from the dispersion curves as they are inefficient as acoustic radiators compared to bending waves [9].The ratio Ra ¼ t p =L is defined as the ratio between the plate thickness t p and the side length of a unit cell L, as shown in Fig. 1.In this work, only square unit cells are considered.The dispersion curves corresponding to the bending modes are not influenced by variations in the ratio Ra, whereas the transverse and longitudinal waves are.
The ratio Ra is chosen small (i.e.Ra < 0:02), such that the inplane waves have a significantly higher gradient than the bending waves and waves are decoupled along the contour 1; 2 # ð0; LÞ; ðL; LÞ.A threshold check is implemented on the gradient between points to remove the in-plane waves from the data.Estimates of the bending waves are created with the use of 4th order polynomial fits.The eigenvalue branches are tracked based on the difference between the estimates and available points.Furthermore, the fits are iteratively updated to improve accuracy.A description of the algorithm is provided in C. The algorithm is able to remove the in-plane waves from the dispersion curves and sort the remaining eigenvalues of the raw data to obtain the dispersion curves even when the bending waves are non-smooth.

Optimization of the unit cell
For the simulations, the same structural properties as used in the work of Ruiz et al. (2016) [6] are considered, namely, density q ¼ 950 kg/m 3 , Young's modulus E ¼ 1750 MPa, and Poisson's ratio m ¼ 0:3.Furthermore, the unscaled unit cells have the same slit size s, plate thickness t p , and plate length L, as the DLR design shown in Fig. 1.A ratio Ra ¼ 0:02 is used for scaling the unit cells.Scaling is performed by multiplying the length of the unit cell L, the slit size s, and the design parameters with a scaling factor.The plate thickness t p is unaffected by scaling.

Genetic algorithms methodology
A genetic algorithm is a search heuristic inspired by the principle of natural selection [30].In this work, genetic algorithms are used to optimize the design of a unit cell design to maximize the size of the first frequency stop-band.The metric for optimization is the stop-band factor (SBF), which is defined as the ratio between the lower bound of a stop-band and its upper bound.Furthermore, two SBFs are considered: the first SBF (SBF1), and the second SBF (SBF2).The index refers to the appearance of the stop-band ranked from low to higher frequencies.In Fig. 4, the dispersion curves and frequency stop-bands are shown for the DLR design as shown in Fig. 1.The frequency range of the first stop-band is around 5 kHz.Note that this range is obtained by using the same design, design parameters, and structural properties as used in the work of Ruiz et al. (2016) [6] and Zielin ´ski et al. ( 2019) [4].A lower range can be realized by choosing a different material with, for instance, a lower Young's modulus.In Fig. 4, it can be seen that the first stopband is located between the first and second curve, as is also found in the work of Claeys et al. (2016) [9], and the second stop-band between the third and fourth curve.SBF1 is then defined as the ratio between the minimum of the second curve and the maximum of the first curve.Likewise, SBF2 is then defined as the ratio between the minimum of the fourth curve and the maximum of the third curve.SBF1 is chosen as the metric for optimization because it is located in the frequency range of application of micro-slit resonant metamaterials.

Influence unit cell design characteristics
In this subsection, the influence of the unit cell design characteristics on the stop-bands is discussed.The design considered, is the DLR design shown in Fig. 1.By linearly increasing the variables d 1 ; d 2 , and d 3 , the relative size of the resonator with respect to the surface of the unit cell increases.An increase in the relative size of the resonator increases both SBF1 and SBF2.By increasing the relative size of the resonant structure, the maximum kinetic energy of the structure increases as well.This leads to a greater fano-typelike interference, and to larger stop-bands.By choosing d 1 and d 2 constant, and varying d 3 , the influence of the resonator mass is investigated.An increase in resonator mass increases SBF1, albeit smaller than the increase observed in SBF1 by increasing the resonator size.By choosing d 1 constant, d 2 þ d 3 constant, and varying d 2 , the influence of the resonator stiffness is investigated.An increase in resonator stiffness increases SBF1 however decreases SBF2.Lastly, an increase in slit size has a small negative effect on SBF1 and SBF2.Note that due to manufacturability and accuracy constraints, a slit size of 0:3 mm is considered.

Designs considered during optimization
The DLR design as shown in Fig. 1, is optimized to increase SBF1.The accompanying dispersion curves are shown in Fig. 4, in which it can be seen that SBF1 ¼ 1:38 and SBF2 ¼ 1:05.For manufacturability, constraints are applied to the design variables to ensure a minimum length of 1 mm for each variable and a minimum distance of 1 mm between the slit and the edges of the plate.The optimized geometry and dispersion curves are shown in Fig. 5.In Fig. 5b, the dispersion curves of the optimized DLR unit cell design are plotted over the IBC.It can be seen that SBF1 ¼ 1:61 and SBF2 ¼ 1:18, an increase of 16.7% and 12.3%, respectively.It can be seen  that the resulting resonant structure is maximized within the given constraints.Other notable designs considered during optimization are shown in Fig. 6.A single-legged resonator design is considered to compare its performance to the DLR design.Double resonator designs are explored to see if the interaction between two resonators at different frequencies can lead to greater stop-bands.Lastly, the implementation of internal slits in the resonant structure is investigated to see how small changes in the mode shapes can alter the stop-band behavior.Again, for manufacturability, constraints are applied to the design variables to ensure a minimum length of 1 mm for each variable and a minimum distance of 1 mm between the slit and the edges of the plate.The resulting SBFs are shown in Table 1.
In Fig. 6 and Table 1, it can be seen that the geometries converge to a configuration where the size of the resonator is maximized within the design constraints.The single-leg design, as shown in Fig. 6a, significantly underperforms the optimized DLR design, as shown in Fig. 5a.The addition of the second leg increases the maximum displacement of the resonant shape, and thereby the maximum kinetic energy of the resonant structure.In Fig. 6, it can be seen that the addition of internal slits, Figs.6d and 6e, significantly improves SBF2.The stop-bands of these designs are located around the same frequencies as the DLR designs considered earlier, see Figs. 4 and 6b.The best performing design is the one with the internal slit in the top right corner of the cell, as shown in Fig. 6e.The motivation for this design is further elucidated in Section 3.4.To improve the manufacturability of the design, a slanted trim (ST) design is proposed.In Fig. 7, the unit cell design is depicted and the corresponding dispersion curves are shown.In Fig. 7, it can be seen that SBF1 ¼ 1:65 and SBF2 ¼ 1:31.All the designs  considered during optimization converge to a configurations where the size of the resonator is maximized within each unit cell.

ST design results
The implementation of the ST design is motivated by comparing the mode shapes of the optimized DLR unit cell design (Fig. 6a) to the ST design (Fig. 7a).In Figs. 8 and 9 the real displacements of the unit cell at the bounds of the first stop-band are shown for the optimized DLR and for the trimmed DLR, respectively.In Figs. 8  and 9, it can be seen that there is a difference between the lower and upper bound mode shapes of the first stop-band.The deflection for the lower bound gradually increases along the diagonal of the cell.For the upper bound, this increase is much steeper.For the design in Fig. 8 and the material properties as shown in Section 3, the resonance frequencies are 4.4 kHz and 7.1 kHz for the lower and upper bounds, respectively.For the design in Fig. 9 and the material properties as shown in Section 3, the resonance frequencies are 4.4 kHz and 7.3 kHz for the lower and upper bounds, respectively.By reducing the mass of the resonant cell at the maximum deflection (the top right corner), the mode shape of the lower bound is not significantly affected.However, the mode shape of the upper bound becomes stiffer and therefore an increase in resonance frequency and SBF1 is realized.The ST design has an increase in SBF1 of 20% and SBF2 of 25% with respect to the DLR design currently used in literature, see Figs. 1 and 4. A small increase (< 1%) in SBF1 can be realized by the implementation of an additional internal and external slit.However, the small increase does not justify the increase in manufacturing complexity.

Absorption
In this section, the absorption curves are compared for the DLR and the ST design.Furthermore, the absorption curves are also compared to the numerical and experimental results obtained in the work of Zielin ´ski et al. (2019) [4].For the simulations, the same properties for the plate are considered as described in Section 3.

Methodology
The two-microphones method, as discussed by Bodén et al. (1986) [31], Jang et al. (1998) [32], and Labašová et al. (2019) [33], is implemented in a numerical model to compute the absorption coefficient at normal incidence for the micro-slitted metamaterial.In Fig. 10, an overview of the numerical model is shown.The linearized Navier-Stokes module is used to consider the viscous and thermal effects caused by the slits [34].The model is composed of four layers, as shown in Fig. 10; a perfectly matched layer (PML), a background pressure field (BPF), the plate with the slits, and a back cavity layer.The PML acts as a perfect absorber, which ensures that no waves are reflected into the BPF.The BPF is used to impose an incident pressure wave with a certain frequency.The structural mechanics module is used to model plate with the slits as an elastic body [35].The height of the BPF is empirically chosen at 30 mm.This height allows the mesh to transition from small elements at the slits of the plate to larger elements at the microphones.Furthermore, taller heights increase the number of total elements in the model but do not significantly improve the accuracy.Similarly, the height of the PML is chosen at 10 mm.Tetrahedral elements are used in the BPF, the plate, and the back cavity layer.A swept mesh is used for the PML and the slits in the plate.To model an infinite plate, symmetric boundary conditions are applied in the x and z directions as displayed in 10.Furthermore, a rigid wall boundary condition is applied at the rightmost plane in Fig. 10.At the microphones, the average gross pressure is computed.The absorption coefficient is then given by where H 12 is the ratio between the average gross pressure at microphone 2 and the average gross pressure at microphone 1 [31] (2010) [24].The JCAPL model is used to take into account the viscous and thermal effects of the slits.

Results
The JCAPL model is used to derive an analytical solution for the absorption curves.The parameters used for the JCAPL model are shown in Table 2.
Absorption curves are computed with the rigid numerical model for two cavity depths: 30 mm and 53 mm.In the rigid numerical model, the plate itself is not modeled; only the slits are.Therefore, the structural effects of the plate, such as the resonance of the resonator, are not considered.Given the configurations, the first stop-band is located between 4.2 and 5.8 kHz for the DLR design, and between 4.4 and 7.3 kHz for the ST design.Since both stop-bands are outside the frequency range of interest, a rigid numerical model is used.The absorption curves for the DLR   In Figs.11 and 12, it can be seen that for both cavity depths, the rigid numerical model closely resembles the analytical solution.Furthermore, it can be seen that the results from the rigid model and the analytical solution for the DLR design are in agreement with the numerical results and impedance tube measurements done in the work of Zielin ´ski et al. ( 2019) [4].It can be seen that the ST design produces a higher sound absorption peak at a lower frequency with respect to the DLR design.The ST design shows a 9% increase in the first peak of absorption coefficient compared to the DLR design at a cavity depth of 30 mm, and an increase of 10% at a cavity of size 53 mm.The lower porosity of the ST design decreases the resonance frequency of the Helmholtz resonator, see Table 2.The smaller characteristic lengths of the ST design lead to a higher sound absorption peak due to the smaller losses in the slits.
To inspect the influence of a stop-band on the resulting absorption curves of the ST design, an elastic numerical model is used.To artificially reduce the frequency bounds of the first stop-band, the Young's modulus is reduced from E ¼ 1750 MPa to E red ¼ 26:25 MPa.This brings the first stop-band to lie between 545 and 900 Hz.In Fig. 13, the absorption curves are shown for the ST design for a cavity depth of 30 mm using the elastic numerical model.
In Fig. 13, it can be seen that there is a negligible difference between the absorption curves for the two different Young's moduli.Also, the difference between the rigid and elastic numerical models is negligible.For the elastic model, it appears that the stop-band behavior and structural properties of the plate have a negligible influence on the resulting absorption curves.From the perspective of the JCAPL model, this is an expected result since there is no dependency on the structural parameters of the plate in this model.For both numerical models, the wall at the back of the cavity is considered to be perfectly rigid, meaning that the transmitted particle velocity is zero.The only losses occur due to the viscous and thermal effects of the slits.When there exists a transmitted pressure wave through the back cavity, for instance due to an absence of a perfectly rigid back cavity wall, the structural properties of the plate will affect the reflected pressure wave and therefore also the absorption coefficient.The elastic numerical model would have to be extended with an elastic back cavity wall to capture these effects.Similarly, the JCAPL model would have to be extended to include a dependency on the effective density and bulk modulus of the plate to capture its structural effects.

Conclusion
An optimized unit cell design of a micro-slit resonant metamaterial is proposed to increase the size of the frequency stop-bands and to enhance the sound absorption at normal incidence.The design is referred to as the ST design.A FEM model is used to derive the dispersion curves of various unit cell designs.To post-process the output of the FEM model, an algorithm is proposed.The algorithm removes the in-plane waves from the dispersion curves and sorts the remaining eigenvalues of the raw data to obtain the dispersion curves even when the bending waves are nonsmooth.The proposed algorithm does not utilize the stiffness and mass matrices, allowing for fast computations.Unit cell designs are optimized to maximize the size of the frequency stop-bands.The first stop-band factor is chosen as the metric for optimization since it is located in the frequency range of application of micro-slit resonant metamaterials.Optimized designs converge to a configuration where the relative size of the resonant  structure with respect to the surface of the unit cell is maximal.By increasing the relative size of the resonant structure, the maximum kinetic energy of the structure increases as well.In turn, this leads to a greater fano-type-like interference, and to larger stop-bands.The ST design reduces the mass of the resonant cell at the maximum deflection.The resonance frequencies of both the lower and upper bound of the mode shapes increase; however, the ratio between the resonance frequencies increases as well.The ST design has an increase in SBF1 of 20% and SBF2 of 25% with respect to the DLR design currently used in literature.A small increase (< 1%) in SBF1 can be realized by the implementation of an additional internal and external slit.However, the small increase does not justify the increase in manufacturing complexity.The ST design shows a 9% increase in the first peak of absorption coefficient compared to the DLR design at a cavity depth of 30 mm, and an increase of 10% at a cavity of size 53 mm.The smaller characteristic lengths of the ST design lead to a higher sound absorption peak due to the smaller losses in the slits.Furthermore, the lower porosity of the ST design reduces the frequency of the first peak.Stop-band behavior does not influence sound absorption at normal incidence of acoustic waves in the frequency range of interest.To capture the structural effects of the unit cell, the elastic numerical model would have to be extended with an elastic back cavity wall.Similarly, the JCAPL model would have to be extended to include a dependency on the effective density and bulk modulus of the plate to capture its structural effects.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.the wave vector k is composed of two components: k = ½k x k y T .By use of superposition, the longitudinal, transverse, and first bending wave, respectively, can be written as functions of k x and k y and are given by where E is the Young's modulus, q the density, m Poisson's ratio, and t p the thickness of the plate.The first frequency stop-band is designed around the maximum value of the first bending wave mode (f B 1 ), therefore bending wave modes up to order 3 are considered and are given by ðB:5Þ where L is the length of the square unit cell.For the dispersion diagram, the frequency is made dimensionless by dividing the wave frequency by which corresponds to the frequency of the first bending wave mode when the applied wavelength is equal to twice the length of the infinite plate k ¼ 2L.The dispersion curves for the analytically and numerically derived wave types for a ratio Ra ¼ 0:002 along the IBC 0; 1; 2; 3 # ð0; 0Þ; ð0; LÞ; ðL; LÞ; ð0; 0Þ, are shown in Fig. 15.In Fig. 15, it can be seen that the numerical solutions converge to the analytical solution.Furthermore, it can also be seen that the numerical solutions jump between wave types.To clarify, near the point ð0; 0Þ on the IBC, the numerical solution jumps from the bending wave to the transverse wave.This effect occurs due to the way the points that make up the curves are connected and is strictly numerical.As a remedy, the decoupling wave types algorithm, as discussed in C, is used.It can be concluded that the used numerical modeling approach to derive dispersion diagrams for square unit cell type plates is valid for the ratio Ra ¼ 0:002.

Appendix C. Decoupling wave types algorithm
The Matlab implementation of this algorithm is released and can be found on [36].The algorithm derives a matrix F that contains the first n c 2 Z >0 decoupled bending waves, in which Z >0 denotes the integer numbers greater than zero.The input data w of the algorithm contains the computed eigenvalues along the IBC, as discussed in A, and the corresponding wavenumbers j spanning the IBC.The number of points in j is given by n k 2 Z >0 .
As a constraint, the following relation must hold: n f P n c þ 2, in which n f 2 Z >0 is the number of eigenvalues computed at each point along the IBC.This constraint is set in place to have sufficient data present in the input data: there exist three eigenvalues at j ¼ 0 (the point (0,0) of the IBC), where the real part of the eigenvalues is zero.These solutions correspond to a single bending and both in-plane waves.Therefore, to include enough information in the input data, the number of eigenvalues computed must be at least 2 higher than the desired number of decoupled bending waves.The assumptions for the algorithm are: Ra is sufficiently small, such that the low-frequency presence of the longitudinal and transverse waves is only visible for j 2 ½0; 1 2 [ ½2 1 2 ; 3. In the input data, the first n c bending waves are decoupled in the range j 2 ½1; 2.
The initial guess for F is the input data w.Note that only the first n c curves of w are taken for the initial guess of F: dim(F)dim(w).The first part of the algorithm (lines 2 to 10) deals with the inplane waves present in the data.For each point in 2 ; 3, the bending curves intertwined with in-plane waves can be found by looking at the gradient.The in-plane waves have a significantly greater gradient than the bending waves when Ra is sufficiently small.Therefore, when the gradient exceeds a certain threshold b 2 R >0 , in which R >0 denotes the real numbers greater than zero, it is most likely that the curve at that j value is an in-plane wave or is intertwined with it and the corresponding entry in F is set to Not a Number (NaN) (lines 2 to 5).The branches that contain NaN gaps, as a result of this threshold check, are then pruned in two directions (lines 6 to 10): backwards for j 2 ½0; 1 2 (i.e.j ¼ 1 2 !0) and forwards for j 2 ½2 1 2 ; 3 (i.e.j ¼ 2 1 2 !3).The matrix F now contains an initial estimate of the decoupled waves as well as empty entries.The second part (lines 11 to 20) allocates points from w into the right position in F for j 2 ½0; 1 and fills in the empty entries.A fit is made for each iteration over j that is extrapolated to create an estimated value for F j;j .This estimated value is compared to all points in w for the same j and is stored in a matrix M. By looking at the minimum of M, the sorting of points in F for w is solved (line 18).Furthermore, to prevent duplication of points in F, used points are immediately replaced by NaN values in M (lines 19 and 20).The third part (lines 21 to 28) allocates points from w into the right position in F for j 2 ½2; 3.An estimate C is made for F with j 2 ½2; 3 by means of superposition of F with j 2 ½1; 2 and F with j 2 ½0; 1 : C ¼ F 1!0 þ F 2!1 .Again, a matrix M is made and the same concept applies as discussed in part two of the algorithm.The output of the algorithm is matrix F which contains the first n c n f 2 Z >0 decoupled bending waves.

13:
Replace all entries in F q min!3 ;j with NaN 14: end for 15: for q ¼ 1 !0 do 16: for j ¼ 1 !n c do 17: Create a fit of order O f on F j;j where j 2 ½0; 1 18: Compute frequency estimates based on the fit, and store results in vector J j nc Â1 19: end for 20: Compute a matrix M n f Ânc that contains the absolute difference between J 1!nc and w q;1!n f 21: for j ¼ 1 !n c do 22: Find the x and y indexes of the minimum of M: minðMÞ ¼ M x;y 23: F q;x ¼ w q;y 24: M y;1!nc ¼ NaN 25: M 1!n f ;x ¼NaN 26: end for 27: end for 28: Create an estimate C n k Ânc for F 29: for q ¼ 2 ! 3 do 30: Compute a matrix M n f Ânc that contains the absolute difference between the C q;1!nc and w q;1!n f 31: for j ¼ 1 !n c do 32: Find the x and y indexes of the minimum of M: minðMÞ ¼ M x;y 33: F q;x ¼ w q;y 34: M y;

Fig. 3 .
Fig. 3. Dispersion curves of the DLR unit cell design.In gray, the dispersion curves plotted from the raw model output data.In black, the dispersion curves after the decoupling algorithm.

Fig. 4 .
Fig. 4. Dispersion curves and frequency stop-bands of the DLR unit cell design.

Fig. 5 .
Fig. 5. Geometry and dispersion curves of the optimized DLR unit cell design with respect to SBF1.

Fig. 8 .
Fig. 8. Modes shapes displaying the out-of-plane real displacement of the first stop-band of the optimized DLR unit cell design.

Fig. 9 .
Fig. 9. Modes shapes displaying the out-of-plane real displacement of the first stop-band of the ST unit cell design.

Fig. 10 .
Fig. 10.Overview of the numerical model.The curly brackets denote the layers of the model, namely, the perfectly matched layer, the background pressure field, the plate with slits, and the back cavity.The planes corresponding to microphones 1 and 2 are denoted by m1 and m2, respectively.The distance between the two microphones is denoted by dmic.Frequency of imposed wave f ¼ 800 Hz.Mesh: 208219 domain elements, 23011 boundary elements, and 2009 edge elements.

Fig. 11 .
Fig. 11.Absorption curves for the DLR design and the ST design.A cavity depth of 30 mm is used.The JCAPL solutions are represented by the solid lines and the rigid numerical solutions by the circles.The results obtained for the DLR design in the work of Zielin ´ski et al. (2019) [4] are displayed by the dashed lines.

Fig. 12 .
Fig. 12. Absorption curves for the DLR design and the ST design.A cavity depth of 53 mm is used.The JCAPL solutions are represented by the solid lines and the rigid numerical solutions by the circles.The results obtained for the DLR design in the work of Zielin ´ski et al. (2019) [4] are displayed by the dashed lines.

Fig. 13 .
Fig.13.Absorption curves the ST design.A cavity depth of 30 mm is used.The JCAPL solution is represented by the solid line, the rigid numerical solutions by the circles, and the elastic numerical solutions by the dotted lines.

Fig. 14 .
Fig. 14.Definition of the source and destination planes for the unit cell.The blue line depicts the source plane and the green line the destination plane.

Table 1
SBF1 and SBF2 of the optimized geometries of notable designs considered during optimization.
Fig. 7. Geometry and dispersion curves of the ST unit cell design.

Table 2
Derived JCAPL parameters for the DLR design and ST design.