Flexural wave bandgap properties of phononic crystal beams with interval parameters

Uncertainties are unavoidable in practical engineering, and phononic crystals are no exception. In this paper, the uncertainties are treated as the interval parameters, and an interval phononic crystal beam model is established. A perturbation-based interval finite element method (P-IFEM) and an affine-based interval finite element method (A-IFEM) are proposed to study the dynamic response of this interval phononic crystal beam, based on which an interval vibration transmission analysis can be easily implemented and the safe bandgap can be defined. Finally, two numerical examples are investigated to demonstrate the effectiveness and accuracy of the P-IFEM and A-IFEM. Results show that the safe bandgap range may even decrease by 10% compared with the deterministic bandgap without considering the uncertainties.


Introduction
Phononic crystal, which is a type of artificial materials that exhibit a periodic arrangement of unit cells, has been of interest for years due to its bandgap properties. Bandgap properties mean that elastic wave is predominantly evanescent (non-propagating) within the bandgap frequency range [1] . Therefore, phononic crystal structures are of application potentials in vibration isolation and sound insulation. For example, Zhao et al. [2] proposed a phononic crystal vibration isolator to attenuate railway-induced vibration. Wang et al. [3] proposed a phononic crystal frame structure for vibration isolation and noise reduction in high-rise buildings. Presently, researchers have developed two formation mechanisms of the bandgap, namely Bragg scattering [4] and local resonance [5] . Compared with the Bragg scattering mechanism, the local resonance mechanism is usually more favorable to the formation of low-frequency bandgap [6] . In addition, common methods for calculating bandgap range include the plane wave expansion method [7] , the finite element method [8] , and the transfer matrix method [9] .
The above numerical methods are all based on deterministic parameters. However, the true values of design parameters are always subject to some extent of uncertainties in the real world. Uncertainties triggered by initial manufacturing errors, harsh environmental factors, and fluctuations in material properties are unavoidable in practical engineering, and phononic crystals are no exception. For example, Li et al. [10] pointed out that, due to manufacturing, measurement and assembling errors, uncertainties mainly exist in the cable forces and node locations in the cable net structure modeling process. In the multibody mechanical systems, Wu et al. [11] pointed out that a variety of uncertainties exist in materials, loads, geometric dimensions, manufacturing tolerances, etc. Jiang et al. [12] showed that, when a vehicle rides over a rough road, the complex road conditions would lead to a high degree of uncertainties in the loads on the tires. In the interval model, only the lower and upper bounds of uncertain parameters are needed, while the specific distribution form is not needed. Therefore, the interval model is used to model these uncertainties in this paper. In recent decades, the interval model has been applied in the field of finite element analysis giving rise to the so-called interval finite element method (IFEM). Most research of IFEM focuses on how to find the solution set of interval linear equations as accurately as possible. To the authors' knowledge, among the numerous IFEMs, two major methods have been proved to be successful, namely the optimization method and the modified interval arithmetic method.
The optimization method is to obtain the global minimum and maximum values of the response by executing the program of searching the extreme value in the interval domain of the uncertain parameters. The sampling-based algorithm, such as the Monte Carlo method (MCM) [13] and the vertex method [14] , is the most popular one due to its nonintrusive property. Other commonly used optimization methods include the genetic algorithm and the particle swarm optimization algorithm. However, with the increase in the number of interval parameters, the calculation amount will increase dramatically. In recent years, the response surface method has been proposed to reduce the complexity of calculation by establishing the approximate response function and optimizing the approximate response function on the response surface [15] .
On the other hand, due to its philosophy of simplicity, interval arithmetic is another way of dealing with interval linear equations. However, such an interval arithmetic strategy often leads to meaningless and overly conservative results, mainly because the dependence between different variables is ignored in the process of interval arithmetic. This is the so-called dependency problem. Therefore, a variety of modified interval arithmetic methods have been developed for dealing with the overestimation phenomenon. By expanding the interval dynamic stiffness matrix and the interval load vector into the first-order Taylor series, and using the first-order Neumann series to approximate the inverse of the interval dynamic stiffness matrix, Xia and Yu [16] proposed a perturbation-based interval finite element method (P-IFEM). Because the series expansion has truncation errors, the above P-IFEM is only applicable for the cases that the uncertain parameters are of a low degree of uncertainties. Therefore, many modified methods are devoted to preserving the high-order terms of the series in the calculation [17] . In addition to Taylor series, affine arithmetic [18] is also a method to automatically keep track of uncertain parameters. On the other hand, Xiang and Shi [19] proposed an element-by-element technique to decompose the uncertain parameters from the system matrix. Through this method, the number of occurrences of uncertain parameters in the calculation process is reduced, thus solving the dependency problem in the interval arithmetic.
At present, a large number of researchers have used optimization algorithms to solve the extreme values of the beginning frequency or ending frequency of the bandgap, when input parameters are limited within the corresponding interval domain. Taking the relative bandgap as the optimization target, Han and Zhang [20] optimized the phononic crystal thin plate structures by using the genetic algorithm. Further, in order to increase the width of the bandgap and reduce the beginning frequency of the bandgap, Chen et al. [21] optimized the hexachiral phononic crystal by integrating the genetic algorithm into the secondary development platform of ISIGHT. By using the response surface method, Zhai et al. [22][23][24] and Liu et al. [25] established the approximate relationship between the beginning frequency or ending frequency of the bandgap and some geometric parameters of two-dimensional locally resonant phononic crystals, and used this relationship to optimize the structures. Li et al. [26] established the inherent relationship between bandgaps and topological features through a deep learning based data-driven method. The trained model can be used to design phononic crystal structures with anticipated bandgaps.
On the other hand, the dynamic response of phononic crystals with interval parameters can also be studied by using modified interval arithmetic methods. Wu et al. [27] obtained the interval band structure and interval frequency response function of acoustic metamaterials based on the P-IFEM. However, the influence of interval parameters on the bandgap was not specifically studied. In contrast, He et al. [28] used the stochastic finite element technique to deal with the material property uncertainty in acoustic metamaterials, and studied the influence of the uncertainty on the beginning frequency of bandgap. To the authors' knowledge, there are few other studies on the dynamic response of phononic crystals with interval parameters. Therefore, further research on the numerical method for phononic crystals with interval parameters is in demand and quite promising. The remainder of this paper is organized as follows. In Section 2, the model of interval phononic crystal beam is established. Then, both the P-IFEM and the affine-based interval finite element method (A-IFEM) are used to study the dynamic response of this interval phononic crystal beam. In Section 3, a phononic crystal beam with interval Young's modulus or thickness is studied to verify the effectiveness and accuracy of P-IFEM and A-IFEM. The conclusions are summarized in Section 4.

