Elliptical Crack Identification in a Nonrotating Shaft

It is known that fatigue cracks are one of the most important problems of the mechanical components, since their propagation can cause severe loss, both personal and economic. So, it is essential to know deeply the behavior of the cracked element to have tools that allow predicting the breakage before it happens. ,e shafts are elements that are specially affected by the described problem, because they are subjected to alternative compression and tension stresses.,is work presents, firstly, an analytical expression that allows determining the first four natural frequencies of bending vibration of a nonrotating cracked shaft, assumed as an Euler–Bernoulli beam, with circular cross section under pinned-pinned conditions, taking into account the elliptical shape of the crack. Second, once the direct problem is known, the inverse problem is approached. Genetic Algorithm technique has been used to estimate the crack parameters assuming known the natural frequencies of the cracked shaft.


Introduction
Shafts are a very common mechanical components that are usually subjected to very hard working conditions, specifically bending loads that originate periodic stresses, changing from tensile to compression stresses for each cycle of the shaft, often accompanied by high temperatures and aggressive environments.All of that frequently causes incipient fatigue cracks which can quickly lead to catastrophic failures [1][2][3].So, it is very important to use diagnosis methodologies that allow detecting possible cracks and also predicting both the shape and size of cracks in order to replace the part before reaching the total breakdown.When a crack appears in a mechanical component it produces a change in the modal properties (natural frequencies and modes shapes) that is accompanied by an increase in its flexibility [1,4].Over the last four decades, from the 70s to the present day, many researches have focused on developing nondestructive maintenance methods founded on the vibration measurement characteristics, clearly the mode shapes and natural frequencies, for instance, some of them are presented in [3,[5][6][7][8][9][10], besides, some extensive reviews can be found in [11][12][13][14].Some of these works are based on the fact that it is possible to estimate the crack parameters by measuring the changes in the natural frequency (inverse problem), since the natural frequencies can be measured accurately and most easily in comparison with other dynamic characteristics of a structure and the experimental errors do not significantly change their values [15,16].Before considering the crack identification inverse problem, the direct problem should be studied, that is, the determination of natural frequencies of the cracked element, assuming the crack parameters are known.
First of all, to establish the relationship between crack parameters and natural frequencies, assuming that, during the test period the shaft does not rotate, countless authors consider shafts as beams with a circular cross section, among them [17][18][19][20].According to Heydari et al. [21], one can find three approaches for the vibration modelling of cracked beams: continuous models with a continuous model for the crack, discrete models with a local flexibility model for the crack, and continuous models with a local flexibility model for the crack.With regard to the third approach, the continuous beam models with a local flexibility crack model many times consider two healthy segments that are linked by a massless rotational spring in order to take into account the flexibility due to the crack [5], so that the spring constant depends on the crack parameters by the fracture mechanics theory.is type of model has been successfully applied to Euler-Bernoulli cracked beams under different boundary conditions and rectangular cross sections [8,22,23], or circular cross sections [18,[24][25][26][27], but also concerning circular section beams or shafts, it is easy to find works that approach the problem by discrete methods [28][29][30][31].
It is known that the increasing flexibility caused by the crack depends on the geometry of the crack, so it is important to estimate in the best possible way the actual crack shape.In the case of cracked shafts, frequently, cracks grow following transversal planes to their longitudinal axis acquiring elliptical shape fronts [32][33][34][35][36][37].However, at present, in the knowledge of the authors, one cannot find in the literature any expression that allows calculating the relationship between crack parameters and natural frequencies in the case of elliptical cracks.In this paper, natural frequencies of bending vibrations of Euler-Bernoulli cracked beam with circular cross section, taking account the elliptical shape of the crack, have been obtained.
e natural frequencies of the damaged beams have been calculated by the perturbation analysis, a method that was developed to calculate the eigenvalue and eigenvector of the damaged structures for decreasing the computation expense [38].e perturbation techniques have been successfully used before by other authors to calculate the variation of natural frequencies induced by cracks in Euler-Bernoulli beams [8,[39][40][41].
Once the direct problem is known, the inverse problem must be approached.It is usually considered as an optimization problem where the crack parameters are regarded as design variables, so that they are estimated by minimizing the differences between the reference values and the ones predicted by the model.Several optimization methods have been used to solve the inverse problem for crack identification, for example, classical optimization technique of minimizing a least square criterion [8,42], the minimization of the mean-square and min-max formulations [43], or evolutionary algorithms [44].Among them, genetic algorithm (GA) optimization techniques have been successfully included in fault detection methods applied to mechanical systems, for instance, in [45][46][47][48][49][50].
In this paper, the perturbation theory has been used to obtain the analytical expressions of the natural frequencies of bending vibrations of nonrotating cracked shafts with circular cross section, taking into account location, size, and shape of the elliptical crack.
e crack is represented by a rotational spring connecting the two parts of the beam produced by the crack, whose constant is evaluated from the expression of flexibility coefficient found in [35].To solve the inverse problem of detecting and identifying the crack, the formulation previously obtained has been used in an optimization genetic algorithm which has been developed within the MatLab environment.Finally, in order to verify the GA, input data for the developed algorithm have been obtained from a modal 3D finite element (FE) model using Commercial Code ABAQUS [51].

Model of the Cracked Shaft Beam.
A cracked Euler-Bernoulli shaft beam that can vibrate in the X-Y plane is considered (Figure 1).It has length L and longitudinal uniform round cross section with diameter D. e beam contains a transversal crack with elliptical front characterized by semiaxes a and b, shown in Figure 1.It remains always open during vibration, and it is situated in an arbitrary location of abscissa x c .e crack characteristic parameters are as follows: (i) α � a/D characteristic depth of the crack (ii) c � a/b shape factor of the crack (c � 1 corresponds to a semicircular crack and c � 0 corresponds to a straight crack) (iii) δ � x c /L dimensionless location of the crack In order to take into account the presence of the crack, the shaft can be modelled with two beams connected by a massless torsional spring.e model leads to a discontinuity in the slope of the beam, Δθ, which is proportional to the bending moment transmitted by the crack section, Λ m (α, c) is the bending flexibility coefficient that depends on the crack characteristics, both depth and shape front.It is given by where E and v are Young's modulus and Poisson's ratio of the material and λ m (α, c) is the nondimensional flexibility coefficient for bending, which can be found in [35].e two beam segments can be treated separately.e vertical displacement of each part of the cracked beam connected by the spring is where subscripts j � 1 and j � 2 mean left and right parts of the cracked beam, respectively, u j (x) is the transverse deflection, and ω c is the frequency of vibration of the cracked shaft.Both can be calculated solving the following eigenvalue problem: where J is the moment of inertia, A is the area of the cross section of the beam, and ρ is the mass density of the material.e general solution for Equation (4) leads to a pair of equations one for each part of the cracked beam: 2 Shock and Vibration where the coordinate x is measured from the left support of the beam and e coefficients A k , B k , C k , and D k (k � 1, 2) can be found using the boundary conditions.In this work, pinnedpinned end condition has been considered.e mathematical relations describing the boundary restraints are as follows: Also, four conditions corresponding to the cracked section are used: (i) Transverse deflection continuity: (ii) Slope discontinuity: (iii) Bending moment continuity: (iv) Shear force continuity: In consonance with other authors [8,39,40], the problem has been solved using the perturbation method.According to this method, the solution for the cracked beam can be expressed from the uncracked beam solution as where ω nc is the frequency of the harmonic vibration of the uncracked shaft and u o (x) is its transverse deflection, ε is a positive parameter of the same order as Λ m and ϕ, and φ and μ are the variables of the problem that have to be obtained as a part of the solution.Substituting ( 13), ( 14), and ( 15) into ( 4) and keeping the first-order terms only, the natural frequencies of the cracked shaft beam can be derived from the natural frequencies of the uncracked one and the flexibility of the spring that models the crack: e displacement u o takes different expressions depending on the support conditions of the shaft beam [52].From (16), the relation between ω c and ω nc (frequency ratio) can be calculated in terms of D, L, and λ m and two parameters called g and f that depend on the support restrains: where β is a value that also depends on the end conditions of the beam [52].Table 1 shows β values for the first four natural frequencies of the shaft beam for pinned-pinned boundary condition.e expressions of g and f until the 4 th natural frequency are shown below (subscript i � 1, 2, 3, and 4 means the order of the natural frequency): As would be expected, the normalized natural frequencies depend on the position of the crack along the beam.If the crack is placed in a mode shape node, the natural frequency of the beam is not affected by the crack, that is, to say, ω * � 1, and also, the normalized natural frequencies decrease with the distance from the crack to the nearby nodes.

Effects of the Crack Depth and Crack Shape.
Figure 3 shows the variations of the first four normalized natural frequencies of the component in terms of crack depth (α) and position (δ) for c � 0.5, and Figure 4 shows the variations of the first four normalized natural frequencies in terms of crack front shape (c) and location (δ) for α � 0.25; the conclusions are also valid for other values of c and α, respectively.Both in Figures 3 and 4, axes have been oriented in order to achieve the best view of the represented surfaces.
With regard to the crack size influence in the frequency values, for all cases, the natural frequency of the cracked shaft decreases as crack is deeper.If one analyzes the variations of the frequency values depending on the crack shape front, it can be seen that, for all given frequencies, the natural frequencies of the element increase as crack shape is less straight, the largest values appear in the case of circular cracks (c � 1).

Genetic Algorithm.
Once the direct problem is known, the crack detection and identification can be considered as an inverse problem.e direct problem, in this case, has consisted in obtaining the first four natural frequencies from the knowledge of the location and the crack properties.e inverse problem will be then, obtaining the location and the flexibility of the cracked section, δ and λ m (α, c), knowing the natural frequencies.In the proposed methodology, the genetic algorithm (GA) optimization technique is applied to solve the inverse problem.e proposed GA is developed under the MatLab environment.To begin the GA, an initial population of 500 individuals are randomly generated.Every individual consists of 2 genes, λ m and δ, encoded as a binary string (chromosome).Among the different selection techniques, the roulette wheel method has been chosen, that gives every member a chance of being selected according to each individual's relative fitness.In relation to the crossover technique, the two-point operator has been employed; in this case, two points of the chromosome are randomly chosen in both individuals, and then the alleles falling in-between the two points are swapped to give two new offspring.On the contrary, the mutation method is used, which randomly creates adaptive rules with respect to the last successful or unsuccessful generation.At last, 85 is the percentage of the population to replace each generation.
e fitness function is the most important element of a optimization genetic algorithm.According to the GA toolbox of MatLab, the fitness function, f, should be minimized.It is computed as the root mean-square (RMS) of the differences between the actual first four dimensionless natural frequencies and the computed dimensionless frequencies as follows: where i denotes the order of the natural frequency and ω * i,a and ω * i,c represent the actual and the computed values of ratio frequencies, respectively.
To decide if the possible optimal values have been calculated, the best fitness value is checked in every step; if the fitness reaches a prefixed level or if the maximum number of iterations is achieved, the algorithm is stopped.As the stopping level, a value of 10 −6 has been chosen, and to guarantee that the algorithm converges to the optimal solution, each case has been run five times during 250 generations.

Numerical Results.
In order to test the proposed algorithm, the same shaft beam previously described has been considered, with cracks in different positions and with different depths and front shapes, that is to say, different flexibilities.e study has been carried out for three couples of α and c and three crack locations δ. e three couples (α, c) have been chosen considering the fact that the front of incipient cracks is practically circular and that when a crack grows its front becomes straight [53] (Figure 5).In summary, 9 cases have been considered.Table 2 shows the description of the considered cases.
Input data for the genetic algorithm have been obtained from a modal 3D finite element (FE) model using Commercial Code ABAQUS [51].
e mesh of the threedimensional model is made by employing 8-node linear brick elements called C3D8, and in order to determine the level of mesh refinement, convergence analysis was carried out.e FE model finally included about 45000 elements and 50000 nodes.In Figure 6, as an example, in the case of healthy beam, both the first four bending vibration modes related to the first four natural frequencies obtained and undeformed beams obtained from the finite element model can be seen.
In Table 3, the values of the reference frequency ratios for 9 cases calculated by the FE model can be seen.4

Shock and Vibration
It is important to take into account that it is possible to nd that cracked shaft beams with di erent depths and front shape crack present identical vibration behavior as long as the λ m value is the same in all cases.In Figure 7, it can be seen the surface that λ m represents values in function of α and c.In addition, Figure 8 shows several curves that have been calculated as the intersection of the aforementioned function and di erent planes characterized by constant λ m values.As shown, there are in nite values of couples (α, c) with which the same value of λ m can be obtained.
After applying the proposed GA, Table 4 shows the obtained results, where the estimated values of λ m and δ have been calculated as the average of the ve repetitions of the same case, and the standard deviation values are given as where x can be λ m or δ and n is the number of repetitions; in this case, n 5.
In all cases the error is calculated as error(%) |estimated value − true value| true value × 100.
Firstly, in the case of the crack location (δ), the estimation can be considered very good, since the error always keeps under 2.5% and the average error, taking into account all the studied cases, is less than 1%.Secondly, regarding the

Shock and Vibration
cracked section exibility (λ m ), one can nd errors ranging from 3.2% to 7.1%.With respect to exibility, as it is explained in previous sections, for a given λ m value, in nite couples (α, c) can be found, so every λ m value can be graphed as a curve composed by the corresponding in nite (α, c) points.Figure 9 shows, for the studied cases that can be found in Table 2, that given a shape front (c), for every λ m constant value, there is a di erence between the real dimensionless depth of the crack obtained from λ m real and the estimated crack size calculated from λ m predicted.In other words, one can see the exibility curves, both real (λ m real) and predicted (λ m predicted), as a function of their   9, the speci c points corresponding to the cases in Table 2 have also been highlighted.
In Table 5, one can see the calculated errors (ε α (%)), according to Equation (20), for the depth crack α corresponding to the previous analyzed cases and for the c limit values, and these are straight front (c 0) and circular front (c 1).
According to the results shown in Table 5, it can be said that the errors for α estimated values, in the limit cases of straight and circular fronts, are less than the errors associated with exibility λ m for the same shape fronts.Ultimately, the proposed algorithm allows estimating well enough the size of the crack for all the studied boundary conditions since in the worst situations, the errors are under 7%.

Conclusions
Analytical novel expressions available to calculate the rst four natural frequencies of an Euler-Bernoulli cracked shaft     8 Shock and Vibration beam have been presented in this work.ese expressions consider that, most of the times, the cracks presented in a shaft have elliptical fronts, so the natural frequencies values of the cracked element depend not only on the location and the size of the crack but also on the shape of the crack front.e method employed is based on linking the two undamaged parts of the shaft by a rotational spring, located in the cracked section, which constantly depends on the crack parameters, that is to say, on the flexibility due to the crack; meanwhile, the formulated problem has been solved by the perturbation technique.Results obtained from the developed expressions show the natural frequency of the cracked shaft decreases as crack is deeper and increases as crack shape is less straight.Once the analytical expressions have been formulated, the inverse problem has been approached.e genetic algorithm technique has been used as an optimization method in order to provide an algorithm that allows to estimate the crack parameters; these are the location and the flexibility due to the crack, which depends on the size and shape of the defect.A finite element model has also been created in order to use its results to check the validity of the proposed algorithm.In regards to the dimensionless crack section flexibility, which depends on the size and shape of the crack, the calculated errors are in a range of values between from 3.2% to 7.1%.However, because that it is possible to find that damaged shaft beams with different depths and front shape crack present identical vibration behavior as long as the λ m value is the same in all cases, the proposed algorithm allows to estimate well enough the size of the crack since in the worst situations, the errors are under 6.5%, and the average error is 2.3%.Last, in the case of the crack location, the estimation can be considered very good, because the error always keeps under 2.5%, and the average error, taking into account all the studied cases, is less than 1%.

Figure 2 :Figure 3 :
Figure 2: Variations of the rst four normalized frequencies in terms of crack position (δ) for α 0.3 and c 0.2.

Figure 4 :
Figure 4: Variations of the rst four normalized frequencies in terms of crack front shape (c) and position (δ) for α 0.25.

Figure 6 :
Figure 6: Both rst four bending vibration mode shapes and undeformed beams, for a healthy beam, obtained by the nite element model.

Figure 7 :Figure 8 :
Figure 7: λ m values in function of α and c.

Table 2 :
Input data: case description.

Table 3 :
Values of the reference frequency ratios.

Table 4 :
Location and exibility errors.