Hardware-In-The-Loop Simulations of Hole/Crack Identification in a Composite Plate

The technology of hardware-in-the-loop simulations (HILS) plays an important role in the design of complex systems, for example, the structural health monitoring (SHM) of aircrafts. Due to the high performance of personal computers, HILS can provide practical solutions to many problems in engineering and sciences, especially in the huge systems, giant dams for civil engineering, and aircraft system. This study addresses the HILS in hole/crack identification in composite laminates. The multiple loading modes method is used for hole/crack identification. The signals of strains measured from the data-acquisition (DAQ) devices are accomplished by the graphical software LabVIEW. The results represent the actual responses of multiple loading mode tests of real specimens. A personal computer is employed to execute the identification work according to the strain data from DAQ devices by using a nonlinear optimization approach. When all the criteria are satisfied, the final identification results will be obtained. HILS will achieve real time identification of hole/crack in the composite plate by using the actual response measured from the sensors. Not only the size, but also the location and orientation of the crack/hole in a composite plate are successfully identified herein.


Introduction
Studies on structural health monitoring (SHM) are popular because the significant development of computer power, data processing capability, and modern sensor technology make real time detection of fractures possible. The various types of sensors are installed onto or embedded in the object's structure and mass data related to the object's certain conditions are automatically collected [1][2][3]. The safety of the object is determined by the analyses of the collected information. Aircraft structural maintenance is one of the critical operations to guarantee continued airworthiness. Reducing the cost and increasing the safety of air transportation are the main goals in the massive aviation business. The characteristics of high strength and low density for composite materials have highly increased the percentage of composites in aircraft structures. However, composite structures, compared to metallic structures, have more complex damage modes because of their anisotropic properties. X. Chen et al. performed a statistical analysis of Chinese airline maintenance departments in 2012. The records of wing structural damages for two types of Boeing aircraft fleet in a 10-year period were obtained. Dents are the most frequent damage mode (38%), followed by paint peeling off (24%). The damage mode for wind erosion is 12%, 2% for hole damage, and 8% for delamination [3]. Therefore, the detection of defects, such as cracks, holes, and delamination in composite structures is increasingly important during maintenance. In order to enhance the safety of structures and reduce the casualties or property losses, the correct and efficient methods to build the monitoring system of the aircraft health are worthy of research topics. sensitivities with respect to hole or crack sizes, locations, or orientations have been shown in [15]. Large variation of sensitive objective function is the key solution to find the size, location, and orientation identification of a crack or hole in the inverse problem. A general solution satisfying the basic equations of strain-displacement, stress-strain, and equilibrium for a two-dimensional anisotropic linear elastic medium has been presented as: where A = a 1 a 2 a 3 (3) z α = x 1 + p α x 2 , α = 1, 2, 3 (6) u and φ are 3 × 1 column vectors denoting the displacements u 1 , u 2 , u 3 and stress functions φ 1 , φ 2 , φ 3 . Strains can be obtained from ε ij = 1 2 u i,j + u j,i . A and B are eigenvectors of the material properties. The material eigenvalues p α and eigenvectors a α , b α are determined by the following eigenrelations: and C ijks are the elastic constants which are assumed to be fully symmetric and positive definite. Detailed descriptions can be found in [26]. The Green's function of an infinite anisotropic plate containing a traction-free hole subjected to a point forceP applied at pointx is expressed as: a, b are the lengths of the semi-axes of the ellipse. For a straight crack, let b = 0. Equation (5) is the special fundamental solution of boundary element formulation. Consider the real structural problem; the finite domain solution can be obtained by using the boundary element method. This part has been completed by Hwu and Liang [15]. For the crack problem, it is always interesting to know the stress intensity factors at the crack tip which are defined as [27] where r is the distance ahead of the crack tip and can be obtained from stress functions shown in Equation (2). The stress intensity factors can be expressed in terms of remote boundary displacements and tractions [28]. Figure 1 shows the profile of a hole in a square composite plate. The conditions of loading and fixed support and the sensor locations are shown in the figure. The position of hole center is (x, y). Hole size are 2a and 2b. θ is the hole orientation. There are three loading modes in the problem: Open, shear, and tear modes, as shown in Figure 2. In practice, the hole or crack located in the plate will not be always horizontal. We cannot predict the orientation of the hole/crack. We cannot produce pure loading modes. If the pure loading mode is not easily actuated, any kind of independent mixed loading modes can be used. It is more reasonable to give the applied loadings in the directions of X, Y, or Z in Figure 1. Therefore, we use the concept of open, shear, and tear modes and apply the mixed loading modes to identify the hole/crack. completed by Hwu and Liang [15]. For the crack problem, it is always interesting to know the stress intensity factors at the crack tip which are defined as [27] where r is the distance ahead of the crack tip and can be obtained from stress functions shown in Equation (2). The stress intensity factors can be expressed in terms of remote boundary displacements and tractions [28]. Figure 1 shows the profile of a hole in a square composite plate. The conditions of loading and fixed support and the sensor locations are shown in the figure. The position of hole center is (x, y). Hole size are 2a and 2b. θ is the hole orientation. There are three loading modes in the problem: Open, shear, and tear modes, as shown in Figure 2. In practice, the hole or crack located in the plate will not be always horizontal. We cannot predict the orientation of the hole/crack. We cannot produce pure loading modes. If the pure loading mode is not easily actuated, any kind of independent mixed loading modes can be used. It is more reasonable to give the applied loadings in the directions of X, Y, or Z in Figure 1. Therefore, we use the concept of open, shear, and tear modes and apply the mixed loading modes to identify the hole/crack.   Open (Mode I)

