Abstract
This paper investigates the scattering behavior of guided shear horizontal (SH) waves in a two-dimensional, isotropic, and linear elastic layered plate with partially debonded interface by analyzing the reflection and transmission coefficients of scattered waves. The partial wave technique is established to form the displacement and stress of guided wave functions, and the boundary element method (BEM) is utilized to handle the numerical calculation with elastodynamic fundamental solutions in the frequency domain. After applying proper boundary conditions including continuity condition on the interface with traction-free debonding, the scattering coefficients can be obtained in terms of boundary element solutions. Two different materials (steel and aluminum) with various debonding lengths and locations in a 1 mm double-layered plate are considered. With several modes of the incident wave over a frequency range up to 4.5 MHz, the variations of scattering coefficients and scattering phenomena are numerically investigated as several parameters such as mode of the incident wave, materials, locations, and length of debonding are changed. The numerical results also suggest the potential of the suitable wave mode for the debonding detection, which can be useful for non-destructive inspection.
1. Introduction
Multi-layered or laminated structures are now widely used in many disciplines of engineering due to their improved mechanical qualities over traditional structures [1, 2]. Carbon fiber sheets or plates, for example, are often employed in infrastructure maintenance to reinforce or repair damaged structures in the civil engineering field [3]. But with more and more complicated layered materials, it is important to develop a non-destructive testing (NDT) method to handle the complexity.
When any structures become aged, a flaw may appear such as a crack, surface damage, and/or debonding along the interface of a layered plate [4, 5]. So, inspection is becoming more necessary to ensure the safety of existing structures as well as new ones. Ultrasonic guided waves have been widely used as NDT for large structures, e.g., plates, pipes, and rail trucks [6–8]. Compared to conventional ultrasonic testing (UT), which can inspect one location at a time, guided wave NDT makes it possible to inspect a large area at the same time due to the waves propagating along the structure geometry. Interaction between incident guided waves and some flaws inside the structure generates scattered waves in all directions, and the scattered waves are used to detect a defect in the hidden zone of the structure before it causes further damage resulting in the failure of the structure. The scattering analysis of guided waves plays an important role for not only the flaw detection but also the flaw evaluation such as the flaw shape reconstruction. Hirose et al. [9–11] proposed an inverse reconstruction method for plate-thinning and inner cavities in plate problems using the numerical simulation of both guide Lamb and SH waves.
In order to calculate the scattering coefficients, many researchers prefer to use numerical techniques due to the complex nature of scattering problems which are too complicated to find the analytical solutions for layered structures. The finite element method (FEM) is the most popular numerical method to this day. Koshiba et al. [12] used the FEM to obtain the solutions of guided Lamb waves for elastic plate problems with internal and surface cracks. Al-Nassar et al. [13] computed the reflection and transmission coefficients by using FEM and Lamb wave modal expression method for a plate with a normal rectangular strip weldment. Hayashi and Rose [14] used a semi-analytical FEM to reduce the computational time and computer resources to simulate the propagation of a guided wave in a plate and pipe. On the other hand, BEM has shown many advantages over the FEM for guided wave scattering analysis since we only need to discretize the boundaries of the domain of interest. The BEM also uses less computational time, memory, and storage than FEM [15]. Rose et al. [16, 17] developed a hybrid BEM with the combination of normal mode expansion technique to calculate the reflection and transmission coefficients of guided Lamb waves for edge reflection and surface breaking problems. Rose et al. [18, 19] further applied the hybrid BEM to both guided Lamb and SH waves for the internal inclusion in a plate and also did the 3-dimensional analysis of guided wave scattering. In these literature reviews, they are mostly concerned about a surface or inner flaw in a single plate. There was little research on debonding problems in a layered plate except for our previous paper on the BEM to analyze the scattering behavior of guided Lamb waves by layered plate debonding [20].
In this paper, the BEM is applied to the scattering analysis of guided SH waves in a layered plate with partially debonded interface. The technical approach is almost the same as in the previous research on the scattering analysis of guided Lamb waves [20]. In Section 2, wave functions of displacement and stress for guided SH waves are formulated using the partial wave technique including the calculation of dispersion curves. Next, in Section 3, the BEM is implemented to formulate the boundary integral equation. After discretization into constant elements and application of proper boundary conditions, scattering coefficients can be obtained in terms of reflection and transmission coefficients defined by ratios between the amplitudes of calculated reflection/transmission waves and given incident waves. In Section 4, numerical results of reflection and transmission coefficients are demonstrated with various configurations such as material constants, mode of the incident wave, and location of the debonding over a frequency range up to 4.5 MHz. This research will provide useful information for guided waves NDT and also serve as fundamental research for more complex waveguide structures such as anisotropic layered plates and pipes.
2. Guided SH Waves in Layered Plate
Consider SH waves in a 2-dimensional, isotropic, and linear elastic double-layered plate with a partially debonding part along the interface, as shown in Figure 1. The thickness of the plate is , and the lateral length is infinite. Material properties are given in terms of the shear modulus and the density , where is the number of layers, i.e., 1 or 2. The top and bottom surfaces of the plate are traction-free, and the interface is perfectly bonded except for the traction-free debonded part with the length . The BEM analysis is carried out for the finite region with the length , which is long enough so that near-field evanescent modes can be neglected on the both ends of the BEM analysis zone. The Cartesian coordinate system is applied with the origin located at the bottom left of the BEM zone.

