A Lamb-wave based SHM for multi-damage localizations in large composite plates by using piezoelectric transducer array

In this paper, a Lamb-wave based structural health monitoring for multi-damage localizations in large composite plates is presented. The Lamb waves are generated and received by piezoelectric transducers, which are arranged in array on the composite plate. In the experiments, three composite plates with various laminate stacking sequences and taper designs were prepared. The damages were created on the specimens by impact testing. In each specimen, 24 piezoelectric transducers were utilized and mounted on the specimen surface. This study proposed an algorithm to identify the damage localizations. The transducer layout is classified by 10 subsets. In each subset, the wave propagation paths can be grouped into path groups pivoted by actuators and that by sensors. Based on the damage index, the mean angle line for each path group in a subset can be obtained. By assuming that the mean angle line passes through the actual damage, the damage localization can be achieved if there exist more than two mean angle lines in one subset. In this study, two exclusion rules are proposed to exclude a path group from the damage localization calculations. The damage localization results show that, for a composite plate with multiple damages, their locations can be identified by using multiple subsets. The damage localization results show that the damage location can be accurately predicted for the case that a damage exists in the interior of a subset. The experiment results also show that the Lamb wave characteristics and the localization results are not affected by the thickness variation of the plate, indicating that the proposed algorithm is available for tapered composite plate.

Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence.Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

Introduction
The Lamb-wave based structural health monitoring (SHM) has attracted attentions in the recent years as it has been applicable in variety of engineering fields, such as aerospace, civil engineering, etc [1,2].Comparing to non-destructive inspection (NDI), the sensors are integrated with or mounted on structures in SHM for real time monitoring.For SHM applications related to the Lamb waves, the waves are generally generated and received by piezoelectric transducers with thin and compact sizes, which are mounted on thin plate structures [3][4][5].For metallic structures, the cracks usually initiate at hot spots such as holes or corners.The piezoelectric transducers are mounted near the hot spots [6][7][8][9].The concept of damage index (DI) has been developed to quantify the change between the base signals and current signals.Based on the correlation between DI and crack size, the DI values can be used to evaluate the crack sizes at hot spots [10][11][12][13][14][15].
The Lamb waves have been also utilized to identify the damage in large plate structures [16].Unlike the damage size monitoring at hot spots, the damage localization in large plate structures needs more sophisticated algorithms.To increase the detectability, multiple piezoelectric transducers are usually layout as an array, which is called as phase array approach.Huan et al utilized the shear horizontal wave generated by piezoelectric transducer to conduct a high-resolution SHM in a large aluminum plate based on phase array approach [17].They proposed a linear phase array and the results showed that it can detect a through-thickness hole with a diameter of 2 mm.It can also detect multiple damages.Yu and Tian presented a generic guided wave phased array beamforming method for anisotropic composite laminates [18].The wave propagation behavior in anisotropic medium was studied in their damage localization model.Senyurek et al proposed a compact phased array for localization of multiple defects in an aluminum plate [19].In their study, the phased array system was composed of three piezoelectric transducers.They proved that their method does not need to apply high voltages to maintain acceptable signal to noise ratios and can detect multiple defects with high accuracy.Leleux et al utilized the phase array approach and developed a Lamb wave-based SHM technique for large composite plates [20].The defects they investigated include delaminations and impact damages.
Another approach for damage localization is the delay-andsum algorithm.This algorithm is based on a fact that a Lamb wave disperses as it propagates through a damage.Wang et al first proposed this approach to detect damages in aluminum plates [21].This approach has been widely utilized to conduct the damage localization in large aluminum and composite plates [22][23][24][25][26].In the delay-and-sum algorithm, the damage localization accuracy relies on the group velocities of Lamb waves.For structures with anisotropy or geometric complexity, it is difficult to develop the Lamb wave propagation model, and the application of delay-and-sum algorithm is limited.
For damage localization and damage imaging, the computerized tomography [27,28], which is inspired by computational tomographic technique used in x-ray imaging, is also a well-developed technique.This technique can offer a damage localization result with high resolution.However, it requires large number of transducers in a relatively small area, i.e. the high transducer density is needed for this technique.The probability-based diagnostic imaging is another approach for presenting damage image by using Lamb waves [29][30][31][32].In this approach, each pixel of the damage imaging corresponds to a calculated probability based on the damage index.The group velocity is no need in the probability calculations so that this approach can be applicable to large complex structures with sparse transducer array.
One of the challenges in some of the previous-mentioned methods in damage localization is that the group velocity analysis model, especially for composite plates for which the group velocity varies with the wave propagation direction.To accurately model the group velocity in composite materials, various approaches have been proposed, including theoretical approaches [33] and experimental measurements [34].For some methods, the accuracy of the damage localization relies on an accurate model for group velocity as the Lamb wave propagates in composite materials.It restricts the applicability to complex structures.
In this paper, an algorithm for multi-damage localizations in large composite plates by using piezoelectric transducer array is developed and presented.In the present approach, the damage indices and propagation path angles are used in the model for damage localization.There is no need to calculate accurate group velocity so that the proposed algorithm can be applicable to anisotropic composite plates.This study also presents the experiments that contain three composite plate specimens with multiple impact damages.The damage localization results calculated by the algorithm are compared with the real damage locations.This paper is organized as follows.
In section 1, the research background and the literature review are presented.The details of the proposed algorithm are given in section 2. In section 3, the experiments, including the specimen preparation, impact testing, piezoelectric transducers and Lamb wave signal acquisitions are presented.The calculations of the algorithm and their comparison to the real damage locations are presented in section 4. The conclusions are given in the last section.

Algorithm for damage localization
The damage localization algorithm developed in this study is to detect the damage locations in a carbon fiber reinforced plastic (CFRP) plate.It should be mentioned that the ultrasonic wave (Lamb wave) used in this study propagates in the inplane directions of the plate.Thus, the algorithm can identify the inplane location of the damages.The depth of the damage along the thickness direction cannot be identified with this algorithm.