Model and method
In this section, the bandgap calculation of a deterministic phononic crystal beam is introduced at first. Then, consider that uncertainties may also exist in this phononic crystal beam. These uncertain parameters are treated as interval parameters, and an interval phononic crystal beam model is proposed. Both the P-IFEM and the A-IFEM are developed to evaluate the bounds of the dynamic response interval vector. Finally, in order to intuitively reflect the impact of uncertainty on the bandgap, an interval vibration transmission analysis is carried out.

Bandgap calculation method of phononic crystal beam
As shown in Fig. 1, by periodically repeating aluminum (see the orange part in Fig. 1) and polyethylene (see the blue part in Fig. 1) in the axial direction, a phononic crystal beam with N u unit cells is formed. The lengths of each aluminum and polyethylene are L 1 and L 2 , respectively, which makes the lattice constant be a = L 1 +L 2 . The cross section of this phononic crystal beam  is rectangular, and the corresponding width and thickness are b and h, respectively. Assume that each beam is divided into N e Euler-Bernoulli beam elements. Then, the whole phononic crystal beam consists of 2N e N u beam elements (each unit cell contains an aluminum beam and a polyethylene beam). The stiffness matrix and mass matrix of the Euler-Bernoulli beam element can be calculated by Here, I (= bh 3 /12) and S (= bh) denote the moment of inertia and area of the cross section, respectively. E and ρ represent Young's modulus and density of the beam element, respectively. L is the length of the beam element. Then, the dynamic equation of the whole phononic crystal beam can be obtained through the standard assembly procedure, i.e., where Here, ω is the circular frequency. L (i) j is the connectivity matrix of the j th element in the i th beam. U and F are the vectors of dynamic responses and forces, respectively. D, K, and M represent the global dynamic stiffness matrix, the global stiffness matrix, and the global mass matrix, respectively. K (i) and M (i) denote the global stiffness matrix and the global mass matrix of the i th beam, respectively.
The vibration transmission analysis of the corresponding finite phononic crystal beam is an effective approach for studying bandgap characteristics, which is defined as follows: where T represents the transmission coefficient of vibration, and u in and u out denote the amplitudes of the displacements of the excitation point and the response point, respectively. As shown in Fig. 1, this phononic crystal beam is only subject to a single frequency harmonic displacement at the left end of the beam while the other end is kept free. The right end of the beam is selected as the response point. In order to facilitate the subsequent study, the numbering of beam elements is carried out in this paper from left to right along the whole phononic crystal beam. Thus, the first row of vector U represents the displacement of the excitation point, and the penultimate row of vector U represents the displacement of the response point. In order to facilitate the vibration transmission analysis, the displacement of the excitation point is set as 1 by the method of setting a big number (U 11 = u in = 1). The specific implementation method is to first modify the (1,1)th entry of the dynamic stiffness matrix and make it equal to a value that is much greater than other elements in the global dynamic stiffness matrix, such as D 11 = 10 40 . Then, define the vector F as follows: When the amplitude of the response point is less than that of the excitation point, that is to say, T < 0, vibration attenuation occurs during the propagation. Accordingly, the frequency range corresponding to T < 0 is called the bandgap in academia.