Nonlinear Optimization
There are some sensitivity analyses before starting the identification of hole or crack parameters. The variations of normal strain and shear strain functions with respect to the hole/crack size location and orientation are shown in [15]. The aim of this paper is to demonstrate an experiment on hole/crack identification in a composite plate by using measured information from static loading. The strains measured from the nonlinear optimization of the problem can be defined as: where S K denotes the sensor position shown in Figure 1 (k = 1-8). The lower index M is the type of loading mode ( M = 1-3, as shown in Figure 2). i is the selected strain when i = 1 means normal strain ε xx ; i = 2 means shear strain ε xy and i = 3 signifies normal strain ε yy for loading modes M = 1 or 2. i = 1 denotes shear strain ε xz and i = 2 denotes shear strain ε yz when the loading mode M = 3. The x, y, a, b, and θ are the parameters of hole or crack which need to be identified. The insensitivities of the static strains to the hole/crack geometry and location may be overcome by the multiple loading modes. We may switch the loading condition to another loading mode and improve the search for the hole/crack geometry and location when the hole/crack cannot be identified by the static strains under a certain single loading condition. Repeat the process until the convergence criterion is satisfied. The multistep nonlinear optimization can be designed as follows. Minimize objective function, Φ M : For loading mode M (= 1-3) subject to: α kM i is the strain value of the reference problem. To avoid any numerical ill conditions, the value of should be normalized to be on the order of unity. L is the number of sensors (L = 8 in Figure 1). e is the error tolerance of convergence which can be decreased in each step, as shown in the flowchart of Figure 3. The upper bounds and the lower bounds of designed variables are: W, H are width and height of the square plate, respectively. The upper bounds of a and b are chosen to be not too large: a u = b u = W/4. The lower bounds of a and b are chosen to be a very small value 10 −6 because the crack is the special case of elliptical hole when b approaches zero. To prevent a multi-value representation, θ u = π/2 and θ l = −π/2. The location bounds are defined as: Equations (28) and (29) express the conditions to enforce the hole or crack inside the plate. This important program is the identification of hole/crack size a, b, location x, y, and orientation θ which will be one system in the hardware-in-the-loop. Figure 3 shows the flowchart of nonlinear optimization with multiple loading modes.
Materials 2020, 13, x FOR PEER REVIEW 6 of 15 The detailed simulations of hole/crack identification by static strains from multiple loading modes can be found in [15]. To save the cost and get the optimal design of HILS, many simulations should be done. According to the simulation results of [15], four, eight, twelve, and sixteen sensors were designed. The range of sensor spacing to the plate width are 0.

Hardware-In-The-Loop Simulations
In this paper, the carbon composite plates with an elliptical hole or a crack need to be made previously. The fiber orientations of six layers of carbon are [0/90/0]s. There are isotropic, orthotropic, anisotropic, and piezoelectric materials designed in our nonlinear optimization programs. Here, we choose the quasi-isotropic plate not only is it easy to manufacture but also it can reduce the sources of errors. The size of the plate is 30 × 30 cm 2 . The heat-compressor machine and a composite plate Optimal values of The detailed simulations of hole/crack identification by static strains from multiple loading modes can be found in [15]. To save the cost and get the optimal design of HILS, many simulations should Materials 2020, 13, 424 7 of 15 be done. According to the simulation results of [15], four, eight, twelve, and sixteen sensors were designed. The range of sensor spacing to the plate width are 0.8, 0.4, 0.267, and 0.2, and the maximum errors of small hole (2a/W = 0.044, 2b/W = 0.022) identified results are 17.1%, 4.7%, 8.86%, and 4.5%, respectively. The proper ratio of sensor spacing to the plate width is 0.4. Hole/crack identifications of various materials, such as isotropic, anisotropic, piezoelectric materials were accomplished and the maximum error is 14.1% for piezoelectric material. Finally, the critical size of the hole/crack is suggested to be 2a/W = 0.03. The experiments of HILS are based on the provided information to manufacture the tested specimen and arrange the locations and number of sensors.