Damage index (DI) within first arrival window (FAW)
An input signal modulated by a Hanning window is used to excite the actuator.A Lamb wave emitted from the actuator is generated and received by sensor.The received signal contains the electromagnetic interference (EMI), FAW and successive wave packets, as shown in figure 1(a).In figure 1(a), S b and S denote the signals before and after the damage exists in the specimen, respectively.The time duration T H of the Hanning window is generally identical to that of the FAW, which is formed by the wave packet of fastest Lamb wave mode.The starting time of the FAW can be estimated by calculating the group velocity of the fastest Lamb wave mode.The velocity can be used to determine the distance between the actuator and sensor.The distance must be large enough to avoid the overlap of EMI and FAW.Thus, it is possible to search the actual starting time of FAW between t = T 0 and t = T 1 .The time T 0 is chosen such that T 0 > T H in order to exclude EMI in the FAW identification.The time T 1 such that the time interval from T 0 to T 1 comprises FAW.In this time interval, the maximum baseline signal is denoted by S max .Note that S max occurs in the FAW, not in the EMI.This value is used to normalize the signals shown in figure 1. Starting at t = T 0 , the jth time stamping is denoted by t j = T 0 + (j − 1) ∆t, where ∆t is the time duration between two adjacent time stampings with a given sampling rate.With the notation S b,k being the baseline signal at time t k , consider the absolute value of moving average Sb of the baseline signal.The value of the moving average at time t j is denoted by Sb,j , which is defined by where N 1/4 denotes the number of data points in one-quarter period.The relation of N 1/4 , excitation frequency f and the sampling rate f s is given by The moving average of the signal is calculated by averaging the signal data from t j to t j + T/4, where T/4 denotes the one-quarter period of the sinusoidal excitation signal with frequency f.During the one-quarter period, there are N 1/4 (= f s /4f) data points to be averaged.The moving average Sb is plotted with the time, which is also shown in figure 1(b).In the beginning of the Sb curve, the moving average is less than a threshold value.As the FAW comes, the Sb curve increases and exceeds the threshold value.By setting this threshold value, the algorithm can identify the starting time of the FAW.In this study, the threshold value is set to be 0.3, i.e. the starting time t s of FAW is to find t j with the lowest j such that Sb,j > 0.3.The final time t f of FAW can be estimated by t f = t s + T H .A typical example to search starting time and final time of FAW is shown in figure 1(b).
The baseline signals before and after impact testing are denoted by S b and S, respectively.Theoretically, a difference between S b and S within FAW can be used indicate the existence of a damage.The difference can be evaluated by a damage index DI which is defined by

Mean angle for a path group
In a subset shown in figure 2, the transducers are arrayed in two rows with horizontal spacing w and vertical spacing h.There are N a transducers that are assigned as actuators in the first row, while there are N s transducers that are assigned as sensors in the second row.In this study, the number of actuators and sensors is N a = N s = 4.A lamb wave path is a straight line from an actuator to a sensor.A path group P (i ) a pivoted at the ith actuator in the subset is formed by four paths emitted from the ith actuator and received by the four sensors, as shown in red lines in figure 2. The jth paths belong to P (i ) a are denoted by p (i,j) a .The propagation angle of path p The path group can be also defined by pivoting at the jth sensor, as the blue dashed line shown in figure 2. A path group P (j ) s pivoted by the jth sensor in the subset is formed by paths from the four actuators that are received by the jth sensor.The ith paths belong to P  Thus, this path from ith actuator to jth sensor can be simplydenoted by p (i,j) .The propagation angle for path p (i,j) is denoted by θ (i,j) .It is obvious that The damage index for path p (i,j) is denoted by DI (i,j) .Figure 3 shows the damage indices versus propagation angles for a path group.If the damage is within the subset, the distribution of the damage indices in figure 3(a) should be similar to a Gaussian distribution.The fitting function ϕ Similarly, the fitting function ϕ s (θ) of the data DI values for paths in path group P (j ) s is given by where m are constants and Φ (i ) a (θ) Φ (j ) s (θ) are the probability density functions of the Gaussian distribution: In equation ( 6), θ(i) a and σ a are mean angle and standard deviation for path group P are mean angle and standard deviation for path group P (j ) s , respectively.These statistical parameters can be calculated by where w a and w s are the weighting according to DI (i,j) : By least square method, the constants m and m (j ) s appeared in equation ( 8) can be determined by By using equations ( 5) and ( 6), the peak value φ (i ) a of the fitting function for path group Similarly, the peak value φ (j ) s of the fitting function for path group P (j )

Damage localization
According to the transducer layout, there are eight path groups in one subset, four of which are pivoted by actuators and the rest four of which are by sensors.For the path group a , the damage is assumed to be located on the straight line L (i ) a passing through the ith actuator with a slope corresponding to the mean angle θ(i) a .Similarly, another straight line L (j ) s passing through the jth sensor can also be determined for the path group P (j ) s .The straight lines are called the mean angle lines.For locations of ith actuator and jth sensor denoting by (x s ), respectively, the mean angle lines can be written by a , for path group pivoted at ith actuator (11a) L (j ) s : y = y (j ) s + x (j ) s − x tan θ(j) s , for path group pivoted at jth sensor. (11b) It should be noted that not every path group can be used for damage localization.Some path groups are excluded from the damage localization calculations.The criterion to exclude a path group will be discussed later.Introduce two parameters f (i ) a and f (j ) s defined by Theoretically, by excluding path groups based on the above two exclusion rules, all the preserved mean angle lines pass through the damage location, i.e. these lines are concurrent at this location.However, due to measurement errors or other uncertain reasons, this concurrent does not occur.Suppose the damage location is estimated at a point (x,ȳ).In order to determine this point, consider the distance d (i ) a between (x,ȳ) and the mean angle line for path group P (i ) a .According to equation (11), a can be written by Similarly, the distance d s for path group P (j ) The geometric meaning of R 2 is the summation for the distance square from the point (x,ȳ) to the preserved mean angle lines.Then the solution for (x,ȳ) can be determined by solving the following simultaneous equations: If all the paths are excluded, R 2 is zero and no damage is detected in this subset.
The maximum φ (i ) a and φ (j ) s values for preserved path groups in a subset can be used as an indicator to evaluate the damage size.This maximum value φ max can be calculated by For the case that no path group is preserved, φ max = 0.The dimensionless distance between the detected location (x,ȳ) and the nearest damage location (x d , y d ) is denoted by D, which is defined as where w denotes the horizontal distance between two adjacent transducers, as indicated in figure 2.

