Localization and Imaging of Micro-Cracks Using Nonlinear Lamb Waves with Imperfect Group-Velocity Matching

Nonlinear Lamb waves have attracted increasing attention for detecting and identifying microstructural changes in structural health monitoring. However, most identification methods that determine the damage locations based on the intersections of the elliptical loci will inevitably cause positioning errors due to the change of the group velocity before and after interaction with the damage. In this work, a method focusing on elliptical rings was proposed for localization and imaging of micro-cracks in a three-dimensional structure using nonlinear Lamb waves with imperfect group-velocity matching. The width of the elliptical rings can be determined by the degree of the group-velocity mismatching of nonlinear S0 modes. The mode pair S0-s0, satisfying approximate group-velocity matching, is mainly introduced by interacting with the micro-crack. The effectiveness of the proposed methodology for damage localization is verified by the experimental testing and numerical simulation. Although the length of the being-tested small crack (about 1 mm) is smaller than the wavelength of the incident fundamental Lamb wave (around 20 mm), it can be well identified and localized using nonlinear Lamb waves. The experimental results show that the proposed method enables more reliable localization of the small crack with the crossover areas, as compared with the intersections based on the ellipse method. Furthermore, a breathing crack not situated in the propagation path can also be well localized by the proposed method in comparison with those by the probability-based diagnostic imaging in the simulation cases.


Introduction
To reduce the risk of catastrophic failure and prolong the lifespan of structures, structural health monitoring (SHM) sensitive to defect-related changes in a real-time manner has been widely applied in aerospace, transportation, civil and mechanical structures, etc. Among various SHM techniques, ultrasonic Lamb wave has been considered as a promising way to detect damages in plate-like structures due to its appealing advantages of high sensitivity to different types of defects and its large inspection range with low attenuation [1]. The damage localization of linear cracks and through holes was achieved by the ellipse method, based on the time-of-flight of damage-scattered waves with a spatially distributed array of sensors attached to a homogeneous plate structure [2,3]. Moreover, efforts have been dedicated to characterizing the features of damage regions, such as the size and shape through various characterization algorithms [4,5]. Qualitative or quantitative detection of multiple damage cases [6,7] in engineering structures has become an increasing concern, although it seems to be a challenging task.
In general, the above-mentioned methods based on the linear features of Lamb waves (e.g., velocity and attenuation) are commonly sensitive to gross defects or open-type ones (e.g., circle holes, notches and cracks), and insensitive to small variations in microstructures (e.g., precipitates, inclusions and micro-voids) or contact-type damages (e.g., disbands, delamination and micro-cracks). If we want to improve the sensitivity of the linear Lamb wave to the small defects, an effective solution is to increase the excitation frequency, which will bring an additional complexity to signal processing because of the multimodal and dispersive properties of Lamb waves at higher frequencies. Recently, an alternative way to overcome the shortcoming is to employ the nonlinear feature of Lamb waves [8]. There are basically two types of acoustic nonlinearities, classical and non-classical, introduced by the coupling of nonlinear ultrasound and damages. The classical nonlinearity is mostly associated with intrinsic nonlinearity of materials such as anharmonicity and imperfection of atomic lattices, whereas the non-classical nonlinearity mainly arises from contact-type damages in structures through various mechanisms such as stress-strain hysteresis, contact nonlinearity, rough surfaces contact, Luxemburg-Gorky effect, etc. [9].
As a typical nonlinear Lamb wave, higher harmonic generation has a close tie with the behavior of the material nonlinearity. In the past two decades, the cumulative second harmonics of Lamb waves have attracted considerable attention because of their high sensitivity to material nonlinearities. Phase-velocity matching and non-zero power flux are regarded as the necessary conditions for cumulative second harmonic of Lamb waves generation [10,11]. Mode pairs that satisfy these conditions (e.g., S1-s2, A2-s4, and S2-s4) have been used for nonlinear ultrasonic detection of plasticity [12], fatigue [13], creep [14], pitting [15], and microstructural damage [16].
Recently, mode pair S0-s0 in a low-frequency range satisfying approximate phasevelocity matching was proposed to quantitatively evaluate the incipient stage of material degradation [17][18][19]. In addition, mode pair S0-s0 was also used to detect closed defects or imperfect bonding interfaces that open and close as the waves pass across them. Mori et al. [20] explored nonlinear interaction of S0 mode Lamb wave with an imperfect joint of plates by perturbation analysis and numerical simulation. Zhao et al. [21] performed numerical simulations on the propagation of S0 Lamb waves in thin plates with randomly distributed micro-cracks to study the behavior of nonlinear Lamb waves. They found that the acoustic nonlinear parameter increased linearly with the micro-crack density and the size of the micro-crack zone. Wang et al. [22] investigated the scattering patterns of interaction between S0 Lamb waves and fatigue cracks from analytical, numerical, and experimental perspectives. They found that an embryonic fatigue crack in an aluminum plate waveguide can be oriented accurately. Yang et al. [23] investigated the influence of the length of the fatigue crack as well as the incident wave angle on the amplitude of second harmonic generated by S0 mode. Although researchers have studied the interaction between nonlinear Lamb waves and interface nonlinearity, few studies have been reported on localization and imaging of micro-damage using nonlinear Lamb waves [24][25][26][27].
In the application of Lamb wave for SHM, appropriate damage identification algorithms have usually been adopted to achieve the damage position by giving an interpretable and intuitive image in a three-dimensional structure. One of the most common damage imaging algorithms is based on the so-called ellipse method [28,29], which determines the damage localization by the intersection of the elliptical loci based on the time of flight (ToF) and assumes that the group velocity remains unchanged before and after interacting with the damage. This will inevitably cause positioning errors due to the change of the group velocity before and after interaction with the damage [30]. When using the mode pair S0-s0 with quasi group-velocity matching, other candidate positioning methods will be implemented for damage localization. A total focus method based on decomposition of the time reversal operator was experimentally investigated by Zhou et al. [24] to locate two fatigue cracks in an aluminum plate. A probabilistic algorithm based on probability scan matrix was also implemented by Li et al. [25] to localize a microcrack in a thin plate, but it was considered that the microcrack had the highest probability of occurrence on the propagation path with a much dense transducer configuration. Furthermore, a damage is generally not a single point, but with a certain dimension [31]. Su et al. [26] proposed a probability-based diagnostic imaging (PDI) based on elliptical loci to determine the most possible intersection areas for damage localization. It was demonstrated that the method of probability distributions can be used to predict or estimate damage area. Hong et al. [27] localized a fatigue crack starting at a rivet hole in an aluminum plate by a probabilistic imaging algorithm using nonlinear Lamb waves, in which the employed mode pair S1-s2 needed to exactly meet the phase-and group-velocity matching.
Therefore, the objective of this study is to propose a modified ellipse method for accurate localization and imaging of a micro-crack in a three-dimensional structure using nonlinear Lamb waves with imperfect group-velocity matching. This study is organized as follows. Section 2 introduces the concept of the modified method and establishes the mathematical formulation for damage localization using nonlinear S0 mode Lamb waves with quasi group-velocity matching. Section 3 provides the details of the experimental setup to estimate the method's performance for localization of a small-crack in an aluminum plate. The experimental results for localization of the small-crack based on linear and nonlinear features are discussed, and the identification results of the proposed method were compared with that of the traditional ellipse method. Section 4 presents finite element simulations of a breathing crack in comparison with the experimental verification. The proposed method is also validated when the breathing crack is not situated in the propagation path. Finally, conclusions are drawn in the Section 5.