Phononic crystal beam with interval parameters
Obviously, uncertainties also exist in phononic crystal models. If only the upper bound and lower bound of uncertain parameters are known, the interval model is one of the most effective methods to deal with uncertainties. Thus, the dynamic equation of the deterministic phononic crystal beam as shown in Eq. (2) can be extended to the following interval dynamic equation: where D(b I ) represents the interval global dynamic stiffness matrix, U I denotes the dynamic response interval vector, F is the structural load vector, and b I represents an n-dimensional closed interval parameter vector The underline and overline represent the lower and upper bounds, respectively. The midpoint and the radius of the interval parameter vector b I are defined as It should be pointed out that the vector F is defined by Eq. (5), so it is not a function of the interval parameter vector b I . In this paper, both P-IFEM and A-IFEM are proposed to solve the interval linear equations as shown in Eq. (6).

P-IFEM
The first-order Taylor series expansion of the dynamic response interval vector U I around the midpoint values of the interval parameters can be expressed as Here, U m and ∆U I represent the midpoint and the radius of the dynamic response interval vector U I , respectively.
can be obtained as Substitute Eq. (10) into Eq. (9). Then, the dynamic response interval vector U I can be expressed in the following form: Obviously, U I is a monotonic function of [ε i ]. Thus, the upper and lower bounds of the dynamic response interval vector U I can be given as