Exclusion rules for path groups
In this study, two exclusion rules for path groups are proposed.
The details for the rules are described as follows.

Exclusion rule 1.
The exclusion rule 1 states that a path group is excluded from damage localization calculations if the damage indices for the group do not fit the shape of Gaussian distribution function.For the four paths in one path group, the mean angle line is excluded in the calculation of damage localization if the largest damage index occurs at the first or the fourth path, as shown in figure 3(b).Such damage index distribution usually occurs in the situation that the damage is not located in the range covered by the paths in the path group, or the damage is located on the border of the subset.For this case, the mean angle θ(i) a (or θ(j) s ) may have large error in estimating the damage location and should be excluded.

Exclusion rule 2.
It is assumed that the damage index of a path is greater than a certain level if the path passes through the damage.It is expected that the damage indices of the paths in a path group are small if no damage exists.For this case, the path group should be excluded.Thus, the exclusion rule 2 states that a path group is excluded from the damage localization calculations if the damage indices of this path group is less than a specific value.By introducing a threshold damage index DI th , a path group is excluded from damage localization calculations if a < DI th for path group P (i ) a φ (j ) s < DI th for path group P (j ) s .(18) In another case that a damage occurs within a subset, the peak value φ (i ) a (or φ (j ) s ) exceeds DI th .However, this path group does not directly interest with the damage.For such case, the damage indices of this path group should be lower than those of other path groups in the subset that are directly affected by the damage.Thus, another reference value DI ref is introduced and the exclusion rule 2 is modified by Equation ( 18) can be considered as a special case of equation ( 19) by letting DI ref = 0.If DI ref and DI th are both set to be zero, i.e.DI ref = DI th = 0, it means that the exclusion rule 2 is not applied.In the above procedure, the values of DI th and DI ref are the key parameters that can affect the damage localization results.The threshold damage index DI th can be estimated by observing the largest φ (i ) a and φ (j ) s for paths in a subset without damage within it.The determination of DI ref will be discussed latter.

Determination of DI ref
As mentioned earlier, the introduction of DI ref is used to exclude the path groups that do not directly interact with the damage.The paths within a subset can be categorized into two types based on whether they hit the damage or not.DI ref can be defined as the minimum damage index from the paths that hit the damage.The following procedure provides a simple way to estimate DI ref .
The reference damage index DI ref is calculated subset by subset.The transducers in a subset are arrayed such that the numbers of actuators and sensors are N a and N s , respectively.Therefore, the total number of paths is N p = N a N s , and there are totally N p damage indices in a subset.In this study, N a = N s = 4 and N p = 16.These damage indices can be ranked according to their values for all paths.Suppose that there exists a damage within a subset.By taking all the damage indices in a subset, they can be sorted in descending order and plotted figure as shown in figure 4. In the descending plot, DI n denotes the nth damage index, which means that DI 1 is the largest damage index, DI 2 is the second largest, etc.It is seen that some of the damage indices are significantly higher than the others, indicating that these paths are likely to hit the damage.It is also seen that the DI n curve can be roughly divided by two segments with different slopes, and the point where the slope changes can be assigned to determine DI ref .However, it is not smooth for each of the two segments and the slope change location is not easy to be identified.In order to smooth the curve, introduce DI (ave)   n which is the average value of the damage indices from DI n to DI Np : ) and last point (N p , DI (ave) Np ) make a line that makes an angle of α with the horizontal line.This angle can be calculated by Define the vector [ n DI (ave)   n ] t that started from the origin point (0, 0) to the point (n, DI (ave) n ), where the superscript t denotes the transpose.The rotation of this vector by angle α counterclockwise about the origin point forms the rotated vector [ n (rot) DI (rot) n ] t , which can be calculated according to In general, the maximum damage index is approximately equal to 1 and DI Np .It gives that α is a small angle and The DI (rot)   n plot is shown in figure 4. With the rotation, it is seen that DI

Specimen preparation
Three specimens made of CFRP were prepared for the experiments in this study.Each specimen has identical area of 27 in by 21 in (685.8 mm by 533.4 mm).In applications to aircraft engineering, the thickness of the composite plates is generally non-uniform.The tapered thickness is a common seen design.It is well known that the Lamb wave propagation depends on the plate thickness.Thus, three specimens labeled A, B and C are considered, for which three thickness arrangements are presented: no taper (uniform thickness), T1 thickness taper and T2 thickness taper, respectively.The taper design is illustrated in figure 5.According to thickness, the T1 and T2 tapers can be divided into five and nine taper zones, respectively.For tapered specimen, the bottom surface is not flat due to the thickness variation in two adjacent taper zones.The sensors are adhered on the top surface.
The reference stacking sequences is symbolized by R 0 , for which the detail is listed in table 1.The R 0 laminate has 24 layers, in which a single lamina has a thickness of 0.005 in (0.127 mm).Table 2 lists the details of the stacking sequences of the three specimens.In table 2, R −n denotes a stacking sequence by deleting the last n laminas of R 0 .For specimen A, there is no taper design and the stacking sequence is R −4 (n = 4).It means that the specimen A is stacked by 20 layers.For tapered specimen, the stacking sequence is denoted according to taper zones.For example, the T1 taper is utilized in specimen B and the stacking sequences for taper zones T1-1 to T1-5 are denoted by R -8 /R -6 /R -4 /R -2 /R 0 .All the specimens were fabricated by Aerospace Industrial Development Corporation (Taichung, Taiwan).