Modified Ellipse Method
When Lamb waves encounter a defect, a propagating wave will be reflected and/or partly transmitted. These damage-scattered waves will carry enough information about the damage. For instance, time of flight (ToF) and amplitude can be used for damage identification. The ToF is defined as a time duration for a specific wave packet travelling a certain distance. Generally, ToF represents, to an extent, a linear correlation with the damage parameters such as its location, size and shape.
To illustrate the principle of damage localization using ToF, a sensing path T I − T J (I, J = 1, 2, . . . N, I = J) from part of a sensor network composed of N pieces of PZT transducers is schematically shown in Figure 1a. Using Cartesian coordinates, the actuator T I , the sensor T J and the center of the damage are supposed to be situated at (x I , y I ), (x J , y J ), (x m , y m ), respectively. The widely known ellipse method can be used to describe their relationship as given by where V g and V d−s are the group velocities of the incident wave and the damage-scattered wave, respectively. t I J is the ToF it takes for the incident wave to travel from the actuator to the damage (t g ), and then to the sensor after damage scattering (t d−s ). l Im is the distance between the actuator and the damage center, and l mJ is the distance between the damage center and the sensor. By repeating the above analysis for at least two different sensing paths in the sensor network, a nonlinear equation set with two unknown parameters (x m , y m ) is created. For each equation described by Equation (1) in the set, the solution, which mathematically depicts an elliptical (provided V g = V d−s ) locus, indicates all possible damage locations perceived by the sensing path T I − T J (Figure 1a). If the nonlinear equation set is solved, the damage locations can be presumably positioned at the intersections of most of the loci. However, the ellipse method will inevitably cause positioning errors when the group velocity changes (V g = V d−s , e.g., mode conversion) before and after interacting with the damage. Furthermore, a damage location should not be simply considered as a single point at the intersection with specific geometric parameters. Instead, it is more likely to be a damage zone. For example, a crack has a certain length, and there may be a plastic-deformation zone at the tip of the crack. Thus, a modified ellipse method focusing on elliptical rings is proposed to solve the problems in accurate damage localization. the damage. For instance, time of flight (ToF) and amplitude can be used for damage identification. The ToF is defined as a time duration for a specific wave packet travelling a certain distance. Generally, ToF represents, to an extent, a linear correlation with the damage parameters such as its location, size and shape.
To illustrate the principle of damage localization using ToF, a sensing path − ( , = 1,2, … N, I ≠ J) from part of a sensor network composed of N pieces of PZT transducers is schematically shown in Figure 1a. Using Cartesian coordinates, the actuator , the sensor and the center of the damage are supposed to be situated at ( , ), ( , ), ( , ), respectively. The widely known ellipse method can be used to describe their relationship as given by Given the velocity difference ∆v = V g − V d−s , the nonlinear equation for damage localization could be expressed as The damage locations are positioned at the intersections of most of the elliptical (when ∆v = 0) or ellipse-like (when ∆v = 0) loci. By neglecting nonzero ∆v, two different elliptical loci are derived from the left-and right-hand side of Equation (2), which are smaller and bigger than the ellipse-like locus, respectively. Moreover, no matter how small the difference between the two elliptical loci is, the ellipse-like locus indicating all the possible damage locations will situate at a region between the two elliptical loci. As shown in Figure 1b, two elliptical loci shown as dashed lines were plotted with the group velocities V g and V d−s in a sensing path, respectively, and the region between the two elliptical loci forms an elliptical ring, described by Then, the elliptical ring, not just one ellipse-like locus, should be defined as the potential region for damage locations. If two or more sensing paths are available, the damage locations can be determined by the crossover area at the intersections of most of the elliptical rings.
To present the identification results of damage localization, a grey-scale image was adopted to describe a two-dimensional structure with exactly the same dimension. The inspection area of the structure was virtually and evenly meshed. At each mesh node, the pixel value D I J (x m , y m ) was quantified with respect to the region between the two elliptical loci. If the mesh node fell within the region, it had the high pixel value of 1 for damage occurrence. If not, it had the low pixel value of 0. By integrating information from all sensing paths in the sensor network, the ultimate identification image could be produced using an image fusion scheme as According to the above analysis, this modified ellipse method makes use of a region between two elliptical loci, instead of an elliptical or ellipse-like locus, for damage localization. Generally, it can improve the accuracy of damage localization, while considering the change in group velocity before and after interacting with the damage, and it has the potential to evaluate quantitatively the size and even the level of structural damage because it makes use of the region between two elliptical loci with respect to the certain dimension of the damage location. Moreover, it is suitable not only for locating micro-damages using nonlinear Lamb waves with quasi group-velocity matching (e.g., mode pair S0-s0) but also for damage localization with mode conversion at the damage.

