Structural health monitoring of 3D frame structures using finite element modal analysis and genetic algorithm

1. Introduction

In the last years, Genetic Algorithm (GA) has been recognized as a promising intelligent artificial technique for difficult optimization problems. GA is considered to be more efficient than classical methods because it is based on precise results. Sazonov et al. [2] used strain energy mode shapes to determine the changes in structural integrity from changes in the vibrational response of a structure using (GA) to produce a sufficiently optimized amplitude characteristics. The detection of structural damage is an inverse problem in structural engineering. In this problem, it is very important to determine the existence, location and extent of the damage. This problem is formulated as an optimization problem, which is then solved using GA [3].

To identify damage using GA, three objective functions were used [2, 4], i.e. a) the changes in natural frequencies, b) Modal Assurance Criterion (MAC) and c) MAC in the frequency domain. In these references, GA based on the strain energy mode shapes was used and a finite element model was constructed to generate training data. The results obtained from GA were efficient and confirmed with the theoretical predictions. Radzieński et al. [5] used only the changes of natural frequencies without mode shapes as the measurements of frequencies were much less time consuming in comparison with the measurements of mode shapes. The benchmark problem on structural health monitoring strategy was applied to vibration data generated with an analytical model of a benchmark structure [6, 7]. The subject of the experimental benchmark problem was the 4-story, 2-bay by 2-bay steel-frame scale-model structure by Dyke et al. [8]

It has been shown that the proposed technique was suitable for damage localization in beam-like structures. The investigations of damage detection and localization based on modal parameters and the combination of natural frequency and mode shape were widely utilized [9-19]. The structural damage detection and localization was defined as a permanent change in the mechanical state of a structural material or component that could potentially affect their performance [20]. Rizos et al. [21] presented a method based on amplitudes of two points of a cantilever beam vibrating at one of its natural modes to identify crack location and depth. Damage detection and localization in thin plates based on vibration analysis using BAT algorithm was presented by Khatir et al. [22]. A comparison between a frequency based and a mode shape based method for damage identification in beam like structures has been published by Kim et al. [23]. A new frequency domain technique was introduced for modal identification of output-only systems, in the case where the modal parameters must be estimated without knowing the input excitation of the system. This technique was applied to two-storey building model [24]. The FE model of a multi-story frame with rotational stiffness at joints under earthquake excitation for structural health monitoring was investigated by Lei et al. [25].

The classical optimization methods are fast, but they suffer from big limitations related to the continuity of the objective function, e.g. the Hessian of the objective function, which might not be positive-defined in all points during optimization, and the substantial possibility of getting a local optimum, which strongly depends on starting point [26]. Genetic Algorithm (GA) [27] and the Particle Swarm Optimization (PSO) [28] are well-known as optimization techniques that are free from the above mentioned restrictions. The new conception of natural frequency vector (NFV) and natural frequency vector assurance criterion (NFVAC), a new damage detection method using the natural frequencies as damage was introduced by Yang [1].

In this study, a GA approach is used for damage detection and localization in 3D complex frame structures by minimizing the objective function, which is based on the changes in natural frequencies. The damage detection algorithm is investigated in case of removed frame elements. The approach is validated using experimental benchmark test from literature.

2. Optimization problem

Genetic algorithm is developed by Holland [29]. In this optimization method, information about a problem, such as variable parameters, is coded into a genetic string known as an individual (chromosome). In this study, the numbers of chromosome presented in two elements are: a) the presence of removed elements and b) the locations of removed elements. Each of these individuals has an associated fitness value, which is usually determined by the objective function to be minimized. Genetic algorithm is used as inverse problem and can be coupled with FEM of 3D structure. It has been shown that to be able to solve the optimization problem through mutation, crossover and selection operation should be applied to individuals in the population (see Fig. 1). In this study, we address the problem for damage detection by removing beam elements from 3D frame structures. As mentioned above, in GA approach, each of the 100 individuals contains two chromosomes representing damage parameters. The maximum number of iteration was set equal to 100. After several applications, a crossover coefficient of 0.8 and mutation of 0.1 were used in the GA parameters.

3. Numerical example

A simulated 3D frame structure, as shown in Fig. 2, is used to verify the proposed technique. The frame model is divided into ten beam elements with 6 DOFs for each node, as explained in section 3.1. The properties of the beam element are listed in Table 1.

Fig. 1Flowchart for damage detection in 3D Structures

Fig. 2The finite element model of the 3D frame structure used in the numerical example