Impact damage
In the experiment, each specimen was placed on the dropweight impact testing machine for creating impact damages, as shown in figure 6.The ASTM D7136 test standard was followed.The composite plate was placed on clamping fixture and impacted by an impactor with a mass of 5988 g from a preset height.The impact energy can be determined by impactor mass and the height.Each specimen was conducted two experiments.For the first experiment, the impact testing was conducted at three selected locations.In the second experiment, a second impact testing was conducted again at the same locations in each specimen.Two additional locations in each specimen were selected for impact testing in the second experiment.Figure 7 shows the specimens after impact testing.Table 3 lists the details of the impact testing experiments.By using ultrasonic NDI (Olympus OmniScan SX), a damage area was marked for each damage spot.This area can be approximated by an ellipse with horizontal semi-axis c and vertical semi-axis a.The damage size for each damage in table 3 is denoted by 2c × 2a.

Piezoelectric transducers
Figure 8 shows the transducer array in one specimen.In each CFRP plate, 24 Smart Layer transducers (Acellent Technologies, Sunnyvale, California) numbered from 1 to 24 are mounted on the top surface.The transducers are arrayed in three rows and eight columns.The distances between two adjacent rows and between two adjacent columns are 63.5 mm (2.5 in) and 152.4 mm (6 in), respectively.The damage locations listed in table 3 are also marked in figure 8 for reference.
where V 0 is the peak voltage of the signal and H (t) is the Hanning window: where m is an integer which defines the width of the Hanning window.In this study, m = 5 is used in the experiments.
A system composed of Smart Layer transducers, hardware ScanGenie III and a laptop was utilized for signal excitation and acquisition, as shown in figure 9(a).The software SHM Composite (version 2.18.1, Acellent Technologies, Sunnyvale, California) installed in the laptop was utilized to control the hardware and transducers, as shown in figure 9(b).
The sampling rate of the receiving signals is f s = 48 MHz.
Figure 10 shows the dispersion curves for various laminate with different layer sequences.The group velocities shown in figure 10 are calculated with the code GUIGUW (Graphical User Interface for Guided Ultrasonic Waves).The stiffness matrix C for a lamina with its principal coordinate is given by The density and thickness of one lamina is 1750 kg m −3 and 0.005 in (0.127 mm), respectively.By taking the path 1-9 and path 1-11 as examples, the taper zones for the three specimens include R −4 , R −6 and R −8 .Note that the plate simulation with GUIGUW must be with uniform thickness.The wave propagates through various thickness along path 1-11 for specimens B and C. For simplification, various thickness conditions are considered for the dispersion curves solved with GUIGUW in figure 10.Also note that the wave propagation characteristics vary with the propagation direction due to the anisotropic behavior of the composite plate.According to the results shown in figure 10, it is seen that the group velocity of the S0 mode approximately keeps constant for frequency less than 300 kHz.For a frequency greater than 300 kHz, the group velocity of S0 mode may be no longer the fastest and the trends of the dispersion curves becomes complex.According to the group velocity shown in figure 10, the shortest time of flight (TOF) between any two transducers is approximately ranged from 30 µs to 33 µs.The time span of the EMI is approximately equal to that for Hanning window, which is inversely proportional to the frequency f.For Hanning window with lower frequency, the time span of the EMI could be long and overlap   24) is applied on the piezoelectric actuator and the output voltage can be monitored as the Lamb waves pass through the sensors.According to a mesh sensitivity test, the element size in the model is approximately equal to 0.8 mm. Figure 11 shows the transient displacement response of the Lamb wave generated by transducer 1 with ABAQUS for specimen C. At the beginning of the wave propagation as shown in figure 11(a), it is seen that the group velocities along vertical and horizontal directions are not the same, indicating the anisotropy of the laminate plates.Figure 11(b) shows the first wave packet, the S0 mode, passing through the transducers 9 and 11 at t = 46.2µs. Figure 11(c) shows the displacement contour at t = 76.5 µs, revealing that there exist multiple complicated wave packets, which comprise higher-order modes and S0 mode reflected from the boundaries.Thus, it is difficult to identify the modes of the signals received by sensors after t = 70 µs.Fortunately, these complicated wave packets can be ignored as the FAW only contains the S0 mode, which has the fastest group velocity.For damage localization, some transducers are grouped to belong to one subset.Note that the cover area of a subset cannot be too large.For example, if transducers 1-16 are grouped in one subset, then the wave propagation distance for path 1-9 differs a lot from that for path 1-16.Thus, the comparison of damage indices between these two paths becomes improper.In addition, the signal of a path with long distance may also be weak and become unavailable for damage localization.A subset with too small covering area is also unavailable because of the limited detection area.In this study, ten subsets are assigned, each of which is composed of 8 transducers, as listed in table 4. It should be also noted that the cover area of a subset overlaps that of the adjacent one.This arrangement can increase the detectability of the damage.