Experimental Verification
The modified ellipse method of nonlinear Lamb waves was then corroborated by ultrasonic experimental testing of a prefabricated small crack in an aluminum plate.

Experimental Setup
The experimental testing system, as shown in Figure 2, consisted of a dual channel arbitrary/function generator (Textronix AFG 3022C), a power amplifier (RITEC GA 2500), a mixed domain oscilloscope (Textronix MDO 3012), and a personal computer (PC) connected with the oscilloscope. The specimen was a 6061-T6 aluminum plate with a dimension of 450 mm × 400 mm × 2 mm. Four PZT piezoelectric transducers, 10 mm in diameter and 1 mm in thickness each, were permanently bonded to the surface of the aluminum plate by means of a thin cyanoacrylate adhesive layer. All the PZT transducers used for exciting or sensing Lamb wave signals had the same polarization direction.
According to the above analysis, this modified ellipse method makes use of a region between two elliptical loci, instead of an elliptical or ellipse-like locus, for damage localization. Generally, it can improve the accuracy of damage localization, while considering the change in group velocity before and after interacting with the damage, and it has the potential to evaluate quantitatively the size and even the level of structural damage because it makes use of the region between two elliptical loci with respect to the certain dimension of the damage location. Moreover, it is suitable not only for locating microdamages using nonlinear Lamb waves with quasi group-velocity matching (e.g., mode pair S0-s0) but also for damage localization with mode conversion at the damage.

Experimental Verification
The modified ellipse method of nonlinear Lamb waves was then corroborated by ultrasonic experimental testing of a prefabricated small crack in an aluminum plate.

Experimental Setup
The experimental testing system, as shown in Figure 2, consisted of a dual channel arbitrary/function generator (Textronix AFG 3022C), a power amplifier (RITEC GA 2500), a mixed domain oscilloscope (Textronix MDO 3012), and a personal computer (PC) connected with the oscilloscope. The specimen was a 6061-T6 aluminum plate with a dimension of 450 mm 400 mm 2 mm. Four PZT piezoelectric transducers, 10 mm in diameter and 1 mm in thickness each, were permanently bonded to the surface of the aluminum plate by means of a thin cyanoacrylate adhesive layer. All the PZT transducers used for exciting or sensing Lamb wave signals had the same polarization direction. A 10-cycle Hanning-windowed sine wave, generated by the arbitrary generator, was applied as the excitation tone burst. The burst was amplified by means of the power amplifier before it was sent to one of the PZT transducers. The rest of the PZT transducers were used to receive Lamb wave signals. The signals were captured by the oscilloscope at a sample rate of 50 MHz with 256-time averaging, and then transmitted to the computer. A 10-cycle Hanning-windowed sine wave, generated by the arbitrary generator, was applied as the excitation tone burst. The burst was amplified by means of the power amplifier before it was sent to one of the PZT transducers. The rest of the PZT transducers were used to receive Lamb wave signals. The signals were captured by the oscilloscope at a sample rate of 50 MHz with 256-time averaging, and then transmitted to the computer. The excitation-acquisition procedure was in turn implemented on each PZT transducer, such that six received signals were collected from these sensing paths. In order to reduce the influence of stochastic error of the experimental measurements, the received signals were averaged with three measurements on each sensing path.
A global cartesian coordinate was introduced, the bottom-left corner of a rectangle with an in-plane dimension of 200 mm × 100 mm was at the origin, and the abscissa axis was parallel to the lower edge of the plate. Based on this reference coordinate system, the PZT transducers labeled as A, B, C, and D were bonded at coordinates (25,10), (25,90), (175, 90), (175, 10), respectively.
Lamb wave signals in the aluminum plate were first received in an intact (crack-free) state. Then, a small crack at coordinates (100, 50) was produced by a milling cutter at the center of the aluminum plate. The physical length of the crack was about 1.025 mm, and the width was about 0.425 mm, as observed by means of an optical microscope ( Figure 2). Afterward, Lamb wave signals in the aluminum plate with the small crack were measured in the current state.