A-IFEM
In fact, a direct solution of Eq. (6) is to use interval arithmetic instead of deterministic arithmetic to solve the interval vector U I . However, such an interval arithmetic strategy often leads to meaningless and overly conservative results, mainly because the dependence between different variables is ignored in the process of interval arithmetic. This is the so-called dependency problem. In order to reduce or eliminate the overestimation phenomenon, affine arithmetic, which can keep track of the dependency between different variables throughout the calculation process, was proposed in the literature [18] .
Affine arithmetic represents uncertain quantities during the analysis by an affine form x, which can be expressed as Here, x 0 is the center value of the interval. x i [ε i ] represents the linear dependency of x on the input interval variables. Specifically, [ε i ], termed noise symbol, represents an independent source of the total uncertainty whose value is unknown but is limited within the range [−1, 1]. Meanwhile, the coefficient x i , called partial deviation, gives the magnitude of that particular uncertainty. Considering the influence of nonlinear dependencies, x e [ε e ] is added at the end of Eq. (14) as an error term. [ε e ] is a new noise symbol, and its range is also an interval [−1, 1]. By definition, x e is a non-negative number. Similarly, if one or more entries of a matrix can be affine variables, the matrix is called an affine matrix A, where m and k represent the numbers of rows and columns of the affine matrix A, respectively. The affine matrix A can be expressed as follows: The numerical examples in this paper mainly consider the uncertainties of Young's modulus and thickness of the beam. First, a phononic crystal beam with interval Young's modulus is introduced. It is assumed that Young's modulus of the 2N u beams shown in Fig. 1 is an independent interval parameter, defined as follows: Here, α I i = [−∆α i , ∆α i ] is a symmetric interval variable, and ∆α i is its radius. The interval variable α I i represents the dimensionless fluctuation around the nominal value E (i) 0 . It is well known that Young's modulus only affects the stiffness matrix. Thus, the interval global stiffness matrix can be expressed as follows: where K 0 represents the global stiffness matrix corresponding to the nominal value of Young's modulus. Similarly, K (i) stands for the global stiffness matrix of the i th beam corresponding to the nominal value of Young's modulus. Further, the interval global dynamic stiffness matrix can be deduced as where Next, a phononic crystal beam with interval thickness is described in detail. Similarly, the thickness of 2N u beams can be described in the affine form, It can be seen from Eq. (1) that the stiffness matrix is proportional to the cubic power of the thickness of the beam. Thus, the interval global stiffness matrix of the i th beam is In affine arithmetic, one common technique for approximating the cubic operation is the Chebyshev approximation, and its result can be given by where Here, [ε ie ] is a new noise symbol. Therefore, the interval global stiffness matrix of the i th beam can be expressed in the affine form, On the other hand, Eq. (1) shows that the mass matrix is proportional to the thickness of the beam. Therefore, it can be deduced that the interval global mass matrix of the i th beam is The interval global dynamic stiffness matrix of the i th beam can be further deduced as Here, D , and K (i) and M (i) represent the global stiffness matrix and the mass matrix of the i th beam corresponding to the nominal value of thickness, respectively. Through the standard assembly procedure, the global interval dynamic stiffness matrix can be obtained as From the analysis of the above two examples, it can be seen that the interval dynamic stiffness matrix can be described in the form of affine matrix shown in Eq. (16). For the convenience of subsequent derivation, the interval global dynamic stiffness matrix is uniformly described as where N ε represents the number of noise symbols. N ε = 2N u when considering the uncertainty of Young's modulus and N ε = 4N u when considering the uncertainty of thickness. It is worth remarking that, because the phononic crystal beam is assumed to be actuated with a single frequency harmonic displacement at the left end of the beam, Eq. (31) needs to be processed as follows: Substituting Eq. (31) into Eq. (6) yields where A i = D −1 0 D i and A e = D −1 0 D e . Degrauwe et al. [18] proposed a novel method to solve the inverse of the affine matrix A, which can be expressed as Through the analysis of the above two examples, the affine form of the interval global dynamic stiffness matrix for the phononic crystal beam with interval Young's modulus or thickness may not contain the error term (D e = A e = 0). Therefore, the dynamic response interval vector can finally be expressed as Then, the upper and lower bounds of the dynamic response interval vector can be obtained according to the affine algorithm rules shown in Ref. [18].