Verification of the FEM results of wave propagation
For the transducer layout shown in figure 8, it should be confirmed that the distance between two transducers of a path is long enough that the starting time of the FAW is longer than the time duration of the EMI. Figure 12 shows the signals for paths 1-9 and 1-11.The excitation frequency is f = 250 kHz.For specimen A, no taper design in the layer sequences.For specimens B and C, the taper definitions are shown in figure 5 and table 2. For paths 1-9 and 1-11, the wave propagation paths pass through different thickness region for the three specimens.The FEM results for specimens A, B and C show that the thickness variations do not significantly affect the group velocity of Lamb waves.Figure 12 also shows the measured baseline signals for the two paths of the specimens A, B and C. It is seen that the measured group velocity of the Lamb wave propagating along different paths agree well with the FEM results.For specimen A, the number of layers is 20.For specimen B, the path 1-9 lies on the taper zone T1-2, which has 18 layers.For specimen C, the path 1-9 lies near the border of the taper zones T2-2 and T2-3, which have 18 layers and 20 layers, respectively.The measured results in figure 12 show that no significant difference exists for the starting time of the FAW of path 1-9 in the three specimens.Due to the longer distance, the TOF of the wave for path 1-11 is larger than that for path It should be noted that the group velocity theoretically varies with different plate thickness.Due to the insignificant thickness variation in the tapered plate, the results in figure 12 show that the group velocity is not significantly affected by the thickness variations in the taper designs.The time duration of EMI for f = 250 kHz is T H = 24µs.The starting time of FAW of path 1-9 for f = 250 kHz is about 40 µs.As mentioned in the discussions on figure 11, the FAW contains the S0 mode.For the succeeded wave packets after FAW, they comprise the higher-order modes or the S0 mode reflected from the boundary, and it is difficult to identify the wave modes.For figure 12 shows the wave signals for an excitation frequency of 250 kHz.For the other two driving frequencies, 275 kHz and 300 kHz, the similar results can be obtained.
Based on the FEM results and measured data, it is found that the starting time t s of FAW is longer than the time duration T H of Hanning window.For the current transducer layout, the FAW is not interfered by the EMI.= 90 • .However, the thickness where the paths located are vary with different specimen.For specimen A, the plate has 20 layers with uniform thickness.For specimen B, the path 6-14 is located in taper zone T1-4 (22 layers), while path 11-19 is in taper zone T1-3 (18 layers).For specimen C, the paths 6-14 and 11-19 are located in taper zones T2-7 (20 layers) and T2-5 (24 layers).Recall that the thickness of each layer is 0.005 in (0.127 mm).The two paths for the three specimens are located in the taper zones ranged from 18 layers to 24 layers.It is seen that the starting time of the FAW (S0 mode) has no significant variations for the two paths in the three specimens.Again, the signals in figure 13 show that the layer sequences for different plates have no significant effect on the group velocity of S0 mode.Recall the damage locations shown in figure 13, there is no damage on the path 6-14 for the three specimens, indicating that the signals with damage are identical to the baseline signal.For path 11-19, there is a damage on this path for specimen A and the signal S has a difference by comparing S b , as shown in figure 13(a).For path 11-19 in other two specimens, no damage exists on the path and signal S has no significant difference from S b , as shown in figures 13(b) and (c).= 50.2• .For specimens B and C, the two paths pass through different tapered zones.Similar to the results shown in figure 13, the signals in figure 14 show that the group velocity of S0 mode is nearly independent of the thickness variations in the tapered composite plate.It is also observed that the signal S has significant change from baseline signal S b if there exists a damage on the path.

Damage index results
The dispersion characteristics for Lamb waves propagating through the damage are quite complicated.However, the arrival times of the S0 mode for the baseline signal are equal to that for the signal with damage, as shown in figures 13 and 14.It suggests that the group velocities from the dispersion curves as shown in figure 10 are also valid for waves propagating through damage.
By analyzing the signals from all paths, the damage indices for these paths can be calculated.Figure 15 shows the damage indices for paths in the experiments.By comparing with the damage location shown in figure 8, it is observed large damage indices of the paths that pass through the damages.For example, damages exist on the path 6-16 for the three specimens, and the damage indices for this path in figure 15 are greater than 0.1.It is also observed that the damage indices for path 6-16 after the second impact are greater than those in the first impact.

Reference damage index and threshold damage index
In figure 8, it is observed that no damage exists within subset 3 for the three specimens.The φ max values are shown in figure 16.In order to determine DI th , the values of φ max in figure 16 are calculated from the path groups thar are preserved by applying the exclusion rule 1 only, i.e. the exclusion rule 2 is not applied in the calculations of φ max .It is seen that the φ max values in subset 3 for the first impact experiment are typically lower than 0.005.After the second impact, the φ max values for specimen A and specimen B are greater than 0.01 because two damages A-4 and B-5 are located near subset 3.In specimen C, no damage is located near subset 3 after the second impact.Thus, the φ max values in subset 3 of specimen C for the second impact are typically unchanged by comparison those for the first impact.According to the φ max values in subset 3 after first impact for the three specimens, the φ max values less than 0.005 can be considered that no damage exists within a subset.Thus, the threshold value of damage index DI th is set as 0.005.For comparisons, two additional DI th values, DI th = 0 and DI th = 0.01 are also considered in this study.
Tables 5 and 6 list the reference damage index DI ref and φ max for the subsets in the three specimens.Two sets of φ max values are presented here for each case.The first set of φ max values determined by applying exclusion rule 1 are presented.By applying equation (19) (the exclusion rule 2) on the first set of φ max values and letting DI th = 0.005, the second set of φ max values are also presented.For the φ max values are zero for some subsets, it means that all the path groups have been excluded by applying the exclusion rules.For a subset that a damage exists within it, the φ max values are higher.A typical example is the subset 5 for the three specimens.Large φ max values in subset 5 for all cases are observed.For such cases, no further path groups can be excluded by applying the exclusion rule 2.
Figures 17-19 show the damage localization results for specimens A, B and C, respectively.Three different threshold damage indices (DI th = 0, 0.005 and 0.01) are considered in the calculations.In each specimen, both the actual damage locations and the damage localization locations are marked for comparison.Each specimen was undergone two impact tests.The damage localization results obtained for the first impact and second impact are presented.In general, if a damage is within the transducer layout, it can be detected by one or more subsets.For example, consider damage B-1 shown in figure 18.It can be detected by the paths in subset 6.For the other damages that are within the transducer layout, most of them are identified by the present method.Conversely, damages A-3, A-5, B-3, C-3 and C5 are not within the transducer layout and are difficult to be identified.
Tables 7-15 list the damage localization results for each subset of specimens under various driving frequencies and threshold damage indices.If two or more mean angle lines are preserved in a subset, a damage location prediction (x,ȳ) is obtained.Then the damage identified by the calculation is assumed to be the nearest actual damage to (x,ȳ).Also, the dimensionless distance D between the predicted location and actual location of the damage is also given in tables 7-15.For example, consider the subset 5 in specimen A, the damage localization results of which are listed in table 7. The damage localization location with 250 kHz after the 2nd impact is (504.16,337.12) in mm.This location remains the same as the threshold damage index ranged from 0 to 0.01.The nearest actual damage is A-2, which is located at (501.65, 337.60) in mm.The distance between the predicted damage location actual damage location is 2.56 mm, resulting a dimensionless distance D = 0.040.For subsets that D > 1, the damage localization results are inaccurate and are marked by 'star'.For the case of DI th = 0, some path groups are preserved due to the low threshold value.These path groups have low φ max values, indicating that no significant damage signals detected.It results in more subsets are marked by 'star' for DI th = 0.By increasing the threshold value, say DI th = 0.005 or DI th = 0.01, the number of the subsets marked by 'star' are significantly reduced.It is concluded that by assigning proper value of DI th , the damage localization results are more accurate.
There are two possible cases that result in no damage detected in a subset, for which the fields for damage localization and D in tables 7-15 are left as blanks.The first case is that no path group or only one path group is preserved in that subset.The second is that two path groups are preserved, but the two corresponding mean angle lines are nearly parallel to each other such that the location (x,ȳ) is far away from the subset.
For the damages identified by the subsets, the relationship between damage sizes and φ max values are shown in figure 20.The φ max values are calculated by using DI th = 0.005.Note that the inaccurate damage localization that are marked by 'star' as listed in tables 7-15 are not included in figure 20.A trend can be confirmed that the φ max value increases with the increase of damage size and can be regressed by a straight line with a slope of 0.001852 mm −1 .
Before more detail discussions for the damage localization results, it is worthy to mention that the damages appeared in the specimens can be categorized by four types: no damage within a subset, a damage within a subset, a damage near border of a subset and two damages within a subset.In the following, some selected examples categorized in the four types are discussed in detail.to the predict damage is C-1.The distance between them is 149.1 mm, which is equivalent to a dimensionless distance of D = 2.349.
The small φ max value after first impact for driving frequency 275 kHz in subset 3 of specimen C agrees with the fact that no damage exists in this subset.However, there is a damage detected in figure 21(b).If the threshold damage index DI th is adjusted to 0.005, all the six path groups are excluded and no damage localization is reported within subset 3.
It should be noted that two of the mean angle lines L (4) a and L (6) a in figure 21(b) toward to damage C-1.It means that some of the preserved mean angle lines may be useful for predicting a damage in the exterior of the subset.Figure 22 shows an example for this.The damage indices and predicted damage shown in figure 22 are for subset 1 of specimen A after first impact with driving frequency 250 kHz.No damage exists within this subset.With DI th = 0 and DI th = 0.002237 as listed in table 5, four mean angle lines are preserved and the predicted damage location in mm is (x,ȳ) = (241.2,473.7).The predicted location is near the actual damage A-3 with a dimensionless distance of D = 0.933.This distance is approximately equal to the horizontal spacing between two adjacent transducers.The accuracy of the damage localization results is not as good as the case that a damage exists within a subset.By applying the exclusion rule 1, the value of φ max in subset 1 of specimen A at the first impact is 0.004892 as listed in table 5.This value is lower than 0.005.Thus, by increasing DI th to 0.005, all the four mean angle lines are excluded by applying the exclusion rule 2 and the calculation results suggest that no damage exists within this subset.
As a summary, if the maximum damage index in a subset is approximately equal to or less than a certain level, say 0.005, it can be concluded that no damage exists within the subset.With DI th = 0, it might result in a false damage location within the subset or a predicted damage in the exterior of the subset at low accuracy.