Mode-Tuning of LAMB Wave
Although Lamb wave signals are complicated because of their multimodal characteristics, two fundamental Lamb modes S0 and A0 can exist and propagate in low frequency range with weak dispersion. S0 mode Lamb wave was selected as an incident wave because this mode is more sensitive to the through-thickness crack than the A0 mode [1].
A mode-tuning experiment [32] was conducted to determine the optimal excitation frequency of the S0 mode. A pair of PZT transducers with a distance of 170 mm was attached to the upper surface of the plate. The tone burst was applied to one of the PZT transducers to generate Lamb wave, and the other was used to receive Lamb wave signals. The excitation frequency of the tone burst was swept from 100 to 500 kHz in steps of 10 kHz. At each excitation frequency, the amplitudes of the S0 and A0 modes were recorded. Figure 3 shows the mode-tuning result. It can be seen that the amplitudes of S0 and A0 modes fluctuate periodically versus the excitation frequency. The amplitude of the S0 mode reached its maximum at an excitation frequency of 280 kHz, while the amplitude of the A0 mode was relatively smaller. Therefore, to obtain a sufficiently strong nonlinear effect and to restrain the A0 mode, the S0-dominated Lamb wave was excited at a frequency of 280 kHz in the experiment. At this frequency, the length of the small crack was much smaller than the wavelength of the incident wave (about 20 mm). The small crack could thus be considered as micro-damage [24].
The excitation-acquisition procedure was in turn implemented on each PZT transducer such that six received signals were collected from these sensing paths. In order to reduce the influence of stochastic error of the experimental measurements, the received signals were averaged with three measurements on each sensing path.
A global cartesian coordinate was introduced, the bottom-left corner of a rectangle with an in-plane dimension of 200 mm × 100 mm was at the origin, and the abscissa axis was parallel to the lower edge of the plate. Based on this reference coordinate system, the PZT transducers labeled as A, B, C, and D were bonded at coordinates (25,10), (25,90) (175, 90), (175, 10), respectively.
Lamb wave signals in the aluminum plate were first received in an intact (crack-free) state. Then, a small crack at coordinates (100, 50) was produced by a milling cutter at the center of the aluminum plate. The physical length of the crack was about 1.025 mm, and the width was about 0.425 mm, as observed by means of an optical microscope ( Figure 2) Afterward, Lamb wave signals in the aluminum plate with the small crack were measured in the current state.

Mode-Tuning of LAMB Wave
Although Lamb wave signals are complicated because of their multimodal characteristics, two fundamental Lamb modes S0 and A0 can exist and propagate in low frequency range with weak dispersion. S0 mode Lamb wave was selected as an incident wave because this mode is more sensitive to the through-thickness crack than the A0 mode [1].
A mode-tuning experiment [32] was conducted to determine the optimal excitation frequency of the S0 mode. A pair of PZT transducers with a distance of 170 mm was attached to the upper surface of the plate. The tone burst was applied to one of the PZT transducers to generate Lamb wave, and the other was used to receive Lamb wave signals The excitation frequency of the tone burst was swept from 100 to 500 kHz in steps of 10 kHz. At each excitation frequency, the amplitudes of the S0 and A0 modes were recorded Figure 3 shows the mode-tuning result. It can be seen that the amplitudes of S0 and A0 modes fluctuate periodically versus the excitation frequency. The amplitude of the S0 mode reached its maximum at an excitation frequency of 280 kHz, while the amplitude of the A0 mode was relatively smaller. Therefore, to obtain a sufficiently strong nonlinear effect and to restrain the A0 mode, the S0-dominated Lamb wave was excited at a frequency of 280 kHz in the experiment. At this frequency, the length of the small crack was much smaller than the wavelength of the incident wave (about 20 mm). The small crack could thus be considered as micro-damage [24].

Comparison of Group Velocities
Before using the ellipse methods, the group velocity of S0 mode was determined. A pair of PZT transducers apart from the same actuator with different distances were used to receive Lamb wave signals in the intact plate. The group velocity was calculated based on the respective distance and arrival time at a frequency. The experimental testing was conducted over a range of frequencies. Following the experimental testing, the group velocity dispersion curves based on dispersion equation were also calculated for Lamb waves in the aluminum plate from 0 to 1 MHz·mm.
As shown in Figure 4, the results of the experimental testing and dispersion curves do not exactly agree with each other, which may result from the inconsistency between the theoretical and actual material parameters. Nevertheless, the results still show that there is effective actuation and reception of Lamb wave signals in the aluminum plate. It was found that the group velocity of S0 mode at a frequency-thickness of 0.56 MHz·mm was about 4978.6 m/s. Meanwhile, the group velocity at a frequency-thickness of 1.12 MHz·mm was about 4688.6 m/s as determined by extrapolation. Figure 3. Lamb wave mode-tuning curve for a 2 mm-thick aluminum plate excited by a 1 mm-thick PZT transducer.

Comparison of Group Velocities
Before using the ellipse methods, the group velocity of S0 mode was determined. A pair of PZT transducers apart from the same actuator with different distances were used to receive Lamb wave signals in the intact plate. The group velocity was calculated based on the respective distance and arrival time at a frequency. The experimental testing was conducted over a range of frequencies. Following the experimental testing, the group velocity dispersion curves based on dispersion equation were also calculated for Lamb waves in the aluminum plate from 0 to 1 MHz.mm.
As shown in Figure 4, the results of the experimental testing and dispersion curves do not exactly agree with each other, which may result from the inconsistency between the theoretical and actual material parameters. Nevertheless, the results still show that there is effective actuation and reception of Lamb wave signals in the aluminum plate. It was found that the group velocity of S0 mode at a frequency-thickness of 0.56 MHz•mm was about 4978.6 m/s. Meanwhile, the group velocity at a frequency-thickness of 1.12 MHz•mm was about 4688.6 m/s as determined by extrapolation.

Extraction of Arrival Time
After the group velocity measurements, the arrival time of the damage-scattered signals was also necessitated for the ellipse methods. The received signals were first processed through a low-pass filter with a frequency of 1 MHz, which eliminated the highfrequency noise from the environment. Figure 5 shows typical waveforms of the signals received from the sensing path AB (A is an actuator and B is a sensor) in the intact and the current (with a small crack) states, respectively. The S0-dominated Lamb wave signal was more obvious, while the damage-scattered signals submerged in multimode and boundary reflection signals could not be directly recognized. Then, signal processing methods were applied to obtain the linear and nonlinear features of the damage-scattered signals.