Guided SH waves have anti-plane displacements in direction. The equations to descript the displacement and stress of guided SH waves can be formulated using the partial wave technique which is derived from Christoffel equations [21, 22]. The modal functions for the displacement of guided SH waves, propagating in the direction in the i-th layer, can be expressed aswhere are the arbitrary amplitudes and is the wavenumber defined as , where and are frequency and phase velocity, respectively. is the ratio of the wavenumbers between and directions, given by
The modal functions of the stresses and can be obtained by applying Hooke’s law for an isotropic elastic solid [23]. Thus, we can get
One important characteristic of guided waves in a plate is a dispersive property, meaning that the phase velocity is dependent on the frequency and the mode of propagation. The exact solution of the phase velocity for the double-layered plate cannot be obtained analytically due to the complexity of the multi-layered plate with different materials. Thus, in order to find the pairs of for each mode, the numerical approach is applied by using the global matrix method combined with a partial wave technique [22]. Applying the following boundary conditions:and solving the eigenvalue problem derived from the system of linear equations (6)–(9), the dispersion curves can be obtained. Note that (6)–(9) are valid for all values of and the exponential term of is extracted as a common item. Thus, equations (6)–(9) give the same equations regardless of the sign .
Figure 2 shows the dispersion curves of a 1 mm steel single plate and a 1 mm double-layered plate made with various thicknesses of steel and aluminum. The material properties are GPa and g/cm3 for steel and GPa and g/cm3 for aluminum. Since the double-layered plate is made of different materials, we cannot separate the mode of propagation into symmetric and anti-symmetric modes. Thus, we named them Mode 1, Mode 2, Mode 3, …, in the order that their cutoff frequencies appear. It is worth noting that Mode 1 of a double-layered plate shown in Figure 2 is also dispersive, unlike a fundamental mode (SH0) in a single-layered regular plate in which the phase velocity is constant.

In order to fix the arbitrariness in the amplitude and make it possible to directly compare the reflection and transmission coefficients of all modes, the modal wave functions are normalized as [16]where is the normalized factor of m-th mode defined as
In (12), and are the average powers of the m-th modal mode and the incident mode passing through out a vertical section of a plate which can be defined aswhere the asterisk means the complex conjugate.
3. Boundary Element Implementation
With use of a standard weighted residual technique with the weight function of the elastodynamic displacement fundamental solution in frequency domain [20], The boundary integral equation for the guided SH wave problem on the surface S is expressed aswhere is the constant and equal to 1/2 in the case of smooth boundary. and are the displacement and traction in the direction at the point , respectively. The traction can be defined as a product of stress and unit normal vector outward the boundary, namely, . The displacement fundamental solution is expressed as
Here, , , , and is the Hankel function of the first kind of zero order. The traction fundamental solution can be obtained by taking the derivative of (15) with respect to a unit outward normal vectorwhere is the Hankel function of the first kind of first order.
Now, boundary integral formulation (14) is applied to two subdomains and surrounded by the boundaries and for each layer, as shown in Figure 3. After discretized into constant elements for the i-th layer, boundary integral equation (14) becomes