Results for a damage within a subset.
Figure 23 shows the DI variations with path angles and the calculation results of damage localization for subset 5 of specimen B. The damage indices are obtained from signals after 1st impact with driving frequency 275 kHz.In subset 5, the four transducers 5, 6, 7 and 8 are assigned as actuators while the four transducers 13, 14, 15, 16 are assigned as sensors.The value of φ max for this case is 0.141268, as listed in table 5.This value is much greater than the threshold value of 0.005.It indicates that there exists a damage within this subset.In figure 23, it is seen that four path groups are preserved to form four mean angle lines and a predicted damage location is calculated.In figure 23, it is seen that the nearest actual damage is B-2.The damage localization is quite accurate with a dimensionless distance of D = 0.077, as seen in table 11.
Figure 24 shows the results subset 5 of specimen B. The damage indices are obtained with driving frequency 275 kHz.The difference of this case from figure 23 is that the results in figure 24 are analyzed after the second impact.The damage size of B-2 is greater after 2nd impact.The value of φ max for this case is 0.499151, which is higher than that for the case of figure 23.By comparing the result for the first impact, the increasing φ max for second impact suggests that this parameter can reflect the damage level.The damage localization in figure 24 is also similar to that in figure 23, showing high accuracy in prediction of damage location.The dimensionless distance for this case is D = 0.084, as seen In table 11.
In one subset with more than one preserved mean angle lines, an intersection point is formed by taking any two of the lines.There may exist more than one intersection points.With the damage localization results shown in figures 23 and 24, it can be concluded that, by applying the two exclusion rules with DI th = 0.005 within one subset, only one damage can be identified in this subset if the intersection points are distributed in a small region, and the damage can be localized within this region.