Extraction of Arrival Time
After the group velocity measurements, the arrival time of the damage-scattered signals was also necessitated for the ellipse methods. The received signals were first processed through a low-pass filter with a frequency of 1 MHz, which eliminated the highfrequency noise from the environment. Figure 5 shows typical waveforms of the signals received from the sensing path AB (A is an actuator and B is a sensor) in the intact and the current (with a small crack) states, respectively. The S0-dominated Lamb wave signal was more obvious, while the damage-scattered signals submerged in multimode and boundary reflection signals could not be directly recognized. Then, signal processing methods were applied to obtain the linear and nonlinear features of the damage-scattered signals.
For processing the linear features, the received signals in the two states were normalized according to the respective maximum value, and then subtraction of the two received signals was performed to obtain the damage-scattered signal, as shown in Figure 6. The damage-scattered signal was further subjected to Hilbert transform to acquire an envelope curve, and the time corresponding to the maximum peak (red dot) of the envelope curve was taken as the arrival time of the damage-scattered signal. For processing the linear features, the received signals in the two states were normalized according to the respective maximum value, and then subtraction of the two received signals was performed to obtain the damage-scattered signal, as shown in Figure 6. The damage-scattered signal was further subjected to Hilbert transform to acquire an envelope curve, and the time corresponding to the maximum peak (red dot) of the envelope curve was taken as the arrival time of the damage-scattered signal. For processing the nonlinear features, short time Fourier transform (STFT) was applied to transform the received signals from the time domain into the time-frequency domain. Spectrogram of a received signal from the sensing path AB in the current state is shown in Figure 7, including amplitude profiles at the respective fundamental and double frequencies for both the intact and the current states. It was found that the amplitude profile at the fundamental frequency in the current state (red dashed line) was essentially the same with that in the intact state (blue line), while the amplitude profile at the double frequency in the current state (black dashed line) was higher and more obvious than that in the intact state. Due to interaction with the small crack, the amplitude profile at the  For processing the linear features, the received signals in the two states were normalized according to the respective maximum value, and then subtraction of the two received signals was performed to obtain the damage-scattered signal, as shown in Figure 6. The damage-scattered signal was further subjected to Hilbert transform to acquire an envelope curve, and the time corresponding to the maximum peak (red dot) of the envelope curve was taken as the arrival time of the damage-scattered signal. For processing the nonlinear features, short time Fourier transform (STFT) was applied to transform the received signals from the time domain into the time-frequency domain. Spectrogram of a received signal from the sensing path AB in the current state is shown in Figure 7, including amplitude profiles at the respective fundamental and double frequencies for both the intact and the current states. It was found that the amplitude profile at the fundamental frequency in the current state (red dashed line) was essentially the same with that in the intact state (blue line), while the amplitude profile at the double frequency in the current state (black dashed line) was higher and more obvious than that in the intact state. Due to interaction with the small crack, the amplitude profile at the For processing the nonlinear features, short time Fourier transform (STFT) was applied to transform the received signals from the time domain into the time-frequency domain. Spectrogram of a received signal from the sensing path AB in the current state is shown in Figure 7, including amplitude profiles at the respective fundamental and double frequencies for both the intact and the current states. It was found that the amplitude profile at the fundamental frequency in the current state (red dashed line) was essentially the same with that in the intact state (blue line), while the amplitude profile at the double frequency in the current state (black dashed line) was higher and more obvious than that in the intact state. Due to interaction with the small crack, the amplitude profile at the double frequency had an extra hump in the current state, considered as the second harmonics of S0 mode. Thus, the time corresponding to the maximum peak (black dot) of the amplitude profile was regarded as the arrival time of the crack-induced S0 mode. The arrival time was also related to the length of the window function in the STFT as well as the absolute time between the excitation and the reception signals.
double frequency had an extra hump in the current state, considered as the second harmonics of S0 mode. Thus, the time corresponding to the maximum peak (black dot) of the amplitude profile was regarded as the arrival time of the crack-induced S0 mode. The arrival time was also related to the length of the window function in the STFT as well as the absolute time between the excitation and the reception signals.

Localization and Imaging
The parameters of the group velocities and the arrival times as detailed above were used to realize localization and imaging of the small crack. By combining individual source images obtained from all sensing paths in the network, the final identification image was produced as shown in Figure 8. Based on the ellipse method, the positioning results of the small crack based on linear and nonlinear features are shown in Figure 8a When the damage is situated in the propagation path or close to the sensor position, the corresponding elliptical loci cannot be constructed [2]. Thus, the signals received from such propagation paths, namely, the sensing paths AC and BD, should be excluded for accurate damage localization. The re-positioning results based on linear and nonlinear features are shown in Figure 8b,d, respectively. It can be seen that the elliptical loci based on nonlinear features fall accurately on the location of the small crack, while the elliptical loci based on linear features still do not.

