A rapid 2D-FDTD method for damage detection in stiffened plate using time-reversed Lamb wave

Abstract

This study proposes an improved two-dimensional finite-difference time-domain method (2D-FDTD) to simulate the propagation of guided waves in three-dimensional (3D) stiffened plates, which can greatly reduce the calculation amount and accelerate the simulation process. This method can be used for damage-detection imaging based on time-reversal with wavefield reconstruction (TR-RW) as a quick numerical simulation of wavefield is available. First, the piezoelectric wafer active sensors excite and collect Lamb wave signals in the damaged stiffening plate, after which these signals are reversed in the time domain. Second, the TR signals are again triggered by each sensor in turn on the reconstructed numerical model of the tested structure with 2D-FDTD. Third, the wavefield animation is obtained through simulation. Finally, the energy of the wavefield can be focused in certain areas, indicating the location of the damage. As an example, TR-RW, combined with 2D-FDTD method, is adopted to detect damage in an integral grid-stiffened plate, showing that the proposed 2D method can decrease the calculation time from dozens of minutes in the 3D model to dozens of seconds, which is of great significance for damage-detection applications. Furthermore, compared with the virtual TR method, the combined TR-RW and 2D-FDTD method can overcome the multipath effect and improve the damage location accuracy of stiffened structures.

Keywords

Lamb wave time-reversal damage detection wavefield simulation numerical simaulation

Introduction

Widely used in spacecraft, stiffened plates can effectively reduce mass while maintaining structural stiffness.^1,2 Spacecraft in outer space are faced with the risk of space debris or micro meteoroid impact, which may cause perforation, leakage, and disintegration.^3,4 The use of guided waves with piezoelectric wafer active sensors (PWAS) for damage detection on stiffened plates has attracted significant research interest.^5–7 Despite the promising features of the abovementioned application in flat plates, the reflection, transmission, and multimode superposition of guided waves in stiffeners will degrade the effect of the method in stiffened plates, or even rule out the feasibility of using it.

The time-reversal (TR) method is effective for damage detection in complex structures.^8,9 It refers to a process in which the received signal is reversed in the time domain and reloaded to the corresponding sensor, and acoustic energy is focused at the source. The damage location can be determined by observing the distribution of the wavefield energy on the surface. Nevertheless, the actual TR process is unfeasible because the surface vibration at each time node and position is difficult to be monitored overall.

In the last two decades, researchers have found that the TR process can be realized based on the principle of delay-and-sum (DAS), such as the virtual TR (VTR)¹⁰ and TR with multiple signal classification algorithms.^11–13 However, when dealing with Lamb waves in complex structures using the VTR method, multimode and multipath inevitably influence the damage imaging results with significant interference. Therefore, some scholars have proposed a method of model compensation to optimize imaging results. For example, in light of the dispersion characteristics of Lamb waves, analytical calculations are performed in the frequency domain to recompress the dispersion signal.¹⁴ This method is effective when the structure is simple, but it is hardly practical in the case of complex structures, such as multiple damages and stiffeners.

The TR with a reconstructed wavefield (TR-RW) realized by numerical simulation is another way to realize damage detection. A rapid Lamb wavefield simulation method must be developed to fully exploit the advantages of the RW-TR method. This method should not only provide the calculation results rapidly and accurately describe the dispersion properties of Lamb waves but also be convenient for damage detection in stiffened structures.

In recent decades, many numerical methods have been proposed to understand the propagation characteristics of Lamb waves in complex structures, as well as their interaction with damage. The most popular multiphysics simulation method for Lamb wave propagation with piezoelectric transducers is the finite element method (FEM),^15–17 which is supported by many commercial software, such as the ABAQUS,^17,18 ANSYS,¹⁹ and COMSOL.²⁰ However, the linear FE need excessive spatial discretization to accurately solve the high-frequency wave propagation problem because their linear interpolation functions cannot precisely represent deformed shapes due to wave propagation within an element. Consequently, the linear FE for three-dimensional (3D) wave propagation analysis requires considerable computation and memory.

In order to accelerate the calculation of Lamb wave propagation simulations, many methods have been proposed, such as the combined analytical finite elements approach (CAFA),²¹ spectral element method (SEM),^22–24 boundary element method,^25,26 local interaction simulation approach,²⁷ and semi-analytical finite element method.²⁸ The CAFA uses a global analytical solution to simulate wave generation, propagation, scattering, mode conversion, and detection, while the wave-damage interaction coefficients are extracted from the harmonic analysis of the local FEM, which speeds up the calculation; however, its modeling is more complex. The SEM is similar to the FEM in that it divides the whole domain into small pieces to solve differential equations, and it has been proven to significantly accelerate the calculation, and well describe the characteristics of the lead zirconate titanate (PZT) and paste layer. However, the 3D-SEM method has not been widely used, mainly because the element division of the complex 3D structure is difficult. It is temporarily limited to the general flat plate structure, or the cross section of the stiffened plate.