3.1. Finite element analysis of 3D frame elements

In the 3D frame element, we consider in each node, 6 degrees of freedom (6 DOFs), as shown in Fig. 3, namely three displacements and three rotations with respect to the three global Cartesian axes; $X$, $Y$, $Z$.

Fig. 33D frame element with 6 DOFs per node

In the local coordinate system, the stiffness matrix, 12×12, of the 3-D frame element is given by:

1

${\left[K\right]}^e=\left[\begin{array}{c}\begin{array}{c}{r}_1{I}_z\\ 0\\ 0\end{array}\\ 0\\ 0\\ 0\\ {-r}_1{I}_z\\ 0\\ 0\\ 0\\ 0\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ 12I_z\\ 0\end{array}\\ 0\\ 0\\ 6l_e{I}_y\\ 0\\ -12I_z\\ 0\\ 0\\ 0\\ 6l_e{I}_z\end{array}\begin{array}{c}\begin{array}{c}\\ \\ 12I_y\end{array}\\ 0\\ 6l_e{I}_y\\ 0\\ 0\\ 0\\ -12I_y\\ 0\\ -6l_e{I}_y\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ r_2{l}_e^2{I}_y\\ 0\\ 0\\ 0\\ 0\\ 0\\ {-r}_2{l}_e^2{I}_y\\ 0\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ 4l_e^2{I}_y\\ 0\\ 0\\ 0\\ 6l_e{I}_y\\ 0\\ 2l_e^2{I}_y\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ \\ 4l_e^2{I}_z\\ 0\\ -6l_e{I}_z\\ 0\\ 0\\ 0\\ 2l_e^2{I}_z\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ \\ \\ r_1{I}_z\\ 0\\ 0\\ 0\\ 0\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ \\ \\ \\ 12I_z\\ 0\\ 0\\ 0\\ -6l_e{I}_z\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ \\ \\ \\ \\ 12I_y\\ 0\\ 6l_e{I}_y\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ Sym\end{array}\\ \\ \\ \\ \\ \\ \\ r_2{l}_e^2{I}_y\\ 0\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ \\ \\ \\ \\ \\ \\ 4l_e^2{I}_y\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ \\ \\ \\ \\ \\ \\ \\ {4l}_e^2{I}_z\end{array}\right],$

where $r_1=A{l}_e^2/I_z$ and $r_2=\left(1+v\right)I_x/{2I}_y$, $A$ is the cross-sectional area, $I_x$ the torsional constant, $I_y$, and $I_z$ are the second moments of inertia with respect to the local axes $y$ and $z$, respectively, and $v$ is Poisson’s ratio.

The mass matrix is evaluated according to the consistent formulation and is given in an explicit form with respect to the element coordinate system as:

2

${\left[M\right]}^e=\frac{\rho A{l}_e}{420}\left[\begin{array}{c}\begin{array}{c}140\\ 0\\ 0\end{array}\\ 0\\ 0\\ 0\\ 70\\ 0\\ 0\\ 0\\ 0\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ 156\\ 0\end{array}\\ 0\\ 0\\ 22l_e\\ 0\\ 54\\ 0\\ 0\\ 0\\ -13l_e\end{array}\begin{array}{c}\begin{array}{c}\\ \\ 156\end{array}\\ 0\\ -22l_e\\ 0\\ 0\\ 0\\ 54\\ 0\\ 13l_e\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ 140r_g^2\\ 0\\ 0\\ 0\\ 0\\ -13l_e\\ 0\\ -3l_e^2\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ 4l_e^2\\ 0\\ -13l_e\\ 0\\ 0\\ 0\\ -3l_e^2\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ \\ 4l_e^2\\ 0\\ 13l_e\\ 0\\ 0\\ 0\\ -3l_e^2\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ \\ \\ 140\\ 0\\ 0\\ 0\\ 0\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ \\ \\ \\ 156\\ 0\\ 0\\ 0\\ -22l_e\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ \\ \\ \\ \\ 156\\ 0\\ 22l_e\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ Sym\end{array}\\ \\ \\ \\ \\ \\ \\ 140r_g^2\\ 0\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ \\ \\ \\ \\ \\ \\ 4l_e^2\\ 0\end{array}\begin{array}{c}\begin{array}{c}\\ \\ \end{array}\\ \\ \\ \\ \\ \\ \\ \\ \\ 4l_e^2\end{array}\right],$

where $\rho$ is the mass density, and $r_g=J/A$ is the radius of gyration. After transformation to the global axes, the stiffness and mass matrices in global coordinates are obtained as:

3

$K=R^T{K}^{e}R,M=R^T{M}^{e}R,$

where the transformation matrix $R$ is defined as:

$R=\left[\begin{array}{cccc}\gamma & 0& 0& 0\\ 0& \gamma & 0& 0\\ 0& 0& \gamma & 0\\ 0& 0& 0& \gamma \end{array}\right],\gamma =\left[\begin{array}{ccc}{C}_{Xx}& C_{Yx}& C_{Zx}\\ C_{Xy}& C_{Yy}& C_{Zy}\\ C_{Xz}& C_{Yz}& C_{Zz}\end{array}\right].$

$C_{Xx}=\cos {\theta }_{Xx},$

where angles ${\theta }_{Xx}$, ${\theta }_{Yx}$, and ${\theta }_{Zx}$, are measured from global axes $X$, $Y$, and $Z$, with respect to the local axis x, respectively.

The dynamic behaviour of a linear mechanical structure is governed by the following equation:

4

$M\ddot{x}+C\dot{x}+Kx=F.$

Ignoring damping and external force terms, Eq. (4) can be written as:

5

$M\ddot{x}+Kx=0,$

where $M$ and $K$ are real symmetric matrices, which are discretized as follows:

6

$M=\sum _{i=1}^{N_e}{M}_i^{\left(e\right)},K=\sum _{i=1}^{N_e}{K}_i^{\left(e\right)}.$

Let:

7

$x\left(t\right)=s\sin \left(\omega t+\varphi \right),$

where vector $s$ is $n\times 1$, the frequency $\omega$, and the phase $\varphi$ can be determined.

Differentiating Eq. (7) twice with respect to time, gives:

8

$\ddot{x}=-{\omega }^{2}s\sin \left(\omega t+\varphi \right).$

Substituting Eqs. (7) and (8) into Eq. (5), gives:

9

$\left(K-{\omega }^{2}M\right)s=0.$

For a nonzero or a nontrivial solution of $s$:

10

$det\left(K-{\omega }^{2}M\right)=0.$

Which will be a polynomial equation of degree $n$ in ${\omega }^2$. Eq. (9) can also be written as:

11

$Ks={\omega }^{2}Ms.$

Or:

12

$M^{-1}Ks={\omega }^{2}s.$

Eq. (12) clearly indicates that ${\omega }^2$ and $s$ are an eigenvalue and eigenvector of the matrix $M^{-1}K$, respectively. In addition, Eq. (11) suggests that ${\omega }^2$ and $s$ are generalized eigenvalues and eigenvectors of the stiffness matrix $K$ with respect to the mass matrix $M$.

3.3. Damage scenarios

In order to validate the proposed damage assessment technique, four damage scenarios are considered, in which single damage, as well as, multiple damage cases are studied as shown in Table 2 (see the 3D frame structure in Fig. 2). In Table 2, element numbers, where the stiffness is reduced, are listed along with and the percentage of reduction (between parenthetic). The naturel frequencies of each damage scenario are presented in the Table 3.

Table 2The 3D frame structure damage scenarios

Damage scenario	Damaged element (% reduction in bending stiffness)
D1	4 (30 %)	–	–	–	–
D2	9 (20 %)	10 (35 %)	–	–	–
D3	7 (25 %)	8 (30 %)	13 (35 %)	15 (20 %)	–
D4	3 (25 %)	5 (30 %)	7 (30 %)	13 (25 %)	15 (35 %)

Table 3Naturel frequencies of undamaged and damaged 3D frame structure

Mode	Undamaged $f$ [Hz]	Damage scenario
		D1		D2		D3		D4
		$f$ [Hz]	%	$f$ [Hz]	%	$f$ [Hz]	%	$f$ [Hz]	%
1	107.723	105.364	2.2	103.629	3.8	107.563	0.1	105.624	1.9
2	142.545	136.840	4.0	142.545	0.0	140.895	1.2	136.232	4.4
3	159.732	155.896	2.4	159.177	.03	158.688	0.7	155.898	2.4
4	435.371	433.476	0.4	416.140	4.4	426.295	2.1	422.780	2.9
5	484.670	475.107	2.0	471.982	2.6	481.979	0.6	474.043	2.2
6	586.661	583.594	0.5	586.654	0.0	548.766	6.5	527.674	11.1