Localization and Imaging
The parameters of the group velocities and the arrival times as detailed above were used to realize localization and imaging of the small crack. By combining individual source images obtained from all sensing paths in the network, the final identification image was produced as shown in Figure 8. Based on the ellipse method, the positioning results of the small crack based on linear and nonlinear features are shown in Figure 8a When the damage is situated in the propagation path or close to the sensor position, the corresponding elliptical loci cannot be constructed [2]. Thus, the signals received from such propagation paths, namely, the sensing paths AC and BD, should be excluded for accurate damage localization. The re-positioning results based on linear and nonlinear features are shown in Figure 8b,d, respectively. It can be seen that the elliptical loci based on nonlinear features fall accurately on the location of the small crack, while the elliptical loci based on linear features still do not.
By seeking the intersections of most of these elliptical loci, the small crack position can be localized and highlighted. The center of the small crack is at coordinates (100, 50), and there are seven intersections of most of the elliptical loci based on linear features, all far away from the coordinates (100, 50). However, there are only two intersections of most of elliptical loci based on nonlinear features at coordinates (95.5, 49.5) and (102, 49.5), which are close to the coordinates (100, 50). It is concluded here that the small crack cannot be identified based on linear feature (Figure 8b), and such a failure can be attributed to the minute dimension of the small crack, which causes not much wave scattering (the wavelength of the incident wave is around 20 mm at the selected excitation frequency, greater than the length of the small crack). While the small crack can be identified and localized based on nonlinear feature (Figure 8d). These elliptical loci are thickened for clarity. The actual location of the small crack is highlighted by a red line, while the blue plus signs indicate the intersections. By seeking the intersections of most of these elliptical loci, the small crack position can be localized and highlighted. The center of the small crack is at coordinates (100, 50), and there are seven intersections of most of the elliptical loci based on linear features, all far away from the coordinates (100, 50). However, there are only two intersections of most of elliptical loci based on nonlinear features at coordinates (95.5, 49.5) and (102, 49.5), which are close to the coordinates (100, 50). It is concluded here that the small crack cannot be identified based on linear feature (Figure 8b), and such a failure can be attributed to the minute dimension of the small crack, which causes not much wave scattering (the wavelength of the incident wave is around 20 mm at the selected excitation frequency, greater than the length of the small crack). While the small crack can be identified and This positioning result seems to be at the location of the small crack, but there is still a certain error from the actual location of the small crack. This error is most likely due to the change in group velocity before and after interaction with the small crack. In order to achieve more accurate damage localization, the modified ellipse method based on nonlinear feature is employed to localize the small crack. The identification results with and without the signals received from the sensing paths AC and BD are shown in Figure 8e,f, respectively. As shown in Figure 8e, the modified ellipse method fails in damage localization when the small crack is situated in the propagation path. While in Figure 8f the crossover area of damage localization determined by the proposed method accurately appears on the actual location of the small crack. This crossover area is enclosed by two vertical blue dashed lines with x = 96.8 mm and 100.8 mm and two horizontal blue dashed lines with y = 47.5 mm and 51.5 mm, respectively. The darker the crossover area is, the more likely it is the location of the damage. The determined crossover area involves the small crack at coordinates (100, 50) and includes but is not limited to the intersections at coordinates (95.5, 49.5) and (102, 49.5). The distance between the two horizontal lines within the crossover area can be used to approximately estimate the length of the small crack (about 1 mm) in the y-axis direction. Thus, the modified ellipse method utilizing the imperfect group-velocity matching of nonlinear S0 modes can achieve the small crack localization, and the result for damage localization becomes more accurate and more convincing compared with that of the ellipse method. It can be further seen from the identification results that the thus-determined location area can also reveal the size and shape, and the damage level of the small crack to a certain extent.

Simulation Validation
In addition to the experimental verification, a three-dimensional finite element (FE) model was developed and realized in conjunction with a commercial FEM software (ABAQUS/CAE) to localize a breathing crack.