The finite-difference time-domain (FDTD) method is one of the most commonly used methods for studying elastodynamics.²⁹ It solves differential equations by using a finite difference method instead of a derivation. The same equation of each grid node in the FDTD method makes it programmer-friendly, and parallel computing can be used to reduce the calculation time. This practical algorithm can easily observe the conversion process of various waveforms during ultrasonic propagation, and even fit for the nonlinear elastic media.²⁹ In addition, a perfectly matched layer can be used to obtain the boundaries of total absorption, semi-absorption, semi-reflection, and complete reflection.³⁰

However, similar to the FEM, the 3D-FDTD method is also inefficient. In many studies, scholars have used two-dimensional finite-difference time-domain method (2D-FDTD) to simulate the acoustic wave propagation of the plate cross section and obtain the dispersion curve, or to study other characteristics of acoustic wave propagation in the structure.³¹ The 2D simulation method enjoys a fast computing speed, and it is much easier to model the tested structure. However, the 2D method can only simulate the wave propagation in the cross section of the plate, or the body wave propagation without dispersion characteristics. For the Lamb wave, the general 2D numerical simulation method can only simulate the wave propagation of the cross section of the plate, but cannot simulate the out-plane wavefield of the whole plate.

The Lamb wave is a unique waveform in the thin plate. It is formed by the repeated reflection and superposition of P wave and SV wave on the upper and lower surfaces of the plate, which has the characteristics of multimode and dispersion. If the general 2D numerical simulation method is used, only the section of the plate can be simulated,³² while the entire 3D structure cannot be simulated simultaneously.^22,23 In order to obtain the wavefield of the plane on the whole structure, the 2D-FDTD method needs to be improved.

The remainder of this paper is organized as follows. In the “Numerical simulation method” section, the 2D-FDTD method for simulating Lamb propagation in stiffened plates is presented. The 2D-FDTD is verified by the scanning laser Doppler vibrometer (SLDV) in the “Method verification” section. Then, Lamb wave propagation in a plate with a single stiffener is simulated using both 2D-FDTD and FEM for comparison. In the “Damage detection with combined TR-RW and 2D-FDTD method” section, simulations are performed to demonstrate how TR-RW method works. The damage detection results of the TR-RW are compared with those of the DAS method to demonstrate the ability of the former method to overcome multipath and multimode interference.

Numerical simulation method

Lamb wave simulation in 3D-stiffened plates using 2D-FDTD method

The FDTD method can be used to simulate the propagation of ultrasonic waves in elastic solids. The wave equation in 2D elastic media can be expressed as follows:

\frac{\partial^{2} u}{\partial t^{2}} = c^{2} (\frac{\partial^{2} u}{\partial x^{2}} + \frac{\partial^{2} u}{\partial y^{2}}),

(1)

where $u$ is the amplitude of the wavefield, $c$ is the wave propagation velocity, $x$ and $y$ are the two Cartesian coordinates, and $t$ is the time. The model is divided into an $N_{x} \times N_{y}$ rectangular grid to obtain a numerical solution. Let $Δ x$ and $Δ y$ be the distance steps on the $x$ and $y$ axes, respectively, and $Δ t$ be the time increment. Then, the second-order central difference method is used to approximate the continuous partial derivative to obtain the discretization equation

\begin{matrix} \frac{u_{i, j}^{t_{n + 1}} - 2 u_{i, j}^{t_{n}} + p_{i, j}^{t_{n - 1}}}{Δ t^{2}} = c_{i, j}^{2} \\ (\frac{u_{i + 1, j}^{t_{n}} - 2 u_{i, j}^{t_{n}} + u_{i - 1, j}^{t_{n}}}{Δ x^{2}} + \frac{u_{i, j + 1}^{t_{n}} - 2 u_{i, j}^{t_{n}} + u_{i, j - 1}^{t_{n}}}{Δ y^{2}}), \end{matrix}

(2)

where $c_{i, j}$ denotes the wave speeds at positions $(i, j)$ . By solving the wavefield values at position $(i, j)$ at time $t_{n + 1}$ , $u_{i, j}^{t_{n + 1}}$ , a finite difference formula to approximate the solution of the 2D acoustic equation can be obtained. Each point in the grid can be updated independently at each iteration according to the explicit formula above. Owing to limitations in the computer memory and calculation speed, FDTD calculations can only be performed in a limited area. To simulate the propagation of waves over a larger area, it is necessary to set up absorbing boundary conditions to absorb the reflected waves and avoid their generation.

The propagation speed of the waves varies with the frequency, that is, $e^{- i ω t} = e^{- i ω d / c (ω)}$ . Based on Huygens’ principle, the dispersion characteristics of Lamb waves can be simulated by superimposing multiple wavefields with different grid velocities. However, the greater the number of superimposed wavefields, the higher the calculation amount. Therefore, choosing a narrowband excitation signal helps to reduce calculation amount. Additionally, it is preferable to select a single-mode wave for damage detection based on the sensors and structures used.