Hardware-In-The-Loop Simulations
In this paper, the carbon composite plates with an elliptical hole or a crack need to be made previously. The fiber orientations of six layers of carbon are [0/90/0] s . There are isotropic, orthotropic, anisotropic, and piezoelectric materials designed in our nonlinear optimization programs. Here, we choose the quasi-isotropic plate not only is it easy to manufacture but also it can reduce the sources of errors. The size of the plate is 30 × 30 cm 2 . The heat-compressor machine and a composite plate with an elliptical hole are shown in Figure 4. The structure of HILS in this paper is shown in Figure 5. There are three steps in the HILS.
the left and the middle of Figure 6. There are eight sensor positions in this plate. The type of strain gauge is KFG-5-120-C1-11 produced by KYOWA electronic instruments. The DAQ card NI9237 with four channels is shown on the right-hand side of Figure 6. The loading time profile is shown in Figure  7. The designed maximum load of the MTS is 5 kN and the stop loading in this paper is 1 kN.
Step 2: Usually, there are eight sensor positions, as shown in Figure 1. In Equation (20), strains , , and for loading mode I or II and strains and for loading mode III are required in each sensor. There will be 24 channels of strain measurement in loading modes I and II, and 16 channels of strain measurement in loading mode III. As it makes no difference whether the strain gauges are embedded into the plate, the strain gauges were set on the top surface of the plate. The sensors embedded in the plate for mode III are a problem. The plate is considered as plane stress for thin thickness. The fractures may generate during embedding the strain gauges into the plate. Therefore, only loading modes I and II is considered in this paper. To reduce the cost of future development, the strain is eliminated to leave only the 16 channels (four DAQ cards) required. Figures 8 and 9 show the execution of strain measurement in the composite plate.
Step 3: The core process in the HILS is the program of nonlinear optimal identification in the personal PC, as shown in Figure 5. The algorithm has been described in Sections 2 and 3. The real response of the composite plate with a hole/crack subjected to loading mode I or II is forwarded into the nonlinear optimization program. The program starts from an initial guess of hole/crack size, location, and orientation (a, b, x, y, θ) under the situation of loading mode I or II and stops at the conditions when all the criteria of Equations (22)(23)(24)(25)(26)(27)(28)(29) are satisfied. If the criteria are not satisfied, the optimal search will go back to another loading mode by using the final result of the previous loading mode to be the initial guess, and execute the loop till the program is convergent.   Step 1: The material test stand of multiple loading modes (loading mode I and loading mode II) is shown on the left and the middle of Figure 6. There are eight sensor positions in this plate. The type of strain gauge is KFG-5-120-C1-11 produced by KYOWA electronic instruments. The DAQ card NI9237 with four channels is shown on the right-hand side of Figure 6. The loading time profile is shown in Figure 7. The designed maximum load of the MTS is 5 kN and the stop loading in this paper is 1 kN.      Step 2: Usually, there are eight sensor positions, as shown in Figure 1. In Equation (20), strains ε xx , ε xy , and ε yy for loading mode I or II and strains ε xz and ε yz for loading mode III are required in each sensor. There will be 24 channels of strain measurement in loading modes I and II, and 16 channels of strain measurement in loading mode III. As it makes no difference whether the strain gauges are embedded into the plate, the strain gauges were set on the top surface of the plate. The sensors embedded in the plate for mode III are a problem. The plate is considered as plane stress for thin thickness. The fractures may generate during embedding the strain gauges into the plate. Therefore, only loading modes I and II is considered in this paper. To reduce the cost of future development, the strain ε xy is eliminated to leave only the 16 channels (four DAQ cards) required. Figures 8 and 9 show the execution of strain measurement in the composite plate.      Step 3: The core process in the HILS is the program of nonlinear optimal identification in the personal PC, as shown in Figure 5. The algorithm has been described in Sections 2 and 3. The real response of the composite plate with a hole/crack subjected to loading mode I or II is forwarded into the nonlinear optimization program. The program starts from an initial guess of hole/crack size, location, and orientation (a, b, x, y, θ) under the situation of loading mode I or II and stops at the conditions when all the criteria of Equations (22)-(29) are satisfied. If the criteria are not satisfied, the optimal search will go back to another loading mode by using the final result of the previous loading mode to be the initial guess, and execute the loop till the program is convergent.