Simulation Setup
The dimensions of the aluminum plate employed for the simulation were 200 mm × 100 mm × 1 mm. In order to have better nonlinear effect, a 3 mm-long crack was set at the center of the plate with an initial clearance of zero between the two interfaces, as shown in Figure 9. Additionally, in order to realize very minimal or no wave reflection from the boundary, a gradually damping artificial boundary was introduced to simulate a non-reflecting boundary condition [33]. The boundary was divided into 25 layers, and the damping coefficients of the individual layers were gradually increased from the inner layer to the finite boundary. Both the boundary length and the damping coefficient of each layer are dependent on the excitation frequency. The intact state without crack and the current state with the crack were studied, separately.  The Young's modulus, Poisson's ratio and density of the aluminum plate were 69 GPa, 0.33 and 2700 kg/m 3 , respectively. The non-reflecting boundary had the same material properties with the aluminum plate except having the damping coefficients. Moreover, the crack was modeled by embedding a seam crack at the center of the plate. Hence, it was initially closed, but it could be opened when Lamb waves traveled through. Hard normal contact and frictionless tangential contact were [22] applied to the interfaces of the seam crack to prevent nodes penetration. Thus, a breathing behavior appeared when Lamb waves interacted with the crack.
Hanning-window modulated five-cycle sinusoidal tone bursts at a central frequency of 250 kHz were excited by applying a pair of symmetric displacement on two points, which were symmetrically positioned on the upper and lower surfaces of the plate. Previous work [34] has described in detail this type of point displacement modeling to gen- The Young's modulus, Poisson's ratio and density of the aluminum plate were 69 GPa, 0.33 and 2700 kg/m 3 , respectively. The non-reflecting boundary had the same material properties with the aluminum plate except having the damping coefficients. Moreover, the crack was modeled by embedding a seam crack at the center of the plate. Hence, it was initially closed, but it could be opened when Lamb waves traveled through. Hard normal contact and frictionless tangential contact were [22] applied to the interfaces of the seam crack to prevent nodes penetration. Thus, a breathing behavior appeared when Lamb waves interacted with the crack.
Hanning-window modulated five-cycle sinusoidal tone bursts at a central frequency of 250 kHz were excited by applying a pair of symmetric displacement on two points, which were symmetrically positioned on the upper and lower surfaces of the plate. Previous work [34] has described in detail this type of point displacement modeling to generate single S0 mode Lamb wave. The displacement U3 at one point on the upper surface was monitored as the received signal. A global cartesian coordinate was also introduced, and the bottom-left corner of the plate was at the origin. The points for exciting and receiving signals were labeled as A, B, C, and D at coordinates (25,10), (25,90), (175, 90), (175, 10), respectively.
In order to obtain sufficient accuracy and high efficiency, the dimension of the elements was around 0.5 mm × 0.5 mm × 0.25 mm, such that there were at least 20 elements per wavelength of the Lamb waves. This spatial resolution was chosen to warrant simulation precision and error convergence. Smaller elements were applied in the vicinity of the crack because of complicated mechanical response at the crack zone, and the elements were meshed automatically. Eight-noded 3D reduced integrated linear brick elements, C3D8R, were used to model the aluminum plate and the boundary. A time step of 2 ns was used to ensure the accuracy of the second harmonic generation [35]. Figure 10 shows the process of the Lamb wave traversing the breathing crack. It can be clearly seen that the stress induced by Lamb wave propagation is continuously transmitted through the breathing crack when it is closed (Figure 10a), while on the contrary, the transmission process is interrupted when the crack is open (Figure 10b). This phenomenon is defined as so-called contact acoustic nonlinearity (CAN) [9]. From the results of the stress transmission, the breathing behavior of the crack was observed; meanwhile, higher harmonics were generated. As a typical result in the intact and the current states, Figure 11a shows the excitation signal and the Lamb wave signal received from the sensing path AB with a distance of 80 mm. In addition to the received signal at fundamental frequency, the received signal at double frequency obtained by means of a band pass filters with cutoff frequencies of 450 and 550 kHz is also shown in Figure 11b. It was found that there was no noticeable difference between the two states at the fundamental frequency, which allowed no acquisition of the crack-induced signal. However, compared to the intact state, a wave packet at the double frequency was distinct in the current state, which represented nonlinear response generated by interacting with the breathing crack. From the results of the stress transmission, the breathing behavior of the crack was observed; meanwhile, higher harmonics were generated. As a typical result in the intact and the current states, Figure 11a shows the excitation signal and the Lamb wave signal received from the sensing path AB with a distance of 80 mm. In addition to the received signal at fundamental frequency, the received signal at double frequency obtained by means of a band pass filters with cutoff frequencies of 450 and 550 kHz is also shown in Figure 11b. It was found that there was no noticeable difference between the two states at the fundamental frequency, which allowed no acquisition of the crack-induced signal. However, compared to the intact state, a wave packet at the double frequency was distinct in the current state, which represented nonlinear response generated by interacting with the breathing crack. means of a band pass filters with cutoff frequencies of 450 and 550 kHz is also shown in Figure 11b. It was found that there was no noticeable difference between the two states at the fundamental frequency, which allowed no acquisition of the crack-induced signal. However, compared to the intact state, a wave packet at the double frequency was distinct in the current state, which represented nonlinear response generated by interacting with the breathing crack.

Signal Processing
The Lamb wave signals were received from all sensing paths in the intact and the current states, separately, using the excitation-acquisition procedure as described in Section 3.1. The group velocity was determined first by the process depicted in Section 3.3. The group velocities at the fundamental and double frequency were 5348.4 and 5319.1 m/s, respectively. The arrival time was collected sequentially, and the process of extracting linear and nonlinear features was essentially the same with that in Section 3.4.

Signal Processing
The Lamb wave signals were received from all sensing paths in the intact and the current states, separately, using the excitation-acquisition procedure as described in Section 3.1. The group velocity was determined first by the process depicted in Section 3.3. The group velocities at the fundamental and double frequency were 5348.4 and 5319.1 m/s, respectively. The arrival time was collected sequentially, and the process of extracting linear and nonlinear features was essentially the same with that in Section 3.4.

Localization and Imaging
With the group velocity and arrival time from signal processing, the breathing crack in the aluminum plate was localized and reconstructed by the ellipse methods.
As shown in Figure 12a,  Figure 12f, the intersection of most of the elliptical rings can be localized and highlighted. The damage location defined as the crossover area accurately situates at the actual location of the breathing crack at coordinates (100, 50). Moreover, the distance between the two horizontal lines within the crossover area can accurately estimate the length of the breathing crack (about 3 mm) in the vertical direction. It can be seen that when the length of the breathing crack is smaller than the wavelength of the incident wave at the selected excitation frequency, linear Lamb waves are not sensitive to the breathing crack and have limitations in identifying and locating micro-damages. In contrast, nonlinear Lamb waves can identify the breathing crack. Moreover, the breathing crack can be localized more precisely by the proposed method based on nonlinear S0 modes. This result is in good agreement with the experimental testing, thus verifying the effectiveness of the proposed method.
that when the length of the breathing crack is smaller than the wavelength of the incident wave at the selected excitation frequency, linear Lamb waves are not sensitive to the breathing crack and have limitations in identifying and locating micro-damages. In contrast, nonlinear Lamb waves can identify the breathing crack. Moreover, the breathing crack can be localized more precisely by the proposed method based on nonlinear S0 modes. This result is in good agreement with the experimental testing, thus verifying the effectiveness of the proposed method. In the simulation, the ellipse method based on nonlinear features for damage localization is almost as accurate as the modified ellipse method. This is because the group velocity at fundamental frequency is not significantly different from that at double frequency, namely, the group velocity remains essentially unchanged before and after interacting with the breathing crack. Thus, both imaging methods can be applied to identify and localize the breathing crack accurately.
In order to further confirm the performance of the modified ellipse method, the breathing crack set at a random position (118.75, 70), not situated in any propagation paths, was identified and localized using the proposed method. With the signals received from all six sensing paths, the identification result of the breathing crack based on nonlinear features is shown in Figure 13a. The intersection of most of the elliptical rings is enclosed by two vertical lines with x = 117.7 mm and 119.1 mm and two horizontal lines with y = 68.3 mm and 72.3 mm (blue dashed lines in Figure 13a). It is found that the crossover area precisely appears on the actual location of the breathing crack at coordinates (118.75, 70), and the distance between the two horizontal lines within the crossover area accurately estimates the length of the breathing crack (about 3 mm) in the vertical direction. In the simulation, the ellipse method based on nonlinear features for damage localization is almost as accurate as the modified ellipse method. This is because the group velocity at fundamental frequency is not significantly different from that at double frequency, namely, the group velocity remains essentially unchanged before and after interacting with the breathing crack. Thus, both imaging methods can be applied to identify and localize the breathing crack accurately.
In order to further confirm the performance of the modified ellipse method, the breathing crack set at a random position (118.75, 70), not situated in any propagation paths, was identified and localized using the proposed method. With the signals received from all six sensing paths, the identification result of the breathing crack based on nonlinear features is shown in Figure 13a. The intersection of most of the elliptical rings is enclosed by two vertical lines with x = 117.7 mm and 119.1 mm and two horizontal lines with y = 68.3 mm and 72.3 mm (blue dashed lines in Figure 13a). It is found that the crossover area precisely appears on the actual location of the breathing crack at coordinates (118.75, 70), and the distance between the two horizontal lines within the crossover area accurately estimates the length of the breathing crack (about 3 mm) in the vertical direction. In addition, a probability-based diagnostic imaging (PDI) based on ToF for damage localization [27] was also analyzed and discussed. In this imaging method based on the elliptical loci, probability distribution value was assigned to the virtualized grid nodes in the structure, and the greater the probability value of this location, the more likely the damage location. Under the same situation as Figure 13a, the identification result using PDI is shown in Figure 13b. As can be seen from Figure 13a,b, both the crossover area determined by the proposed method and the highlighted damage region defined by the PDI with a threshold fall on the actual location of the breathing crack with a certain dimension. It was found that these two imaging methods based on nonlinear S0 modes can achieve localization and imaging of the breathing crack more precisely. Thus, the modified ellipse method using nonlinear S0 mode Lamb waves may be a suitable candidate for identifying and locating micro-damages.