Equation (17) gives the total equations with terms of and ( and ). Combining all of equation (17) yields a matrix form:where and are the vectors of the total displacement and traction, respectively, on the boundaries and . is the integral of the displacement fundamental solution, and is the traction counterpart defined as . To make the number of unknown terms equal to , which is the same as the number of the equations in (18), the following continuity and boundary conditions are introduced.
Here, equations (19) and (20) are the continuity conditions of displacement and traction, respectively, on the interface on the assumption that two layers are perfectly bonded together except for the debonding zone. Equations (21) and (22) are the traction-free boundary conditions on the bottom and top surface and , respectively, and also the debonding part on the interface. Now further required condition to solve (18) is the relation between displacement and traction on the artificial boundaries , , , and which will be discussed as follows.
The relation between displacement and traction on the artificial boundaries, and , is derived by using the normal mode expansion technique, in which wave fields are represented as the linear combination of normalized wave function for possible propagating modes. For example, on the left boundary , the total displacement and traction can be expressed as the summation of the incident and reflected waveswhere k and l are the numbers of elements on and , respectively, and . M denotes the maximum number of propagating modes at a certain frequency, e.g., at MHz (see Figure 2). is the given amplitude of the incident wave and is the unknown amplitude of m-th mode of the reflected waves R. For the right boundary, similar expression can be written aswith T indicating the transmitted waves. Next, rearranging (23) and (25) in terms of the unknown scattered amplitudes,where is the generalized complex inverse matrix of (see [24]). It is worth noting that the unknown scattered amplitudes, and , are common for both displacement and tractions, e.g., in (23) and (24) and in (25) and (26). Consequently, we can substitute (27) and (28) into (24) and (26), respectively, to obtain
Now total tractions (29) and (30) are expressed in terms of the total displacements. Then, substituting (29) and (30) into (18), the total displacements are now the only remaining unknowns on the artificial boundaries, and , and hence (18) can be solved by any linear solver.
Finally, substituting the calculated total displacements into (27) and (28) yields the scattered amplitudes of reflection and transmission waves, respectively. Then, the reflection and transmission coefficients can be obtained as the ratios of absolute values between scattered and incident waves:where is the reflection coefficient and is the transmission coefficient for m-th mode with the incident wave of inc mode.
4. Numerical Results
The numerical results of reflection/transmission coefficients are shown to discuss the influences of incident mode, material, location of the interface, and length of the debonding for bi-material plate structures made of steel and/or aluminum. As for the numerical accuracy of calculations, the power balance of scattered waves is validated as shown in the last part of this section. Figure 4 shows the modeling configuration with various combinations of incident mode of Mode 1 or Mode 2, location of interface mm, and length of the debonding mm. The bottom layer is always steel with the thickness , and the top layer with the thickness can be either steel or aluminum depending on the considered problem. Note that the total thickness is fixed as 1 mm. The material properties for steel and aluminum are the same as before. In order to apply the far-field approximation, the lateral length of the model is chosen to be L = 10 mm which is longer than two times of the incident wavelength and the near-field evanescent modes can be neglected [12]. The total number of elements is chosen to be more than ten per incident wavelength to guarantee the convergence of the results. It is also mentioned that the crack-tip singularity on the interface can be neglected for the calculation of scattered coefficients in far field, and no special element is used at the crack tip.

Figures 5 and 6 show the reflection and transmission coefficients versus frequency range from 0–4.5 MHz of a steel-steel plate with 1 mm debonding at the center of the plate subjected to the incident waves of Mode 1 and Mode 2, respectively. “R” and “T” denote the reflection and transmission, respectively, and the number after the symbols “R” and “T” indicates the mode number.

(a)

(b)

(a)

(b)
Figure 5 shows that there is no scattering occurring in the case of Mode 1 incidence for all frequency ranges. This behavior is because Mode 1 in the steel-steel plate, which is equivalent to a single steel plate, is a fundamental mode (SH0) with no dispersive properties and the amplitude is constant across the vertical section, and hence no stress disturbance occurs at the horizontal debonding.
Figure 6 shows the same results as Figure 5, but for the frequency range from 2–4.5 MHz and Mode 2 incidence. In the case of Mode 2, scattered behavior can be seen for the larger frequency than the cutoff frequency around 2 MHz. The resonance phenomena can also be seen at 3.5 and 4.15 MHz where the reflection coefficient becomes nearly unity and most of incident energy is reflected from the debonding. Since material, geometry, and debonding configurations are symmetric and symmetric and anti-symmetric modes in a homogenous single plate are uncoupled, no mode conversion occurs where Mode 1 and Mode 3 remain zero. Furthermore, the results confirmed the power balance so that reflection and transmission coefficients reverse the trend to each other.
Figure 7 shows the reflection and transmission coefficients versus frequency range from 0–4.5 MHz of a steel-aluminum plate with 1 mm debonding at the center of the plate subjected to Mode 1 incidence. Since the uncouple property for a single plate cannot be held for the bi-material plate, all possible modes are generated as the scattered waves by the debonding. Mode conversions from Mode 1 to Modes 2 and 3 occur for the frequencies beyond the cutoff frequencies 1.5 MHz and 3.2 MHz of Modes 2 and 3, respectively. In the case of Mode 1 incidence, however, mode conversions are mild, and there is no large resonance phenomenon in Figure 7.