Hole Identification
The parameters of the hole are a = 9 mm, b = 6.5 mm, x = 200 mm, y = 153 mm, and θ = 45 • . The strains calculated from the numerical solutions of the boundary element method are shown in Table 1. Table 2 shows the strains measured in the HILS. In the bottom of Table 2, the value of Φ M indicates the differences between the numerical data ( Table 1) and experiment data ( Table 2). It implicates that Φ M is small. The values of ε xx and ε yy in Tables 1 and 2 are close. The objective function Φ 1 M is calculated by a different hole. The large value of Φ 1 M means the objective function of the different hole is sensitive. The identified results of the hole by using the multiple loading modes are presented in Table 3. The errors increase when identifications do not use only pure simulations. The source of errors may be the loosened jigs of MTS. The jigs were sometimes slightly loosened when the plate was removed. The other error source may be the installations of the strain gauges. The identified procedures are plotted in Figure 10. IG means initial guess. TG means target. BEM means identifying a hole by using pure simulations. HILS is the identification results of a hole by using HILS. No. 1 is the identified results by setting the error tolerance e in Equation (6) 1 (mode I). The error tolerance e = 0.5 (mode II) for No.2 and e = 0.1 (mode I) for No.3. The sequence convergence criterion of error tolerance e in Equation (6) is 1 (mode I)→0.5 (mode II)→0.1 (mode I)→0.05 (mode II)→0.01 (mode I)→0.005 (mode II)→0.001 (mode I)→0.0005 (mode II). The final result of each loading mode is the initial guess of the next loading mode until the convergence criterion is satisfied. The number marked in Figure 10 signifies these seven identified processes.

Crack Identification
It is well known that crack identification is more difficult than hole identification because of the nonsensitive response of crack parameters to the strains around the plate [15]. Tables 4 and 5 show the strains calculated from BEM and the strains obtained from the experiment, respectively. The objective function of another crack is presented in the last line of Table 5. The identified results by using the numerical method or HILS are shown in Table 6 and plotted in Figure 11. IG means initial Figure 10. Procedures of hole identification by using multiple loading modes. Red points are sensor positions. The marked number signifies these seven identified processes.

Crack Identification
It is well known that crack identification is more difficult than hole identification because of the nonsensitive response of crack parameters to the strains around the plate [15]. Tables 4 and 5 show the strains calculated from BEM and the strains obtained from the experiment, respectively. The objective function of another crack is presented in the last line of Table 5. The identified results by using the numerical method or HILS are shown in Table 6 and plotted in Figure 11. IG means initial guess. TG means target. BEM means identifying a crack by using pure BEM. HILS is the identification results of a crack by using HILS. The crack identification by using HILS is not good. The errors increase to 22.9% for identified orientation, 21.48% for size identification, and 18.97% for vertical position identification; this may result from the open mode effect of a larger crack with larger orientation being the same as the one of a smaller crack with smaller orientation. In shear mode, a vertical crack is more sensitive than a horizontal crack. On the other hand, the design of jig support may be improved by some fixtures to completely clamp the plate. The sequence convergence criterion of error tolerance e is 1 (mode I)→1 (mode II)→0.5 (mode I)→0.5 (mode II)→0.1 (mode I)→0.05 (mode II)→0.01 (mode I)→0.005 (mode II)→0.001 (mode I)→0.0005 (mode II)→0.0005 (mode I).

Conclusions
Strains are the most convenient data to acquire in a static state. This paper demonstrates the feasibility of identifying a hole/crack by using only the strains measured from the composite plate. The experiment method is the hardware-in-the loop simulations (HILS). HILS include: A material test stand to give multiple loadings and fixed support to the composite plate having a hole/crack inside, four DAQ cards to acquire the strains around the inner plate boundary, a well-developed nonlinear optimization program to stably search the global minimum results. The experiment of hole/crack identification by using multiple loading modes is performed and the identified results of hole/crack size, location, and orientation are satisfied. The errors of HILS (from the experiments) are larger than the ones of BEM (from the simulations). Moreover, any kind of independent mixed loading mode can be applied if the pure loading mode is not easily actuated. This conforms to the real loading situations of structures. The accuracy may be increased when the independent mixed loading modes are added into HILS.

Conclusions
Strains are the most convenient data to acquire in a static state. This paper demonstrates the feasibility of identifying a hole/crack by using only the strains measured from the composite plate. The experiment method is the hardware-in-the loop simulations (HILS). HILS include: A material test stand to give multiple loadings and fixed support to the composite plate having a hole/crack inside, four DAQ cards to acquire the strains around the inner plate boundary, a well-developed nonlinear optimization program to stably search the global minimum results. The experiment of hole/crack identification by using multiple loading modes is performed and the identified results of hole/crack size, location, and orientation are satisfied. The errors of HILS (from the experiments) are larger than the ones of BEM (from the simulations). Moreover, any kind of independent mixed loading mode can be applied if the pure loading mode is not easily actuated. This conforms to the real loading situations of structures. The accuracy may be increased when the independent mixed loading modes are added into HILS.