Conclusions
A modified ellipse method using nonlinear S0 mode Lamb waves with approximate group-velocity matching is proposed for localization and imaging of micro-cracks. Compared with the traditional ellipse method based on the intersections of the elliptical loci, the modified ellipse method utilizes the elliptical rings between two elliptical loci to localize damages accurately.
In the experimental testing, a rectangle-like in-situ array of PZT piezoelectric transducers was used to excite and receive Lamb wave signals in an aluminum plate. When the length of a small crack (about 1 mm) is smaller than the wavelength of the incident wave (around 20 mm), linear Lamb waves are insensitive to the small crack, while the small crack can be identified and localized more precisely by the proposed method with nonlinear Lamb waves. For the small crack centered in the aluminum plate, the crossover area determined by the proposed method includes, but is not limited to, the intersections based on the ellipse method. Moreover, the dimension of the crossover area can be used to estimate the length of the small crack in the y-axis direction, while an intersection point is defined by the ellipse method. Extensive FE simulations were conducted to localize a breathing crack in the aluminum plate. The numerical results are in good agreement with the experimental findings. In addition, the breathing crack not situated in any of the propagation paths can also be localized by the proposed method as well as by the PDI.
In total, the results presented here demonstrate that the modified ellipse method based on nonlinear S0 mode Lamb waves can accurately realize the localization and imaging of micro-cracks, and it has the potential to evaluate quantitatively the size of microcracks and even to monitor crack growth. In addition, a probability-based diagnostic imaging (PDI) based on ToF for damage localization [27] was also analyzed and discussed. In this imaging method based on the elliptical loci, probability distribution value was assigned to the virtualized grid nodes in the structure, and the greater the probability value of this location, the more likely the damage location. Under the same situation as Figure 13a, the identification result using PDI is shown in Figure 13b. As can be seen from Figure 13a,b, both the crossover area determined by the proposed method and the highlighted damage region defined by the PDI with a threshold fall on the actual location of the breathing crack with a certain dimension. It was found that these two imaging methods based on nonlinear S0 modes can achieve localization and imaging of the breathing crack more precisely. Thus, the modified ellipse method using nonlinear S0 mode Lamb waves may be a suitable candidate for identifying and locating micro-damages.

Conclusions
A modified ellipse method using nonlinear S0 mode Lamb waves with approximate group-velocity matching is proposed for localization and imaging of micro-cracks. Compared with the traditional ellipse method based on the intersections of the elliptical loci, the modified ellipse method utilizes the elliptical rings between two elliptical loci to localize damages accurately.
In the experimental testing, a rectangle-like in-situ array of PZT piezoelectric transducers was used to excite and receive Lamb wave signals in an aluminum plate. When the length of a small crack (about 1 mm) is smaller than the wavelength of the incident wave (around 20 mm), linear Lamb waves are insensitive to the small crack, while the small crack can be identified and localized more precisely by the proposed method with nonlinear Lamb waves. For the small crack centered in the aluminum plate, the crossover area determined by the proposed method includes, but is not limited to, the intersections based on the ellipse method. Moreover, the dimension of the crossover area can be used to estimate the length of the small crack in the y-axis direction, while an intersection point is defined by the ellipse method. Extensive FE simulations were conducted to localize a breathing crack in the aluminum plate. The numerical results are in good agreement with the experimental findings. In addition, the breathing crack not situated in any of the propagation paths can also be localized by the proposed method as well as by the PDI.
In total, the results presented here demonstrate that the modified ellipse method based on nonlinear S0 mode Lamb waves can accurately realize the localization and imaging of micro-cracks, and it has the potential to evaluate quantitatively the size of micro-cracks and even to monitor crack growth.