The first six natural frequencies, listed in Table 3, are utilized in this case for all damage scenarios to calculate the fitness, and consequently the damage locations and their severities. The identified results are shown in Figs. 4 to 7, for the different damage scenarios.

Fig. 4Damage scenario D1

a) Convergence of fitness

b) Damage identification of location and severity

c) Convergence of damaged element

In Figs. 4 to 7, three graphs are plotted: a) Convergence of fitness, b) Damage identification of locations and severities and c) Convergence of damaged elements. For the single damage scenario, D1, it can be seen in Fig. 4 that the algorithm is converged after few iterations and the damage location and severity are correctly identified. Similarly, for multiple damage scenarios, scenarios D2 to D4, all damage locations and severities are correctly found. However, as the number of damage elements increases, the number of iterations increases and the convergence becomes slower. For damage scenario D4, more than 70 iterations are required in order to reach convergence, as it can be seen in Fig. 7.

Fig. 5Damage scenario D2

a) Convergence of fitness

b) Damage identification of locations and severities

c) Convergence of damaged elements

Fig. 6Damage scenario D3

a) Convergence of fitness

b) Damage identification of locations and severities

c) Convergence of damaged elements

Fig. 7Damage scenario D4

a) Convergence of fitness

b) Damage identification of locations and severities

c) Convergence of damaged elements

3.4. Effect of noise

In order to investigate the effect of noise on our damage detection technique, White Gaussian noise was added into damage scenario D2 and D4 in the first six modes with 5 % and 10 %. The $i$th noisy response $N_{di}$ (noise), is simulated by [15]:

14

${N_{d_i}\left(Noise\right)=(1+\sigma \gamma )N}_{d_i},$

where $\sigma$ is the noise level and $\gamma$ is a random number in the interval [−1, 1]. From the results shown in Fig. 8, we can observe that, when the noise is included in the problem of fault detection, our approach based on GA can detect damage with high accuracy. However, the severity of damage is affected by the level of noise.

Fig. 8Damage location for scenarios D2 and D4 with noise level 5 % and 10 %

a) Damage scenario D2

b) Damage scenario D4

4. Experimental validation

In order to validate the proposed technique, experimental benchmark test of an eight-story shear frame is used from literature [1]. The eight-story frame and the experimental set-up are shown in Fig. 9. In this experimental setup, the frame was excited using a MB Modal 50A electromagnetic exciter driven by a sweeping sine signal generated by a GW GFG-8019G signal generator and amplified by a MB SS250VCF amplifier. The vibration responses were measured using PCB 3330B accelerometers and sampled using a DIFA/S CADAS acquisition system. The first six frequencies were identified in this case. Two local damages were assumed to be located at the 54th and 63rd elements. The damage was introduced by removing elements. The finite element model of the eight-story shear frame is shown in Fig. 10. The calculated natural frequencies are compared with the ones obtained from the experimental model [1] in Table 4. The dimension of each beam is 139 mm × 27 mm × 1 mm. The damage is introduced by removing one beam element. The frequencies of undamaged 3D frame model are presented in the Table 4. The loss of rigidity was obtained by multiplying the rigidity matrix by rigidity coefficients [0-1] in Section 3. However, in this section, the damage is simulated by removing elements, therefore the rigidity coefficient is equal to 0.

The results are shown in Figs. 11-13 for the considered damage cases. Damages have been identified accurately in both locations. In this analysis, the maximum number of iteration is set equal to 500 with 1000 population. The proposed technique based on GA and FEM is efficient to determine the removed element. The results show that our approach using the natural frequencies of damaged 3D frame structure can detect damage with high accuracy.

Fig. 9The experimental eight-story frame [1]

Fig. 10The FE model of the 8-story shear frame

Fig. 11Convergence and loss of rigidity for the case of removing element 54

a)

b)

Fig. 12Convergence and loss of rigidity for the case of removing element 63

a)

b)

Fig. 13Convergence and loss of rigidity for the case of removing elements 54 and 64

a)

b)

5. Conclusions

In this article, a method for inverse problem is proposed in order to quantify damage in 3D complex frame structure. The proposed technique based on GA coupled with FEM and the objective function is defined as the difference between calculated by GA and measured natural frequencies. The results show clearly that the proposed methodology can be used to quantify damage in case of single, as well as, multiple damage scenarios. From the numerical example, in which noise was considered, the comparison between the estimated and real damage illustrated the efficiency of the algorithm in damage detection. In the last section, we validated our technique using experimental benchmark test from literature, in which damage was introduced by removing elements. Good results were obtained.