Interval vibration transmission analysis
According to the previous description, the penultimate row of vector U I represents the displacement of the response point, whose the lower bound u out and the upper bound u out can be obtained by the P-IFEM or the A-IFEM. Thus, the interval transmission coefficient at a given frequency point can be defined as

Numerical examples
In this section, an example of bandgap calculation for a deterministic phononic crystal beam is first introduced. Then, two numerical examples are studied to validate the effectiveness and accuracy of P-IFEM and A-IFEM to analyze the phononic crystal beam with interval Young's modulus or thickness. Meanwhile, an MCM with one million samples is regarded as a reference solution for the phononic crystal beam with interval parameters.

Phononic crystal beam with interval Young's modulus
In this subsection, an interval vibration transmission analysis of the phononic crystal beam with interval Young's modulus is carried out. During the analysis, Young's moduli of all beams are set as interval parameters, and all the other parameters are the same as those in Subsection 3.1. First, let us assume that the radius of the interval variable α I i of each beam is ∆α i = ∆α Al = ∆α Po = 0.05 (i = 1, 2, · · · , 2N u ). ∆α Al and ∆α Po represent the radii of the interval variable α I i of each aluminum and polyethylene, respectively. Figure 3(a) shows the results of the interval vibration transmission analysis. The black line in Fig. 3(a) denotes the vibration transmission corresponding to the nominal value of Young's modulus, namely the calculation results in Fig. 2. The other three lines represent the upper bounds of the interval transmission coefficient predicted by the P-IFEM, MCM, and A-IFEM. In order to clearly show the influence of uncertainty on the bandgap, Figs. 4(a) and 4(b) show the interval vibration transmission analysis results around the bandgap. As can be seen from Fig. 4(a), the  Table 1. It can be seen from Fig. 4(a) and Table 1 that the first safe bandgap predicted by the P-IFEM is slightly wider than that calculated by the MCM, while the first safe bandgap predicted by the A-IFEM is significantly narrower than that calculated by the MCM. Figure 4(b) shows the second safe bandgap, and the bandgap ranges are also listed in Table  1. Similar conclusions can be drawn, i.e., the safe bandgap range calculated by the P-IFEM is the widest, followed by that calculated by the MCM, while the safe bandgap calculated by the A-IFEM is the narrowest. In addition, in order to further quantify the impact of uncertainty on the bandgap, let η denote the proportion of width of the safe bandgap, where f ws and f wd represent the widths of the safe bandgap and deterministic bandgap, re- spectively. Obviously, η is a value between 0 and 1, and the closer the value of η is to 1, the smaller effects the uncertainty has on the bandgap. The proportions of width of the safe bandgaps calculated by the three methods are also listed in Table 1 to show the influence of uncertainties of Young's modulus on the bandgap. It can be seen that for the first bandgap, the uncertainty of Young's modulus reduces the reliable bandgap range by about 10%, while for the second bandgap, the reliable bandgap range decreases by about 5%. Obviously, the uncertainty of Young's modulus has a certain effect on the bandgap. Meanwhile, by comparing the calculation results of P-IFEM and A-IFEM, it can be inferred that the P-IFEM has a higher accuracy in predicting the safe bandgap. For comparison purpose, the safe bandgap considering only the uncertainty of Young's modulus of aluminum, i.e., ∆α i = ∆α Al = 0.05 (i = 1, 3, 5, · · · , 2N u − 1) and ∆α i = ∆α Po = 0 (i = 2, 4, 6, · · · , 2N u ), is shown in Fig. 3(b). The interval vibration transmission analysis results around the bandgap are also shown in Figs. 4(c) and 4(d), and the ranges and proportions of width of the safe bandgaps are listed in Table 1. As can be seen from Fig. 4(c), Fig. 4(d), and Table 1, the ranges of safe bandgaps calculated by the P-IFEM and the MCM are very close, while the range of safe bandgaps calculated by the A-IFEM is narrower. On the other hand, it can be found that the safe bandgap ranges obtained by considering only the uncertainty of Young's modulus of aluminum are all wider than those obtained by considering both the uncertainty of Young's modulus of aluminum and polyethylene. A similar phenomenon can be found that, compared with the upper bounds of the interval transmission coefficient in Fig. 3(a), the upper bounds of the interval transmission coefficient in Fig. 3 Similarly, an interval vibration transmission analysis considering only the uncertainty of Young's modulus of polyethylene, i.e., ∆α i = ∆α Al = 0 (i = 1, 3, 5, · · · , 2N u − 1) and ∆α i = ∆α Po = 0.05 (i = 2, 4, 6, · · · , 2N u ), can also be implemented, and its results are shown in Fig. 3(c). Figures 4(e) and 4(f) show the interval vibration transmission analysis results around the bandgaps, and the corresponding ranges and proportions of width of the safe bandgaps are listed in Table 1. It can be seen from Fig. 4(e), Fig. 4(f), and Table 1 that the safe bandgaps calculated by the MCM are wider than those calculated by the A-IFEM, but narrower than those calculated by the P-IFEM. By comparing the data in Table 1, it can be found that the proportions of width of the safe bandgaps obtained by only considering the uncertainty of Young's modulus of polyethylene are smaller than those obtained by only considering the uncertainty of Young's modulus of aluminum. In other words, compared with Young's modulus of aluminum, Young's modulus of polyethylene has a greater effect on the bandgap.