Results for a damage near boarder of a subset.
There are two cases for a damage existing near border of a subset.The first is that a damage exists near a transducer of a subset.Figure 25 shows the DI variations with path angles and the calculation results of damage localization for subset 9 of specimen A. The damage indices are obtained from signals after second impact with driving frequency 300 kHz.In figure 25(a), it is seen that the damage indices for path group p (13) a are greater than those for path groups pivoted by the other three actuators since the damage A-4 is near actuator 13.The damage indices for paths 13-22 and 13-23 are the largest as these two paths directly hit the damage.However, distributions of the DI's for path groups pivoted by actuators do not fit the Gaussian distribution curves.For the four path groups pivoted by the sensors, the damage index variations fit the Gaussian distributions and are preserved by applying the two exclusion rules.The damage localization results in an accurate damage location prediction with a dimensionless distance of D = 0.215.
Figure 26(a) shows the DI variations with path angles for subset 7 of specimen B after the second impact with driving frequency 250 kHz.This is another demonstration for a   The second case is that a damage does not exist near a transducer, but near the right or left border of the subset.Two examples will be demonstrated for this case: subset 1 in specimen A at second impact with 275 kHz and subset 4 in specimen B at first impact with 300 kHz. Figure 27 shows the DI variations with path angles for these two examples.For subset 1 in specimen A at second impact with 275 kHz, the damage A-5 is near the left border of the subset, as indicated in figure 8.In figure 27(a), it is seen that no path group is preserved by  applying the two exclusion rules.It is obvious that the path 1-9 has the highest damage index of 0.049 97 as it is the nearest path to the damage.For the other path, the damage indices are much less than that for path 1.The second largest damage index in this subset is 0.009 854.The subset 1 is the most left one in the transducer layout, and the path 1-9 is the most left one in the subset 1.Thus, the information shown in figure 27(a) can conclude that a damage exists one the left side area of the subset.This conclusion coincides with the actual location of damage A-5.Similar results can be seen in figure 27(b), where the DI variations with path angles for subset 4 in specimen B at first impact with 300 kHz are presented.The path 7-15, the most left path in this subset, has the highest damage index of 0.1510, while the damage indices for the other paths are much smaller.The variations of the damage index can be accepted since damage B-2 is located one the right border of this subset.Again, no path group in this case can be preserved.The damage B-2 cannot be identified in subset 4, however, it can be identified by the next subset, which encloses this damage.By summary, the criterion to identify the case with only one damage near boarder of a subset is similar to the case with only one damage within the subset.The intersection points by applying the two exclusion rules with DI th = 0.005 are located in a small region.Figure 28 shows the damage index variations and damage localization results for this subset.The driving frequency and threshold damage index are 250 kHz and 0.005, respectively.By looking at the damage index variations with the path angles as shown in figure 28(a), it is observed that the damage index for path hitting the damage is higher.For example, consider the path group p pivoted at sensor 22.The four paths in this path group are 12-22, 13-22, 14-22 and 15-22.The second path of this path group (13)(14)(15)(16)(17)(18)(19)(20)(21)(22) does not hit the damage.The path 12-22 passes near the damage C-1.Both the other two paths, paths 14-22 and 15-22, hit damage C-4.Thus, the distribution of damage index for path group p  does not follow the Gaussian distribution and has been excluded in the damage localization calculations.For another path group p (12) a , the two damages affect the path-angle dependent DI variation in a similar  manner.However, it follows the Gaussian distribution and is preserved after applying the exclusion rule 1.As a result, there are five path groups are preserved and the intersection points formed by the mean angle lines are distributed in a relatively large region.Based on the present algorithm, only one damage location can be reported.In this case, a predicted damage location is obtained at the middle of the two real damages.
Nevertheless, the predicted damage location is more near to damage C-4 with D = 0.805.
The damages C-1 and C-4 can be separately predicted by subset 7 and subset 10, as shown in figures 29 and 30.The subset 7 encloses only one damage C-1 and has a predicted damage location with D = 0.281.Similarly, the predicted damage location by subset 10 is near to damage C-4 with D = 0.113.It is concluded that for the case that one subset encloses two damages, these locations of the two damages can be separately identified by the neighborhood subsets.
With the observations from figures 28 to 30, a subset is identified that there are two damages within it if the intersection points are distributed in a relatively large region by applying the two exclusion rules with DI th = 0.005.The predicted damage may be located between the two actual damages.The localizations for the two actual damages can be obtained from the results of the two adjacent subsets.

Conclusions
In this paper, the multi-damage localizations for laminated composite plate is conducted by Lamb waves, which are excited and received by piezoelectric transducers.Three specimens were prepared for testing, each of which has a transducer layout composed of 24 piezoelectric transducers.There are twenty layers for specimen A without taper, while taper designs are used in specimens B and C. The three specimens were put in the drop testing machine to create actual damages.By using the algorithm developed in this study, the locations of multiple damages are identified.The conclusions are as follows.
1. Accurate damage localization: The damage localization results show that the damage location can be accurately predicted for the case that a damage exists in the interior of a subset.2. Multiple-damage localization: For a composite plate with multiple damages, their locations can be identified by using multiple subsets.3. Available for tapered composite plate: In theory, the dispersion characteristics depend on the plate thickness.
According to the results from experiments and finite element analysis, no significant wave propagation characteristic is found by comparing the three specimens.It means that the same transducer layout can be applied on specimens A (without taper) and the other two specimens (with tapers).
Except for the above conclusions, additional remarks are also given as follows: 1.In the experiment, the baseline signals are needed in the damage index calculations.It may increase the measurement complexity.However, the proposed algorithm can be available for baseline-free signals if the distributions of such signals in a path group fit a Gaussian distribution for the case that a damage exists in a subset.2. The signals may be affected by environments such as temperatures or loadings.However, the environmental uncertainties can be included in the proposed model.By taking the temperature variation as the example, a database of the baseline signals measured from various temperature environments can be pre-established.Then this database can be used to compensate the damage index calculation if the temperature varies.

Figure 1 .
Figure 1.Signals received by a piezoelectric transducer: (a) signals before and after the damage; (b) signals within the FAW.

s
are denoted by p (j,i) s .The propagation angle of path p (j,i) s is denoted by θ (j,i) s .It should be noted that p (i,j) a and p (j,i) s are actually the same path.They both refer the path from ith actuator to jth sensor.

Figure 2 .
Figure 2. Illustrations for path groups in one subset.
of the data DI values for paths in path group P (i ) a is given by

Figure 3 .
Figure 3. Damage indices versus propagation angles for paths belong to path group: (a) DI distribution similar to Gaussian distribution and (b) DI distribution not similar to Gaussian distribution.

Figure 4 .
Figure 4. Damage index ranking plot in one subset.