(a)

(b)
Figure 8 shows the same results as Figure 7, but for the frequency range from 2–4.5 MHz and Mode 2 incidence. In this case, large variations in reflection and transmission coefficients are clearly seen by resonance phenomena at 3.4 and 4 MHz. It is noted that because of mode conversion from incident Mode 2 to other modes occurring for the bi-material plate, the peak values of Mode 2 in Figure 8 are slightly smaller than the values for the single plate as seen in Figure 6.

(a)

(b)
From the results shown in Figures 5–8, it can be said that the incident mode plays an important role to detect the debonding when an SH plate wave is used in the non-destructive testing. Mode 1 incidence cannot be applied to detection of any horizontal flaw in a homogenous plates and is very difficult to detect horizontal debonding even in a double-layered plate, although it is known that in the case of in-plane Lamb waves, the first mode A0 is very effective for the detection of the debonding [20]. Instead, in the case of SH plate waves, Mode 2 incidence is much suitable, since a large amount of scattered waves is generated. Material differences that produce unsymmetric deformation in a plate also induce the mode conversion significantly. In the following, therefore, the influence of location of the interface and length of the debonding will be discussed only for bi-material plates with Mode 2 incidence.
Figures 9 and 10 show the reflection and transmission coefficients versus frequency in the bi-material plate with different interface positions of mm and mm, respectively. The incident wave is Mode 2. From Figures 8–10, it is said that the peaks of the resonance phenomena are shifted to the lower frequency zone as decreases, namely, the aluminum portion becomes large, which is related to the dispersion curves as the phase velocity is changed depending on material thickness (see Figure 2). It is also worth noting that for unsymmetric cases in geometry of mm and mm, the mode conversion from the incident wave of Mode 2 into the transmission of Mode 3 becomes strong, resulting in the large transmission amplitude of Mode 3 as shown in Figures 9(b) and 10(b), compared with the case of shown in Figure 8.

(a)

(b)

(a)

(b)
Figure 11 demonstrates the reflection and transmission coefficients of Mode 2 versus frequency range from 3.6–4.2 MHz where the peak and valley in the coefficients are clearly visible. The lengths of the debonding are varied from to mm. With the increase of the length of the debonding, the peak of the resonance trend shifted to the lower frequency zone. Similar behavior has been also observed for Lamb waves in [25] which is related to the phase interaction from two sides of the debonding tips.

(a)

(b)
Finally, the numerical accuracy of calculations is validated by the concept of power balance. Since all wave functions are normalized by the power of amplitude ratio in (10) and (11), the sum of the squares of the absolute scattering coefficients of all scattered modes should equal 1 at each frequency [26]. Figure 12 shows the error in the power balance for the example shown in Figure 8. It is observed that the error tends to increase as the frequency is increased, and that several error spikes occur in the local frequency ranges around the cutoff frequency at 3.2 MHz and the resonance phenomena at 3.4 and 4.0 MHz. As a whole, however, the errors are less than 2% for all frequency ranges with the average of 0.76%. Hence, it can be said that the numerical results obtained by the BEM analysis are sufficiently accurate for investigating the scattering of guided SH waves.

5. Conclusion
The scattering analysis of SH waves in a double-layered plate was carried out by means of the boundary element method together with the partial wave technique, and the reflection and transmission coefficients were obtained. The validity was checked by the concept of power balance. From the calculation, it was shown that the incident mode, materials, location, and length of the debonding affected the scattered coefficients which can be used as reference information in ultrasonic guided wave inspection. Particularly, it was suggested that Mode 1 incidence was not suitable for detection of debonding due to the lack of scattered waves or mild scattering, while Mode 2 could be used for the NDT of debonding where the mode conversion and resonance phenomena occurred throughout a whole frequency range.
In the future, three-dimensional analysis of layered plate structures will be carried out by combining Lamb and SH waves in a full 3-D space.
Data Availability
The data used to support the findings of this study are included within the article.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of the paper.
Acknowledgments
This study was supported by JSPS KAKENHI (grant no. JP20H02230).