Phononic crystal beam with interval thickness
The second example in this paper concerns the phononic crystal beam with interval thickness. In this example, all the parameters except for the thickness are the same as those in Subsection 3.1, and the thickness of each beam is treated as an interval parameter. First, let us assume that the radius of the interval variable α I i of each beam's thickness is ∆α i = ∆α Al = ∆α Po = 0.02 (i = 1, 2, · · · , 2N u ). The corresponding interval vibration transmission analysis results are shown in Fig. 5(a). Figures 6(a) and 6(b) show the interval vibration transmission analysis results around the bandgaps, and Table 2 lists the ranges and proportions of width of the safe bandgaps in detail. It can be seen that, under the influence of the uncertainty of thickness, the upper bounds of the interval transmission coefficient deviate significantly from the midpoint line, making the reliable bandgap ranges for vibration reduction gradually decrease. The safe bandgaps calculated by the three methods are slightly different. For the first safe bandgap, the beginning frequency calculated by the A-IFEM is slightly higher than those calculated by the P-IFEM and the MCM, and the ending frequencies predicted by the three methods are in perfect agreement. Moreover, for the second safe bandgap, the safe bandgap range calculated by the P-IFEM is the widest, and the safe bandgap range predicted by the A-IFEM is the narrowest, while the safe bandgap range calculated by the MCM is in between. It can also be obtained that, for the first bandgap, the uncertainty of thickness reduces the reliable bandgap range by about 10%, while for the second bandgap, the reliable bandgap range decreases by about 5%. Obviously, the uncertainty of thickness has a certain effect on the bandgap. Meanwhile, by comparing the calculation results of P-IFEM and A-IFEM, it can be inferred that the P-IFEM has a higher accuracy in predicting the safe bandgap.
For comparison, Fig. 5(b) displays the interval vibration transmission analysis results considering only the uncertainty of thickness of aluminum obtained by setting ∆α i = ∆α Al = 0.02 (i = 1, 3, 5, · · · , 2N u − 1) and ∆α i = ∆α Po = 0 (i = 2, 4, 6, · · · , 2N u ). In Figs. 6(c) and 6(d), the results of the interval vibration transmission analysis around the bandgap are dis-played. In Table 2, the ranges and proportions of width of the safe bandgaps calculated by the P-IFEM and the A-IFEM are compared with those calculated by the MCM. It can be seen from Table 2 that the safe bandgap ranges predicted by the P-IFEM agree well with those predicted by the MCM. On the other hand, for the A-IFEM, the calculation results except for the beginning frequency of the second safe bandgap are basically consistent with those of the MCM. This point can also be seen from Figs. 6(c) and 6(d). Because in Figs. 6(c) and 6(d), not only the upper bounds of the interval transmission coefficient calculated by the three methods are basically consistent, but even the upper bounds are in very good agreement with the midpoint line. In addition, when ∆α Al = 0.02 and ∆α Po = 0, all the proportions of width of the safe bandgaps listed in Table 2 are very close to 1. These phenomena all indicate that the thickness of material aluminum has little influence on bandgap characteristics.
Similarly, the interval vibration transmission analysis can also be performed when only the uncertainty of thickness of polyethylene is considered. Figure 5(c) displays the interval vibration transmission analysis results obtained by setting ∆α i = ∆α Al = 0 (i = 1, 3, 5, · · · , 2N u − 1) and ∆α i = ∆α Po = 0.02 (i = 2, 4, 6, · · · , 2N u ). Figures 6(e) and 6(f) show the interval vibration transmission analysis results around the bandgap, and Table 2 lists the ranges and proportions of width of the safe bandgaps. By inspection of Fig. 6(e) and Table 2, it is observed that the first safe bandgap ranges calculated by the three methods match well. On the other hand, it can be found from Fig. 6(f) and Table 2 that, for the second safe bandgap, the beginning frequency calculated by the P-IFEM is slightly lower than those calculated by the A-IFEM and the MCM, while the ending frequency predicted by the P-IFEM is the highest, the ending frequency predicted by the A-IFEM is the lowest, and the ending frequency predicted by the