Figure 4
Figure 4 also shows the DI (ave) n plot, which shows two segments with different slopes.An apex can be identified at the point with slope change.In the DI (ave) n plot, the coordinate of the nth point can be written as (n, DI (ave) n ).The first point (1, DI (ave) 1

≪ 1 .
With N p = 16, we have N p ≫ DI (ave) Np .It is also seen that the index n = n ref corresponding to the slope change in DI (ave) n plot can be obtained by seeking the minimum of the DI (rot) n plot.Then DI ref is assigned to be the damage index DI n ref corresponding to this minimum.

Figure 5 .
Figure 5. Taper design for the composite plate: (a) T1; (b) T2.The unit for dimensions indicated in the figure is in inch and the equivalent values in millimeters are also provided in the parentheses.

Figure 8 .
Figure 8. Transducer layout and actual damage locations.

Figure 9 .
Figure 9. (a) The CFRP plate and the data acquisition system; (b) the software SHM composite.

Figure 13
Figure13shows the comparisons of baseline signal S b and signal S with damage for paths 6-14 and 11-19.Note that the

Figure 14 shows
Figure 14 shows the comparisons of baseline signal S b and signal S with damage for paths 6-16 and 11-21.The signals are extracted with f = 300 kHz.These two paths have the same path angle, i.e. θ (6,16) a = θ (11,21) a

Figure 12 .
Figure 12.Comparisons of measured signals and computed signals for f = 250 kHz.

Figure 21
(a)   shows the DI variations with path angles and the calculation results of damage localizations for subset 3 of specimen C. The damage indices are obtained from signals after first impact with driving frequency 275 kHz.In subset 3, the four transducers 3, 4, 5 and 6 are assigned as actuators while the four transducers 11, 12, 13, 14 are assigned as sensors.It is seen that the path groups p utilizing exclusion rule 1.In table 5, the φ max values for the six path groups is 0.002039 and the reference damage index DI ref is 0.001109 for this subset.Thus, all the six path groups are preserved for damage localization calculations by applying exclusion rule 2 with DI th = 0.As shown in figure21(b), the red sold line and blue dashed line indicate the mean angle lines L (i ) a and L (j ) s , respectively.The damage localization in mm is (x,ȳ) = (343.7,343.4), which is marked in figure 21(b).The actual damages are also marked in figure 21(b) for comparison.The nearest actual damage

Figure 15 .
Figure 15.Variations of damage index for various paths, driving frequencies and specimen: (a) damage index for paths in subsets 1-5 at first impact; (b) damage index for paths in subsets 6-10 at second impact; (c) damage index for paths in subsets 1-5 at first impact; (d) damage index for paths in subsets 6-10 at second impact.

Figure 17 .
Figure 17.Damage localization results for specimen A.

Figure 18 .
Figure 18.Damage localization results for specimen B.

Figure 19 .
Figure 19.Damage localization results for specimen C.

4. 4 . 4 .
Results for two damages within a subset.For second impact in specimen C, there are two damages in subset 9.

Figure 21 .
Figure 21.Damage localization results for subset 3 in specimen C under first impact.The driving frequency is 275 kHz and DI th = 0. (a) damage index for path groups; (b) actual damage locations, preserved mean angle lines and damage localization.

Figure 22 .
Figure 22.Damage localization results for subset 1 in specimen A under first impact.The driving frequency is 250 kHz and DI th = 0. (a) damage index for path groups; (b) actual damage locations, preserved mean angle lines and damage localization.
13-22.The damage indices for the other three paths are relatively higher.The path-angle dependent DI variation for path group p

Figure 23 .
Figure 23.Damage localization results for subset 5 in specimen B under first impact.The driving frequency is 275 kHz and DI th = 0.005.(a) damage index for path groups; (b) actual damage locations, preserved mean angle lines and damage localization.

Figure 24 .
Figure 24.Damage localization results for subset 5 in specimen B under second impact.The driving frequency is 275 kHz and DI th = 0.005.(a) Damage index for path groups; (b) actual damage locations, preserved mean angle lines and damage localization.

Figure 25 .
Figure 25.Damage localization results for subset 9 in specimen A under second impact.The driving frequency is 300 kHz and DI th = 0.005.(a) Damage index for path groups; (b) actual damage locations, preserved mean angle lines and damage localization.

Figure 26 .
Figure 26.Damage localization results for subset 7 in specimen B under second impact.The driving frequency is 250 kHz and DI th = 0.005.(a) Damage index for path groups; (b) actual damage locations, preserved mean angle lines and damage localization.

Figure 27 .
Figure 27.Damage index distributions versus path angles for (a) subset 1 in specimen A at second impact with 275 kHz and (b) subset 4 in specimen B at first impact with 300 kHz.

Figure 28 .
Figure 28.Damage localization results for subset 9 in specimen C under second impact.The driving frequency is 250 kHz and DI th = 0.005.(a) Damage index for path groups; (b) actual damage locations, preserved mean angle lines and damage localization.

Figure 29 .
Figure 29.Damage localization results for subset 7 in specimen C under second impact.The driving frequency is 250 kHz and DI th = 0.005.(a) Damage index for path groups; (b) actual damage locations, preserved mean angle lines and damage localization.

Figure 30 .
Figure 30.Damage localization results for subset 10 in specimen C under second impact.The driving frequency is 250 kHz and DI th = 0.005.(a) Damage index for path groups; (b) actual damage locations, preserved mean angle lines and damage localization.

Table 1 .
The reference stacking sequence R 0 .

Table 2 .
Definitions of layer stacking sequences and transducer locations for specimens.

Table 3 .
Damage locations and sizes measured by NDI for each damage.

Table 5 .
DI ref and φ max at first impact.

Table 6 .
DI ref and φ max at second impact.The reference damage index for this subset is DI ref = 0.027451.The damage index ranking plot that is used to calculate this DI ref value has been shown in figure 4. By further applying the exclusion rule 2 with DI th = 0.005 and DI ref = 0.027451, it is seen that p

Table 7 .
Damage localization results for specimen A at 250 kHz.

Table 8 .
Damage localization results for specimen A 275 kHz.

Table 9 .
Damage localization results for specimen A 300 kHz.

Table 10 .
Damage localization results for specimen B 250 kHz.

Table 11 .
Damage localization results for specimen B 275 kHz.

Table 12 .
Damage localization results for specimen B 300 kHz.

Table 13 .
Damage localization results for specimen C 250 kHz.

Table 14 .
Damage localization results for specimen C 275 kHz.

Table 15 .
Damage localization results for specimen C 300 kHz.