In a given material and plate size, the phase velocity of the Lamb wave, $c (ω, d)$ , is determined by the frequency and thickness, and can be solved by the Rayleigh-Lamb wave equation.³³ When the frequency, structural material, and size remain unchanged, the velocity is determined only by thickness. The relationship between thickness and wave velocity can be uniquely obtained from the Lamb wave dispersion curve. In contrast, the difference in the wave velocity can be used to simulate changes in the structural thickness.

In this study, only the integral-stiffened plate is considered; that is, the stiffener and the main structural plate are made of the same material, and there is no welding or bonding discontinuous faces. This type of stiffener plate, which is made of the same material, can be regarded as the splicing of two plates with different thickness. Figure 1 shows a typical stiffened plate structure and transmission path of guided waves passing through a stiffener. The position matrix of the stiffeners is constructed as follows:

T_{i, j} = {\begin{matrix} 1, & with stiffener \\ 0, & without stiffener \end{matrix}

(3)

where the positions with and without stiffeners were 1 and 0, respectively.

For the sake of clarity, the cross section part in the blue dotted box in Figure 1 is redrawn in Figure 2. The thickness of the plate is $2 d$ , the height of the stiffener is $h$ , and the width is $b$ . A transmitting sensor at coordinate (0, 0) converts the electrical signal into a stress wave which simultaneously propagates to the left and right in the plate structure. Only the main waveforms propagating from left to right are considered. The position of the stiffener is $(x_{s}, 0)$ , and the two observation points are located on both sides of the stiffener at the coordinates of $(x_{1}, 0)$ and $(x_{2}, 0)$ , respectively.

Note that the wave velocity in the plate is $c (ω, 2 d)$ , and the wave velocity in the stiffener is $c (ω, 2 d + h)$ for the same mode. This assumption can be only made at a low frequency, otherwise the Lamb wave may change to Rayleigh wave. For the stiffener plate several centimeters high, the inquiring frequency should be within several hundred kilohertz to ensure that only zero-order mode Lamb waves exist. Thus, we can obtain the wave velocity matrix:

c_{i, j} = c (ω, 2 d + h) T_{i, j} + c (ω, 2 d) (1 - T_{i, j}) .

(4)

As shown in Figure 2, the stiffener divides the incident wave into three parts when it propagates from left to right. The first part of the wave is reflected to the left at the bottom of the stiffener, the second part continues to propagate to the right through the stiffener, and the third part propagates upwards. When a wave propagating upward meets the stiffener top, it is reflected to the bottom of the stiffener and then divided into two parts: one to the left and the other to the right. If observation point 1 is set on the left side of the stiffener, its time-domain waveform appears as a superposition of three wave packets. When observation point 2 is set on the right side of the stiffener, its time-domain waveform appears as a superposition of two wave packets.

Reflection, transmission, and absorption occur when guided waves propagate on different media or geometrically discontinuous surfaces. Without the consideration of environment radiation, there won’t be any energy loss, then

σ_{R} + σ_{I} + σ_{A} = 1,

(5)

where $σ_{R}, σ_{I}, σ_{A}$ represent the reflection, transmission, and absorption coefficients, respectively. The reflection coefficient at normal incidence is defined as³³

R_{12} = \frac{σ_{R}}{σ_{I}} = \frac{W_{2} - W_{1}}{W_{1} + W_{2}},

(6)

Transmission coefficient is

T_{12} = \frac{2 W_{2}}{W_{1} + W_{2}},

(7)

where the subscripts 1 and 2 represent different materials and geometric areas, respectively. For stiffened plates of isotropic materials, the acoustic impedance is related only to the wave velocity. The stiffeners change the thickness of the medium and then the velocity of the Lamb wave, defining the acoustic impedance of media with different thickness:

W = ρ c (ω, d) .

(8)

Then the reflection and transmission coefficient can be simplified as

R_{12} = \frac{c_{2} - c_{1}}{c_{1} + c_{2}}, T_{12} = \frac{2 c_{2}}{c_{1} + c_{2}} .

(9)

From the above equation, the relationships among the reflection coefficient, transmission coefficient, frequency, and thickness can be obtained. Because the height of the stiffener is significantly different from that of the flat plate, the wave propagation speed changes evidently. When a numerical method is used to simulate the Lamb wave propagation, the stiffener can be simulated by providing different velocities to nodes in different areas.

The new wavefield is the sum of the original wavefield and the wavefield scattered by the top surface of stiffeners:

U_{total} = A_{0} U_{0} + A_{T} U_{T} .

(10)

The amplitude of the wavefield scattered by the top surface should be determined. The energy of the wavefield scattered by the top surface originates from the incident part, which enters the stiffener upward. After being reflected by the top surface, it is split into two parts at the root of the stiffener.

For a typical system of the excitation source structure to be detected by the receiving sensor, the spectrum of the excitation signal is expressed as $S (ω)$ , and the spectrum of the received signal can be described as

H (ω) = S (ω) T_{x} (ω) A (ω) X (ω) Δ (ω) R_{x} (ω),

(11)

where $T_{x} (ω)$ is the output signal of the excitation source, $A (ω)$ is the attenuation, $X (ω)$ is the reflection or transmission coefficient, and $Δ (ω)$ is the time delay caused by the propagation distance:

Δ (ω) = e^{- i ω τ} = e^{- i ω \frac{d}{c (ω)}} = e^{- ik (ω) d},

(12)

where $τ = \frac{d}{c (ω)}$ represents the transmission delay between two points and the wave velocity $c (ω)$ is related to the signal frequency. When the propagation speed changes with the frequency, the signal waveform is deformed, that is the dispersion phenomenon. With respect to sensor characteristics,

T_{x} (ω) = I (ω) D_{F} (ω),

(13)

where $I (ω)$ indicates the frequency characteristics of the sensor, and $D_{F} (ω)$ is the directional function of the sensor. For circular sensors, it is the same in all directions. The signal at observation point 1 can then be expressed as the sum of the direct signal, reflection signal from the bottom of stiffeners, and reflection signal from the top surface.

Δ_{1} (ω) = A_{I} e^{- ik x_{1}} + A_{R} e^{- ik (2 x_{s} + x_{1}) + π} + A_{R, T} e^{- ik (2 x_{s} + x_{1} + 2 h) + π} .

(14)

In this formula, the factor is introduced into the last two reflection terms because the reflection causes a phase reversal of 180°. The signal at observation point 2 is expressed as the sum of the transmission signal and reflection signal from the top surface:

Δ_{2} (ω) = A_{I} e^{- ik x_{2}} + A_{R, T} e^{- ik (x_{2} + 2 h) + π} .

(15)

Step-discontinuous geometries lie in the left and right sides of the stiffener, and they produce their own reflection and transmission waves. If the width of the stiffener is $b$ and the signal wavelength is $λ$ , the phase difference of the reflection wave generated by the left and right sides of the stiffener is

Δ θ = θ_{l} - θ_{r} = 2 π \frac{b}{λ},

(16)

To simplify the analysis, the wavefield of the hole area is marked as $(x_{d}, y_{d})$ , that is,

u_{i = xd, j = yd} \equiv 0,

(17)

With this definition, the obstacle will be like a complete reflector.

Simulation of Lamb wave emission sensor

To obtain the waveforms of different modes, it is necessary to use multiple source functions to initiate waves of different modes. Then, the multimodal Lamb wavefield is obtained by wavefield superposition:

u (x, y) = \sum_{m = 1}^{M} [α^{A_{m}} u^{A_{m}} (x, y) + α^{S_{m}} u^{S_{m}} (x, y)],

(18)

where $M$ is the number of Lamb wave modes considered, $α^{A_{m}}$ and $α^{S_{m}}$ are the amplitudes of the mth antisymmetric and symmetric modes, respectively.

The velocities used in different modes can be calculated using the Lamb wave equation. The amplitudes are calculated using an analytical method. The upper surface of the left side of the plate is pasted with a piezoelectric sheet with a thickness of $d_{a}$ , with radius $a$ and an elastic modulus of $E_{a}$ , as the excitation source. By applying voltage $v_{T} (t)$ to it, the Lamb wave can be excited in the plate, and its spectrum is recorded as $V_{T} (ω)$ . When the excitation frequency is given, multiple symmetric (S) and asymmetric (A) wavenumbers can be obtained using the Rayleigh Lamb wave equation. The structure transfer function is multiplied by the spectrum of the excitation signal to obtain the spectrum of the received signal,

U (ω, r) = V_{T} (ω) [G^{S} (ω, r) + G^{A} (ω, r)] .

(19)

When a circular piezoelectric patch is used, the frequency characteristics of sensors of different sizes and materials can be calculated according to Victor’s research,^21,34

\begin{matrix} G^{S} (ω, r) = - π i \frac{a^{2} κ_{PZT} (ω)}{2 μ} \sum_{ξ^{S}} \frac{J_{1} (ξ^{S} a) N_{S} (ξ^{S})}{{D'}_{S} (ξ^{S})} H_{1}^{(1)} (ξ^{S} r) \\ G^{A} (ω, r) = - π i \frac{a^{2} κ_{PZT} (ω)}{2 μ} \sum_{ξ^{A}} \frac{J_{1} (ξ^{A} a) N_{A} (ξ^{A})}{{D'}_{A} (ξ^{A})} H_{1}^{(1)} (ξ^{A} r), \end{matrix}

(20)

where

\begin{matrix} N_{S} (k) = k η_{S} (k^{2} + η_{S}^{2}) \cos η_{P} d \cos η_{S} d \\ N_{A} (k) = - k η_{S} (k^{2} + η_{S}^{2}) \sin η_{P} d \sin η_{S} d \\ D_{S} = {(k^{2} - η_{S}^{2})}^{2} \cos η_{P} d \sin η_{S} d + \\ 4 k^{2} η_{P} η_{S} \sin η_{P} d \cos η_{S} d \\ D_{A} = {(k^{2} - η_{S}^{2})}^{2} \sin η_{P} d \cos η_{S} d + \\ 4 k^{2} η_{P} η_{S} \cos η_{P} d \sin η_{S} d \end{matrix},

(21)

where $r$ is the propagation distance of the wave, $τ_{a}$ is the shear stress between the transducer and the host structure, and $μ$ is the shear modulus of the plate structure. $J_{1}$ is a first-order Bessel function that dominates the modulation effect between the piezoelectric sheet and host structure. $H_{1}^{(1)}$ is a first-order Hankel function of the first type that represents the outward propagation of the 2D wavefield. $κ_{PZT} (ω)$ is the transfer function between the piezoelectric sheet and host structure, that is, the coefficient converting voltage into shear stress, $τ_{a}$ , which can be given by

κ_{PZT} (ω) = \frac{d_{31}}{S_{11}^{E}} \frac{χ (ω)}{1 - χ (ω)},

(22)

where $χ (ω)$ is the stiffness ratio between the host structure and excited PZT. $d_{31}$ is the coupling coefficient of the piezoelectric materials, which can be available in relevant manuals. $S_{11}^{E}$ is the in-plane transverse electric field of the piezoelectric wafer.

Method verification

Verification by SLDV

An experiment is conducted on a 600-mm long, 600-mm wide, and 2.5-mm thick aluminum-alloy plate with a single stiffener to validate the proposed 2D-FDTD model, as shown in Figure 3. The elastic wave is generated by a piezoceramic wafer, and an LV-SC500 SLDV is used to measure the vibration velocity of plane surface. An integral stiffener with the height of 20 mm and width of 4 mm is placed in the center of the plate.

The piezoceramic sensors used are a disk with the thickness of 0.5 mm and the diameter of 8 mm. The peak voltage of the excitation is kept at 100 Vpp. In the experiment, the PZT actuator excites a five-peak tone-burst wave with a center frequency ranging from 100 to 350 kHz. According to the Lamb wave dispersion equation, with this frequency and thickness, there are only fundamental symmetric and antisymmetric modes, as the excitation frequency is lower than the cutoff frequency of the first antisymmetric mode.

As shown in Figure 3, the collected analog signals are converted into digital signals and transmitted to the system control module for data analysis and processing. An Aigtek ATA-4052 power amplifier is used to provide an excitation of 100 Vpp. The line scanning mode is used to collect the vibration signal along the propagation path at a distance of 200 mm. The scanning point was set every 1 mm for a total of 200 points.

Figure 4 presents the time-space diagram obtained by the SLDV at 200 kHz. The magnitude of the normal component of the displacement is plotted in time from t = 0 to t = 150 µs and in space from the source at x = 0 to the stiffener at x = 96.5 mm and beyond to 193 mm from the source. The S0 and A0 modes are both visible in the plotted data. When stronger but slower incident A0 waves meet the stiffener, there is significant reflection and transmission of energy with little conversion. Through the experiment, the transmission and reflection coefficients of each mode are obtained, which are consistent with the simulation.

Figure 5 shows the velocity waveforms of the vertical plates on the left and right sides of the stiffeners measured by the laser vibrometer and gives the simulation results for comparison. The two observation points are, respectively, located at 50 mm on the left and right sides of the stiffener, as shown in Figure 3(b). The waveforms are normalized to better observe the fitting degree of the two waveforms. It can be seen from Figure 5(a) that S0 mode and A0 mode overlap because of the relatively short propagation distances. In addition, the amplitude of A0 mode is much larger than that of S0 mode, because the laser vibrometer mainly measures the out-plane displacement of the vertical plate, and the out-plane displacement is dominated by the antisymmetric mode.

It can be observed that the A0 mode of the two signals matches well, but there is a small error in the amplitude and phase of the A0 mode. The error is mainly caused by the fact that the linear element cannot ensure an accurate representation of the A0 mode behavior considering the large mesh size of the mode. Although the waveforms are not completely identical, the positions and widths of the wave packet are fairly the same, which is important for damage imaging.

A comparison with FE method

As a benchmark numerical simulation method, finite element analysis (FEA) has been proved effective in analyzing piezoelectric components and guided wave propagation.^15–17 In order to evaluate the performance of the proposed method, a plate with a single-integral stiffener is modeled using both the FEA and the 2D-FDTD methods. The size and material are consistent with those of the experimental model illustrated in the previous section.

The overall structure of the numerical simulation model is illustrated in Figure 6. The material properties are as follows: Young’s modulus of 68.9 GPa, mass density of $2700 kg / m^{3}$ , and Poisson’s ratio of 0.33.

The stiffener and plate are meshed using a 3D block element with a reduced integration (C3D8R), with a size of 1 mm along the length and width dimensions. By contrast, the thickness direction is meshed using a more refined element (0.5 mm) to ensure at least five elements. The boundaries are meshed with 8-node linear one-way infinite brick (CIN3D8) to eliminate unwanted side reflections. Six PWAS with a diameter of 8 mm are located on the top surface of the plate and can excite or receive Lamb wave signals radially outward. The model is excited by a five-cycle Hanning windowed pulse of 200 kHz. According to the dispersion curve of Lamb wave, only the S0 and A0 modes exist under these excitation frequencies. The total duration of the simulation is 500 µs. An explicit dynamic analysis with a fixed step size of $10^{- 7}$ s is used for all numerical cases.

Figure 6 shows the snapshots of the wavefield animation captured by both methods. When the guided wave meets the stiffener, the energy will be reflected, transmitted, and trapped in the stiffener. It can be seen that there are reflection waves from the top of the stiffener, but it is weaker than the reflections at the bottom of the stiffener because of the energy dissipation of stiffener. On the other hand, the superposition of these two waveforms results in the broadening of waveforms, and this can be used as another proof for the existence of top surface reflection wave.

Figure 7 shows the normalized waveforms of the pristine plate, damage plate, and their difference signals. It is obvious that the waveforms obtained by the two methods are in good consistency, which shows that the proposed method is accurate and efficient in simulating the integral forming of stiffened plates.

Table 1 gives a detailed comparison of these two methods. The 2D-FDTD method significantly reduces the computational time from 51 min to 23 s. It should be noted that 51 min is unacceptable for a rapid damage-detection applications.

Table 1.

Comparison of simulation performance between the FEA and the proposed 2D-FDTD methods.

	3D-FEA	2D-FDTD
Element type	C3D8R and CIN3D8	2D nodes
Element size	1 $\times$ 1 $\times$ 0.5 mm	1.5 $\times$ 1.5 mm
Solver	ABAQUS explicit	MATLAB
Time step	$0.1 μ s$	$0.1 μ s$
Total duration	$500 μ s$	$500 μ s$
Number of nodes	2,367,642	160,000
Time cost	51 min	23 s

2D-FDTD: two-dimensional finite-difference time-domain method; 3D-FEA: three-dimensional finite element analysis.

Damage detection with combined TR-RW and 2D-FDTD method

Method description

The use of the TR method in detecting structural damage requires placement of sensors on the structure to be tested, and each sensor should be capable of sending and receiving signals. The damage-detection system is developed using a pitch–catch setup, in which one of the transducers acts as an actuator to generate Lamb waves and other transducers receive the propagating Lamb waves.

The damage-detection method TR-RW mainly includes two parts: collecting the damage signals in the structure to be tested, and reconstructing the wavefield as well as detecting the damage via numerical simulation. As shown in Figures 8 and 9, this process mainly includes five steps. First, the excitation sensor excites the Lamb wave in the plate. Second, the receiving sensor collects the guided wave on the surface of the plate, which contains the scattering signal of the damage. Third, the damage scattering signal is obtained by subtraction between the baseline signal before damage and after damage, and the wavefield can be reconstructed with TR signals via the 2D-FDTD method. Fourth, the image with damage focus which results from the continuously changing animation is located. Finally, the image is post-processed, filtered, and the location and shape of the damage are output. The first two steps are carried out on an actual plate to obtain the desired signal. The signals are then post-processed to achieve Lamb wavefield reconstruction via numerical simulation, as described in the last three steps.

Multipath effect caused by stiffener with TR

During propagation in all directions, the Lamb wave reflects, diffracts, and refracts in the structure when encountering boundaries, stiffeners, and damages, which leads to the generation of multiple paths when propagating from the actuator to the receivers, also the so-called multipath effect of the guided waves.

Based on a simplified model with a single stiffener as shown in Figure 10, the signal sent by sensor $i$ propagates along two paths, and two wave packets are received at $j$ .

Figure 10 shows the influence of the multipath characteristics of Lamb waves on the TR process. Because the propagation distance of the direct path is shorter than that of the along-stiff path, these two wave packets, direct wavepacket $D$ and reflected wavepacket $R$ , are separated in the time domain when they arrive at the receiver position. After the signal is time-reversed and re-excited, each wave packet generates signals of the $D$ and $R$ packets for a total of four wave packets. Two wave packets along the original path are superimposed to form the main lobe, and they simultaneously arrive at the receiving sensor after the TR to form the main lobe: $RR + DD$ . The main lobe is superimposed by two wave packets in the same phase, and its modal dispersion is well compensated. However, the part that does not return according to the original path forms two symmetrical side lobes.

Affected by multimode and multipath effects, the reconstructed signal of the TR contains a main lobe with a large amplitude and several weak side lobes. The shape of the main lobe is consistent with that of the original excitation signal, and the TR method can still realize signal reconstruction.

Damage detection of grid-stiffened plate using the TR and 2D-FDTD method

In this section, the damage to the grid-stiffened plate is simulated and detected. The plate is a 600-mm square aluminum plate with a thickness of 2.5 mm. This plate was made of a large aluminum block by cutting to ensure the strength of the structure and its parameters are the same as those of real aerospace structures. The height and width of the stiffeners are 20 and 4 mm, respectively. There are two stiffeners, placed horizontally and vertically, with a spacing of 200 mm. The entire panel is divided into nine grids. A sensor is pasted at the center of each grid, and signals can be sent and received. In the way of working in turn, when one sensor sends a signal, other sensors receive it. Thus, a total of nine experiments are performed.

Figure 11 shows the damage-detection results for the grid-stiffened panel obtained using FEA. The TR wavefield produces an obvious focus on the position of a 5 mm × 5 mm square hole. The wave used belongs to the A0 mode. The wavelength is 15 mm, and the defect size is only one-third of the wavelength. Although the FEM can be used to show the damage wavefield image, this method takes too much time to model and calculate, and it is very inconvenient to extract the wavefield image of each frame for subsequent processing.

Then, the wavefield of the 3D-stiffened plate is simulated by 2D-FDTD numerical method. The wave velocities at the flat plate and stiffener are separately calculated using the Lamb wave equation. The reflection and transmission coefficients are corrected based on the experiments. Finally, the focused TR damage-detection results are obtained.

The wavefield image here is filtered, which is easy to achieve for regular grids. In contrast, the results obtained using commercial FE software, like the wavefield cloud diagram in Figure 11, consist of too many nodes and disordered numbers, making post-processing difficult.

Figure 12 presents the baseline, damage, difference, and TR waveforms. As an example, only the waveforms when sensor 1 transmits and sensor 2 receives are displayed here. The baseline wave is collected in the pristine plate, and the damage wave in the plate with hole-type damage. It can be observed that there is little difference between the baseline and damage waves. Therefore, they were subtracted to get the difference signal and then conduct damage detection. This method is called baseline damage imaging. This study does not use a non-baseline method, because it requires high-quality signals. However, excessive number of scattering signals caused by stiffeners in plates submerge the scattering signals of damage, making them undetectable.

Figure 13(a) presents the damage-detection results with conventional DAS or VTR, while Figure 13(b) shows the results using the proposed 2D-FDTD combined with TR-RW methods. The black cross marks the real damage location. It can be seen that the red area representing the probable damage location does not appear in the VTR imaging results. Instead, it is relatively scattered in a large area. The DAS method cannot locate the damage owing to the interference of the stiffeners, but the TR method works well. The location of the maximum probability value is consistent with the actual damage location. This shows that the TR-RW method overcomes the multipath scattering effect caused by stiffeners.

Conclusion

In this study, an improved 2D finite-difference time-domain (2D-FDTD) method is proposed to simulate the propagation of Lamb waves in 3D-stiffened plates. This method can significantly reduce the calculation amount and accelerate the simulation speed. As an application of the fast numerical simulation of the wavefield, 2D-FDTD can be applied to TR Lamb wave damage-detection imaging. Also, the TR damage detection of a plate structure with grid stiffeners is performed. Compared to the VTR method, the results show that the combined TR-RW and 2D-FDTD method can effectively overcome the multipath effect caused by complex structures and improve the damage location accuracy of stiffened plates.

Footnotes

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This study was funded by the high-level innovative talent project of the National University of Defense Technology.

ORCID iD

Sunquan Yu

References

Reusser

Chimenti

Roberts

, et al. Guided plate wave scattering at vertical stiffeners and its effect on source location. Ultrasonics 2012; 52(6): 687–693.

Reusser

Holland

Chimenti

, et al. Reflection and transmission of guided ultrasonic plate waves by vertical stiffeners. J Acous Soc Am 2014; 136(1): 170–182.

Mespoulet

Héreil

Abdulhamid

, et al. Experimental study of hypervelocity impacts on space shields above 8 km/s. Procedia Eng 2017; 204: 508–515.

Chen

Giambene

Yang

, et al. Analysis of inter-satellite link paths for leo mega-constellation networks. IEEE Trans Veh Technol 2021; 70: 2743–2755.

Sharif-Khodaei

Ghajari

Aliabadi

. Determination of impact location on composite stiffened panels. Smart Mater Struct 2012; 77(10): 662–686.

Han

Yang

Wang

, et al. The research on propagation characteristics of acoustic emission signals in stiffened plates based on the multipath propagation model. Ultrasonics 2020; 108: 106177.

Radzieński

Doliński

Krawczuk

, et al. Damage localisation in a stiffened plate structure using a propagating wave. Mech Syst Signal Process 2013; 39(1): 388–395.

Fink

. Time reversal of ultrasonic fields. I. Basic principles. IEEE Trans Ultrason Ferroelectr Freq Control 1992; 39(5): 555–566.

Thomas

Fink

. Time reversal of ultrasonic fields. Il. Experimental results. IEEE Trans Ultrason Ferroelectr Freq Control 1992; 39(5): 567–578.

10.

Wang

Shen

. An enhanced lamb wave virtual time reversal technique for damage detection with transducer transfer function compensation. Smart Mater Struct 2019; 28(8): 085017.

11.

Yuan

. Lamb wave-based subwavelength damage imaging using the dort-music technique in metallic plates. Struct Health Monit 2016; 15(1): 65–80.

12.

Labyed

Huang

. Ultrasound time-reversal music imaging with diffraction and attenuation compensation. IEEE Trans Ultrason Ferroelectr Freq Control 2012; 59(10): 2186–200.

13.

Zuo

Yang

, et al. Damage identification for plate-like structures using ultrasonic guided wave based on improved music method. Compos Struct 2018; 203: 164–171.

14.

Wang

Yin

, et al. A focusing music algorithm for baseline-free lamb wave damage localization. Mech Syst Signal Process 2022; 164: 108242.

15.

Terrien

Osmont

Royer

, et al. A combined finite element and modal decomposition method to study the interaction of lamb modes with micro-defects. Ultrasonics 2007; 46(1): 74–88.

16.

Willberg

Duczek

Perez

, et al. Comparison of different higher order finite element schemes for the simulation of lamb waves. Comput Methods Appl Mech Eng 2012; 241: 246–261.

17.

S-q

Fan

C-g

Gao

, et al. A finite element simulation method of piezoelectric wafer active sensors with mechanical elements. In: 2020 15th symposium on piezoelectrcity, acoustic waves and device applications (SPAWDA), Zhengzhou, China, 2021, pp. 348–352.

18.

Song

Yang

. Noncontact super-resolution guided wave array imaging of subwavelength defects using a multiscale deep learning approach. Struct Health Monit 2021; 20(4): 1904–1923.

19.

Haider

Poddar

Giurgiutiu

. Experimental validation of an analytical method to predict lamb wave scattering from a discontinuity. Smart Mater Struct 2018; 28(1): 015012.

20.

Qiu

Yan

Lin

, et al. Multiphysics simulation method of lamb wave propagation with piezoelectric transducers under load condition. Chin J Aeronaut 2019; 32(5): 1071–1086.

21.

Shen

Giurgiutiu

. Combined analytical fem approach for efficient simulation of lamb wave damage detection. Ultrasonics 2016; 69: 116–128.

22.

Kim

Chang

. Time-domain spectral element method for built-in piezoelectric-actuator-induced lamb wave propagation analysis. AIAA J 2008; 46(3): 591–600.

23.

Chang

. Optimizing a spectral element for modeling PZT-induced lamb wave propagation in thin plates. Smart Mater Struct 2009; 19(1): 015015.

24.

Wang

. Accurate modeling of PZT-induced lamb wave propagation in structures by using a novel spectral finite element method. Smart Mater Struct 2014; 23(9): 095018.

25.

Cho

Rose

. An elastodynamic hybrid boundary element study for elastic guided wave interactions with a surface breaking defect. Int J Solids Struct 2000; 37(30): 4103–4124.

26.

Galan

Abascal

. Lamb mode conversion at edges. A hybrid boundary element-finite-element solution. J Acoust Soc Am 2005; 117(4): 1777–1784.

27.

Shen

Cesnik

CES

. Local interaction simulation approach for efficient modeling of linear and nonlinear ultrasonic guided wave active sensing of complex structures. J Nondestruct Eval Diagn Progn Eng Syst 2017; 1(1): 011008.

28.

Rose

. The semi-analytical finite element method. Cambridge: Cambridge University Press, 2014, pp. 135–154.

29.

Matsuda

Biwa

. A finite-difference time-domain technique for nonlinear elastic media and its application to nonlinear lamb wave propagation. Jpn J Appl Phys 2012; 51(7S): 07GB14.

30.

Chagla

Cabani

Smith

. Perfectly matched layer for FDTD computations in piezoelectric crystals. In: IEEE ultrasonics symposium, 2004, Montreal, QC, Canada, 2005, vol. 1, pp. 517–520. IEEE.

31.

, et al. Wideband dispersion reversal of lamb waves. IEEE Trans Ultrason Ferroelectr Freq Control 2014; 61(6): 997–1005.

32.

Sun

, et al. Guided wave propagation analysis in stiffened panel using time-domain spectral finite element method. Chin J Aeronaut 2022; 35(10): 208–221.

33.

LRose

. Ultrasonic guided waves in solid media. Cambridge: Cambridge University Press, 2014.

34.

Giurgiu

. Structural health monitoring with piezoelectric wafer active sensors. 2nd ed. Amsterdam: Elsevier, 2014.