Computational efficiency of P-IFEM, A-IFEM, and MCM
To illustrate the computational efficiency of P-IFEM, A-IFEM, and MCM, the time consumption of phononic crystal beam with interval thickness, for ∆α Al = ∆α Po = 0.02, is listed in Table 3. Computations of the three methods are implemented by MATLAB version R2022a on the 3.70 GHz Intel Core CPU i9-10900K. By observing Table 3, it can be found that both the P-IFEM and the A-IFEM spend much less time than the MCM. In other words, the computational efficiency of P-IFEM and A-IFEM is much higher than that of the MCM.

Conclusions
Because uncertainties are unavoidable in phononic crystals, the reliable bandgap range for vibration reduction is actually smaller than the theoretical bandgap range without considering the uncertainties. In this paper, these uncertainties are treated as interval parameters, based on which an interval phononic crystal beam model is established. Meanwhile, the authors of this paper call the bandgap that can always achieve the purpose of vibration reduction no matter how the interval parameters change as safe bandgap. Both the P-IFEM and the A-IFEM are used to obtain the lower and upper bounds of the dynamic response interval vector of this interval phononic crystal beam, based on which an interval vibration transmission analysis can be performed. Finally, two numerical examples are investigated to demonstrate the effectiveness and accuracy of P-IFEM and A-IFEM.
Results show that the reliable bandgap (safe bandgap) range may even decrease by 10%, compared with the deterministic bandgap range without considering the uncertainty. In most cases, the safe bandgap predicted by the P-IFEM is wider than the reference solution (see the MCM result), while the safe bandgap predicted by the A-IFEM is narrower than the reference solution. In other words, the A-IFEM can usually provide a conservative solution to the safe bandgap. It is well known that although the MCM is supposed to provide the correct solution to the safe bandgap, it is very inefficient. Therefore, the A-IFEM is a good choice when both security and computational efficiency must be considered. However, the P-IFEM usually has a higher computational accuracy while ensuring the same computational efficiency as the A-IFEM. Finally, it should be pointed out that only the undamped phononic crystal beam is studied in this paper. Next, the authors of this paper will further study the influence of interval damping on the bandgap characteristics.

Open Access This article is licensed under a Creative Commons Attribution 4.0 International
License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.